dcdfort
Modern Fortran library for analyzing DCD trajectories
Data Types | Functions/Subroutines | Variables
dcdfort_reader Module Reference

Data Types

type  dcdfile
 

Functions/Subroutines

subroutine dcdfile_open (this, filename)
 
subroutine dcdfile_read_header (this, nframes, istart, nevery, iend, timestep, natoms)
 
subroutine dcdfile_close (this)
 
subroutine dcdfile_read_next (this, xyz, box)
 
subroutine dcdfile_skip_next (this)
 

Variables

real(8), parameter pi = 2.0d0*dacos(0.0d0)
 

Function/Subroutine Documentation

◆ dcdfile_close()

subroutine dcdfort_reader::dcdfile_close ( class(dcdfile), intent(inout)  this)
143 
144  implicit none
145  class(dcdfile), intent(inout) :: this
146 
147  close(this%u)
148 

◆ dcdfile_open()

subroutine dcdfort_reader::dcdfile_open ( class(dcdfile), intent(inout)  this,
character (len=*)  filename 
)
44 
45  implicit none
46  character (len=*) :: filename
47  class(dcdfile), intent(inout) :: this
48  integer :: filesize
49  logical :: ex
50 
51  inquire(file=trim(filename), exist=ex, size=this%filesize)
52 
53  if (ex .eqv. .false.) then
54  call error_stop_program(trim(filename)//" does not exist.")
55  end if
56 
57  open(newunit=this%u, file=trim(filename), form="unformatted", access="stream")
58 

◆ dcdfile_read_header()

subroutine dcdfort_reader::dcdfile_read_header ( class(dcdfile), intent(inout)  this,
integer  nframes,
integer  istart,
integer  nevery,
integer  iend,
real  timestep,
integer  natoms 
)
62 
63  implicit none
64 
65  integer :: dummy, nframes, istart, nevery, iend, natoms, i, ntitle, n, framesize, nframes2
66  integer :: curr_pos
67  character (len=4) :: cord_string
68  character (len=80) :: title_string
69  real :: timestep
70  class(dcdfile), intent(inout) :: this
71 
72  read(this%u) dummy
73  if (dummy .ne. 84) then
74  call error_stop_program("This dcd file format is not supported, or the file header is corrupt.")
75  end if
76 
77  read(this%u) cord_string
78  if (cord_string .ne. "CORD") then
79  call error_stop_program("This dcd file format is not supported, or the file header is corrupt.")
80  end if
81 
82  ! Number of snapshots in file
83  read(this%u) nframes
84 
85  ! Timestep of first snapshot
86  read(this%u) istart
87 
88  ! Save snapshots every this many steps
89  read(this%u) nevery
90 
91  ! Timestep of last snapshot
92  read(this%u) iend
93 
94  do i = 1, 5
95  read(this%u) dummy
96  end do
97 
98  ! Simulation timestep
99  read(this%u) timestep
100 
101  do i = 1, 12
102  read(this%u) dummy
103  end do
104 
105  read(this%u) ntitle
106  do i = 1, ntitle
107  read(this%u) title_string
108 
109  n = 1
110  do while (n .le. 80 .and. title_string(n:n) .ne. c_null_char)
111  n = n + 1
112  end do
113 
114  write(error_unit,'(a)') trim(title_string(1:n))
115  end do
116 
117  read(this%u) dummy
118  read(this%u) dummy
119  if (dummy .ne. 4) then
120  call error_stop_program("This dcd file format is not supported, or the file header is corrupt.")
121  end if
122 
123  ! Number of atoms in each snapshot
124  read(this%u) natoms
125 
126  read(this%u) dummy
127 
128  inquire(unit=this%u, pos=curr_pos)
129 
130  ! Each frame has natoms*3 (4 bytes each), plus 6 box dimensions (8 bytes each)
131  framesize = natoms*12 + 48
132  ! Just an estimate
133  nframes2 = (this%filesize-curr_pos)/framesize
134  if ( nframes2 .ne. nframes) then
135  write(error_unit,'(a,i0,a,i0,a)') "WARNING: Header indicates ", nframes, &
136  &" frames, but file size indicates ", nframes2, "."
137  nframes = nframes2
138  end if
139 

◆ dcdfile_read_next()

subroutine dcdfort_reader::dcdfile_read_next ( class(dcdfile), intent(inout)  this,
real, dimension(:,:), intent(inout), allocatable  xyz,
real(8), dimension(6), intent(inout)  box 
)
152 
153  implicit none
154  real, allocatable, intent(inout) :: xyz(:,:)
155  real(8), intent(inout) :: box(6)
156  integer :: dummy
157  class(dcdfile), intent(inout) :: this
158 
159  ! Should be 48
160  read(this%u) dummy
161 
162  read(this%u) box(1) ! A
163  read(this%u) box(6) ! gamma
164  read(this%u) box(2) ! B
165  read(this%u) box(5) ! beta
166  read(this%u) box(4) ! alpha
167  read(this%u) box(3) ! C
168  if (box(4) >= -1.0 .and. box(4) <= 1.0 .and. box(5) >= -1.0 .and. box(5) <= 1.0 .and. &
169  box(6) >= -1.0 .and. box(6) <= 1.0) then
170  box(4) = 90.0 - asin(box(4)) * 180.0 / pi
171  box(5) = 90.0 - asin(box(5)) * 180.0 / pi
172  box(6) = 90.0 - asin(box(6)) * 180.0 / pi
173  end if
174 
175  ! 48 again
176  read(this%u) dummy
177  read(this%u) dummy
178 
179  read(this%u) xyz(1,:)
180 
181  ! 48 again
182  read(this%u) dummy
183  read(this%u) dummy
184 
185  read(this%u) xyz(2,:)
186 
187  ! 48 again
188  read(this%u) dummy
189  read(this%u) dummy
190 
191  read(this%u) xyz(3,:)
192 
193  read(this%u) dummy
194 

◆ dcdfile_skip_next()

subroutine dcdfort_reader::dcdfile_skip_next ( class(dcdfile), intent(inout)  this)
198 
199  implicit none
200  real :: real_dummy
201  integer :: dummy, i
202  real(8) :: box_dummy(6)
203  class(dcdfile), intent(inout) :: this
204 
205  ! Should be 48
206  read(this%u) dummy
207 
208  read(this%u) box_dummy
209 
210  ! 48 again
211  read(this%u) dummy
212  read(this%u) dummy
213 
214  do i = 1, dummy/4
215  read(this%u) real_dummy
216  end do
217 
218  ! 48 again
219  read(this%u) dummy
220  read(this%u) dummy
221 
222  do i = 1, dummy/4
223  read(this%u) real_dummy
224  end do
225 
226  ! 48 again
227  read(this%u) dummy
228  read(this%u) dummy
229 
230  do i = 1, dummy/4
231  read(this%u) real_dummy
232  end do
233 
234  read(this%u) dummy
235 

Variable Documentation

◆ pi

real(8), parameter dcdfort_reader::pi = 2.0d0*dacos(0.0d0)
28  real(8), parameter :: pi = 2.0d0*dacos(0.0d0)