File: | gromacs/gmxana/gmx_vanhove.c |
Location: | line 192, column 5 |
Description: | Value stored to 'natom' is never read |
1 | /* |
2 | * This file is part of the GROMACS molecular simulation package. |
3 | * |
4 | * Copyright (c) 1991-2000, University of Groningen, The Netherlands. |
5 | * Copyright (c) 2001-2004, The GROMACS development team. |
6 | * Copyright (c) 2013,2014, by the GROMACS development team, led by |
7 | * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, |
8 | * and including many others, as listed in the AUTHORS file in the |
9 | * top-level source directory and at http://www.gromacs.org. |
10 | * |
11 | * GROMACS is free software; you can redistribute it and/or |
12 | * modify it under the terms of the GNU Lesser General Public License |
13 | * as published by the Free Software Foundation; either version 2.1 |
14 | * of the License, or (at your option) any later version. |
15 | * |
16 | * GROMACS is distributed in the hope that it will be useful, |
17 | * but WITHOUT ANY WARRANTY; without even the implied warranty of |
18 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
19 | * Lesser General Public License for more details. |
20 | * |
21 | * You should have received a copy of the GNU Lesser General Public |
22 | * License along with GROMACS; if not, see |
23 | * http://www.gnu.org/licenses, or write to the Free Software Foundation, |
24 | * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. |
25 | * |
26 | * If you want to redistribute modifications to GROMACS, please |
27 | * consider that scientific software is very special. Version |
28 | * control is crucial - bugs must be traceable. We will be happy to |
29 | * consider code for inclusion in the official distribution, but |
30 | * derived work must not be called official GROMACS. Details are found |
31 | * in the README & COPYING files - if they are missing, get the |
32 | * official version at http://www.gromacs.org. |
33 | * |
34 | * To help us fund GROMACS development, we humbly ask that you cite |
35 | * the research papers on the package. Check out http://www.gromacs.org. |
36 | */ |
37 | #ifdef HAVE_CONFIG_H1 |
38 | #include <config.h> |
39 | #endif |
40 | |
41 | #include <math.h> |
42 | #include <stdlib.h> |
43 | #include <string.h> |
44 | |
45 | #include "gromacs/utility/smalloc.h" |
46 | #include "macros.h" |
47 | #include "gromacs/commandline/pargs.h" |
48 | #include "gromacs/math/utilities.h" |
49 | #include "gromacs/utility/futil.h" |
50 | #include "index.h" |
51 | #include "typedefs.h" |
52 | #include "gromacs/fileio/xvgr.h" |
53 | #include "viewit.h" |
54 | #include "gstat.h" |
55 | #include "gromacs/fileio/tpxio.h" |
56 | #include "gromacs/fileio/trxio.h" |
57 | #include "gromacs/math/vec.h" |
58 | #include "gromacs/fileio/matio.h" |
59 | #include "gmx_ana.h" |
60 | |
61 | |
62 | int gmx_vanhove(int argc, char *argv[]) |
63 | { |
64 | const char *desc[] = { |
65 | "[THISMODULE] computes the Van Hove correlation function.", |
66 | "The Van Hove G(r,t) is the probability that a particle that is at r[SUB]0[sub]", |
67 | "at time zero can be found at position r[SUB]0[sub]+r at time t.", |
68 | "[THISMODULE] determines G not for a vector r, but for the length of r.", |
69 | "Thus it gives the probability that a particle moves a distance of r", |
70 | "in time t.", |
71 | "Jumps across the periodic boundaries are removed.", |
72 | "Corrections are made for scaling due to isotropic", |
73 | "or anisotropic pressure coupling.", |
74 | "[PAR]", |
75 | "With option [TT]-om[tt] the whole matrix can be written as a function", |
76 | "of t and r or as a function of [SQRT]t[sqrt] and r (option [TT]-sqrt[tt]).", |
77 | "[PAR]", |
78 | "With option [TT]-or[tt] the Van Hove function is plotted for one", |
79 | "or more values of t. Option [TT]-nr[tt] sets the number of times,", |
80 | "option [TT]-fr[tt] the number spacing between the times.", |
81 | "The binwidth is set with option [TT]-rbin[tt]. The number of bins", |
82 | "is determined automatically.", |
83 | "[PAR]", |
84 | "With option [TT]-ot[tt] the integral up to a certain distance", |
85 | "(option [TT]-rt[tt]) is plotted as a function of time.", |
86 | "[PAR]", |
87 | "For all frames that are read the coordinates of the selected particles", |
88 | "are stored in memory. Therefore the program may use a lot of memory.", |
89 | "For options [TT]-om[tt] and [TT]-ot[tt] the program may be slow.", |
90 | "This is because the calculation scales as the number of frames times", |
91 | "[TT]-fm[tt] or [TT]-ft[tt].", |
92 | "Note that with the [TT]-dt[tt] option the memory usage and calculation", |
93 | "time can be reduced." |
94 | }; |
95 | static int fmmax = 0, ftmax = 0, nlev = 81, nr = 1, fshift = 0; |
96 | static real sbin = 0, rmax = 2, rbin = 0.01, mmax = 0, rint = 0; |
97 | t_pargs pa[] = { |
98 | { "-sqrt", FALSE0, etREAL, {&sbin}, |
99 | "Use [SQRT]t[sqrt] on the matrix axis which binspacing # in [SQRT]ps[sqrt]" }, |
100 | { "-fm", FALSE0, etINT, {&fmmax}, |
101 | "Number of frames in the matrix, 0 is plot all" }, |
102 | { "-rmax", FALSE0, etREAL, {&rmax}, |
103 | "Maximum r in the matrix (nm)" }, |
104 | { "-rbin", FALSE0, etREAL, {&rbin}, |
105 | "Binwidth in the matrix and for [TT]-or[tt] (nm)" }, |
106 | { "-mmax", FALSE0, etREAL, {&mmax}, |
107 | "Maximum density in the matrix, 0 is calculate (1/nm)" }, |
108 | { "-nlevels", FALSE0, etINT, {&nlev}, |
109 | "Number of levels in the matrix" }, |
110 | { "-nr", FALSE0, etINT, {&nr}, |
111 | "Number of curves for the [TT]-or[tt] output" }, |
112 | { "-fr", FALSE0, etINT, {&fshift}, |
113 | "Frame spacing for the [TT]-or[tt] output" }, |
114 | { "-rt", FALSE0, etREAL, {&rint}, |
115 | "Integration limit for the [TT]-ot[tt] output (nm)" }, |
116 | { "-ft", FALSE0, etINT, {&ftmax}, |
117 | "Number of frames in the [TT]-ot[tt] output, 0 is plot all" } |
118 | }; |
119 | #define NPA((int)(sizeof(pa)/sizeof((pa)[0]))) asize(pa)((int)(sizeof(pa)/sizeof((pa)[0]))) |
120 | |
121 | t_filenm fnm[] = { |
122 | { efTRX, NULL((void*)0), NULL((void*)0), ffREAD1<<1 }, |
123 | { efTPS, NULL((void*)0), NULL((void*)0), ffREAD1<<1 }, |
124 | { efNDX, NULL((void*)0), NULL((void*)0), ffOPTRD(1<<1 | 1<<3) }, |
125 | { efXPM, "-om", "vanhove", ffOPTWR(1<<2| 1<<3) }, |
126 | { efXVG, "-or", "vanhove_r", ffOPTWR(1<<2| 1<<3) }, |
127 | { efXVG, "-ot", "vanhove_t", ffOPTWR(1<<2| 1<<3) } |
128 | }; |
129 | #define NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))) asize(fnm)((int)(sizeof(fnm)/sizeof((fnm)[0]))) |
130 | |
131 | output_env_t oenv; |
132 | const char *matfile, *otfile, *orfile; |
133 | char title[256]; |
134 | t_topology top; |
135 | int ePBC; |
136 | matrix boxtop, box, *sbox, avbox, corr; |
137 | rvec *xtop, *x, **sx; |
138 | int isize, nalloc, nallocn, natom; |
139 | t_trxstatus *status; |
140 | atom_id *index; |
141 | char *grpname; |
142 | int nfr, f, ff, i, m, mat_nx = 0, nbin = 0, bin, mbin, fbin; |
143 | real *time, t, invbin = 0, rmax2 = 0, rint2 = 0, d2; |
144 | real invsbin = 0, matmax, normfac, dt, *tickx, *ticky; |
145 | char buf[STRLEN4096], **legend; |
146 | real **mat = NULL((void*)0); |
147 | int *pt = NULL((void*)0), **pr = NULL((void*)0), *mcount = NULL((void*)0), *tcount = NULL((void*)0), *rcount = NULL((void*)0); |
148 | FILE *fp; |
149 | t_rgb rlo = {1, 1, 1}, rhi = {0, 0, 0}; |
150 | |
151 | if (!parse_common_args(&argc, argv, PCA_CAN_VIEW(1<<5) | PCA_CAN_TIME((1<<6) | (1<<7) | (1<<14)) | PCA_BE_NICE(1<<13), |
152 | NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm, asize(pa)((int)(sizeof(pa)/sizeof((pa)[0]))), pa, asize(desc)((int)(sizeof(desc)/sizeof((desc)[0]))), desc, 0, NULL((void*)0), &oenv)) |
153 | { |
154 | return 0; |
155 | } |
156 | |
157 | matfile = opt2fn_null("-om", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm); |
158 | if (opt2parg_bSet("-fr", NPA((int)(sizeof(pa)/sizeof((pa)[0]))), pa)) |
159 | { |
160 | orfile = opt2fn("-or", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm); |
161 | } |
162 | else |
163 | { |
164 | orfile = opt2fn_null("-or", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm); |
165 | } |
166 | if (opt2parg_bSet("-rt", NPA((int)(sizeof(pa)/sizeof((pa)[0]))), pa)) |
167 | { |
168 | otfile = opt2fn("-ot", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm); |
169 | } |
170 | else |
171 | { |
172 | otfile = opt2fn_null("-ot", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm); |
173 | } |
174 | |
175 | if (!matfile && !otfile && !orfile) |
176 | { |
177 | fprintf(stderrstderr, |
178 | "For output set one (or more) of the output file options\n"); |
179 | exit(0); |
180 | } |
181 | |
182 | read_tps_conf(ftp2fn(efTPS, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), title, &top, &ePBC, &xtop, NULL((void*)0), boxtop, |
183 | FALSE0); |
184 | get_index(&top.atoms, ftp2fn_null(efNDX, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), 1, &isize, &index, &grpname); |
185 | |
186 | nalloc = 0; |
187 | time = NULL((void*)0); |
188 | sbox = NULL((void*)0); |
189 | sx = NULL((void*)0); |
190 | clear_mat(avbox); |
191 | |
192 | natom = read_first_x(oenv, &status, ftp2fn(efTRX, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), &t, &x, box); |
Value stored to 'natom' is never read | |
193 | nfr = 0; |
194 | do |
195 | { |
196 | if (nfr >= nalloc) |
197 | { |
198 | nalloc += 100; |
199 | srenew(time, nalloc)(time) = save_realloc("time", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 199, (time), (nalloc), sizeof(*(time))); |
200 | srenew(sbox, nalloc)(sbox) = save_realloc("sbox", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 200, (sbox), (nalloc), sizeof(*(sbox))); |
201 | srenew(sx, nalloc)(sx) = save_realloc("sx", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 201, (sx), (nalloc), sizeof(*(sx))); |
202 | } |
203 | |
204 | time[nfr] = t; |
205 | copy_mat(box, sbox[nfr]); |
206 | /* This assumes that the off-diagonal box elements |
207 | * are not affected by jumps across the periodic boundaries. |
208 | */ |
209 | m_add(avbox, box, avbox); |
210 | snew(sx[nfr], isize)(sx[nfr]) = save_calloc("sx[nfr]", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 210, (isize), sizeof(*(sx[nfr]))); |
211 | for (i = 0; i < isize; i++) |
212 | { |
213 | copy_rvec(x[index[i]], sx[nfr][i]); |
214 | } |
215 | |
216 | nfr++; |
217 | } |
218 | while (read_next_x(oenv, status, &t, x, box)); |
219 | |
220 | /* clean up */ |
221 | sfree(x)save_free("x", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 221, (x)); |
222 | close_trj(status); |
223 | |
224 | fprintf(stderrstderr, "Read %d frames\n", nfr); |
225 | |
226 | dt = (time[nfr-1] - time[0])/(nfr - 1); |
227 | /* Some ugly rounding to get nice nice times in the output */ |
228 | dt = (int)(10000.0*dt + 0.5)/10000.0; |
229 | |
230 | invbin = 1.0/rbin; |
231 | |
232 | if (matfile) |
233 | { |
234 | if (fmmax <= 0 || fmmax >= nfr) |
235 | { |
236 | fmmax = nfr - 1; |
237 | } |
238 | snew(mcount, fmmax)(mcount) = save_calloc("mcount", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 238, (fmmax), sizeof(*(mcount))); |
239 | nbin = (int)(rmax*invbin + 0.5); |
240 | if (sbin == 0) |
241 | { |
242 | mat_nx = fmmax + 1; |
243 | } |
244 | else |
245 | { |
246 | invsbin = 1.0/sbin; |
247 | mat_nx = sqrt(fmmax*dt)*invsbin + 1; |
248 | } |
249 | snew(mat, mat_nx)(mat) = save_calloc("mat", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 249, (mat_nx), sizeof(*(mat))); |
250 | for (f = 0; f < mat_nx; f++) |
251 | { |
252 | snew(mat[f], nbin)(mat[f]) = save_calloc("mat[f]", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 252, (nbin), sizeof(*(mat[f]))); |
253 | } |
254 | rmax2 = sqr(nbin*rbin); |
255 | /* Initialize time zero */ |
256 | mat[0][0] = nfr*isize; |
257 | mcount[0] += nfr; |
258 | } |
259 | else |
260 | { |
261 | fmmax = 0; |
262 | } |
263 | |
264 | if (orfile) |
265 | { |
266 | snew(pr, nr)(pr) = save_calloc("pr", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 266, (nr), sizeof(*(pr))); |
267 | nalloc = 0; |
268 | snew(rcount, nr)(rcount) = save_calloc("rcount", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 268, (nr), sizeof(*(rcount))); |
269 | } |
270 | |
271 | if (otfile) |
272 | { |
273 | if (ftmax <= 0) |
274 | { |
275 | ftmax = nfr - 1; |
276 | } |
277 | snew(tcount, ftmax)(tcount) = save_calloc("tcount", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 277, (ftmax), sizeof(*(tcount))); |
278 | snew(pt, nfr)(pt) = save_calloc("pt", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 278, (nfr), sizeof(*(pt))); |
279 | rint2 = rint*rint; |
280 | /* Initialize time zero */ |
281 | pt[0] = nfr*isize; |
282 | tcount[0] += nfr; |
283 | } |
284 | else |
285 | { |
286 | ftmax = 0; |
287 | } |
288 | |
289 | msmul(avbox, 1.0/nfr, avbox); |
290 | for (f = 0; f < nfr; f++) |
291 | { |
292 | if (f % 100 == 0) |
293 | { |
294 | fprintf(stderrstderr, "\rProcessing frame %d", f); |
295 | } |
296 | /* Scale all the configuration to the average box */ |
297 | m_inv_ur0(sbox[f], corr); |
298 | mmul_ur0(avbox, corr, corr); |
299 | for (i = 0; i < isize; i++) |
300 | { |
301 | mvmul_ur0(corr, sx[f][i], sx[f][i]); |
302 | if (f > 0) |
303 | { |
304 | /* Correct for periodic jumps */ |
305 | for (m = DIM3-1; m >= 0; m--) |
306 | { |
307 | while (sx[f][i][m] - sx[f-1][i][m] > 0.5*avbox[m][m]) |
308 | { |
309 | rvec_dec(sx[f][i], avbox[m]); |
310 | } |
311 | while (sx[f][i][m] - sx[f-1][i][m] <= -0.5*avbox[m][m]) |
312 | { |
313 | rvec_inc(sx[f][i], avbox[m]); |
314 | } |
315 | } |
316 | } |
317 | } |
318 | for (ff = 0; ff < f; ff++) |
319 | { |
320 | fbin = f - ff; |
321 | if (fbin <= fmmax || fbin <= ftmax) |
322 | { |
323 | if (sbin == 0) |
324 | { |
325 | mbin = fbin; |
326 | } |
327 | else |
328 | { |
329 | mbin = (int)(sqrt(fbin*dt)*invsbin + 0.5); |
330 | } |
331 | for (i = 0; i < isize; i++) |
332 | { |
333 | d2 = distance2(sx[f][i], sx[ff][i]); |
334 | if (mbin < mat_nx && d2 < rmax2) |
335 | { |
336 | bin = (int)(sqrt(d2)*invbin + 0.5); |
337 | if (bin < nbin) |
338 | { |
339 | mat[mbin][bin] += 1; |
340 | } |
341 | } |
342 | if (fbin <= ftmax && d2 <= rint2) |
343 | { |
344 | pt[fbin]++; |
345 | } |
346 | } |
347 | if (matfile) |
348 | { |
349 | mcount[mbin]++; |
350 | } |
351 | if (otfile) |
352 | { |
353 | tcount[fbin]++; |
354 | } |
355 | } |
356 | } |
357 | if (orfile) |
358 | { |
359 | for (fbin = 0; fbin < nr; fbin++) |
360 | { |
361 | ff = f - (fbin + 1)*fshift; |
362 | if (ff >= 0) |
363 | { |
364 | for (i = 0; i < isize; i++) |
365 | { |
366 | d2 = distance2(sx[f][i], sx[ff][i]); |
367 | bin = (int)(sqrt(d2)*invbin + 0.5); |
368 | if (bin >= nalloc) |
369 | { |
370 | nallocn = 10*(bin/10) + 11; |
371 | for (m = 0; m < nr; m++) |
372 | { |
373 | srenew(pr[m], nallocn)(pr[m]) = save_realloc("pr[m]", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 373, (pr[m]), (nallocn), sizeof(*(pr[m]))); |
374 | for (i = nalloc; i < nallocn; i++) |
375 | { |
376 | pr[m][i] = 0; |
377 | } |
378 | } |
379 | nalloc = nallocn; |
380 | } |
381 | pr[fbin][bin]++; |
382 | } |
383 | rcount[fbin]++; |
384 | } |
385 | } |
386 | } |
387 | } |
388 | fprintf(stderrstderr, "\n"); |
389 | |
390 | if (matfile) |
391 | { |
392 | matmax = 0; |
393 | for (f = 0; f < mat_nx; f++) |
394 | { |
395 | normfac = 1.0/(mcount[f]*isize*rbin); |
396 | for (i = 0; i < nbin; i++) |
397 | { |
398 | mat[f][i] *= normfac; |
399 | if (mat[f][i] > matmax && (f != 0 || i != 0)) |
400 | { |
401 | matmax = mat[f][i]; |
402 | } |
403 | } |
404 | } |
405 | fprintf(stdoutstdout, "Value at (0,0): %.3f, maximum of the rest %.3f\n", |
406 | mat[0][0], matmax); |
407 | if (mmax > 0) |
408 | { |
409 | matmax = mmax; |
410 | } |
411 | snew(tickx, mat_nx)(tickx) = save_calloc("tickx", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 411, (mat_nx), sizeof(*(tickx))); |
412 | for (f = 0; f < mat_nx; f++) |
413 | { |
414 | if (sbin == 0) |
415 | { |
416 | tickx[f] = f*dt; |
417 | } |
418 | else |
419 | { |
420 | tickx[f] = f*sbin; |
421 | } |
422 | } |
423 | snew(ticky, nbin+1)(ticky) = save_calloc("ticky", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 423, (nbin+1), sizeof(*(ticky))); |
424 | for (i = 0; i <= nbin; i++) |
425 | { |
426 | ticky[i] = i*rbin; |
427 | } |
428 | fp = gmx_ffopen(matfile, "w"); |
429 | write_xpm(fp, MAT_SPATIAL_Y(1<<1), "Van Hove function", "G (1/nm)", |
430 | sbin == 0 ? "time (ps)" : "sqrt(time) (ps^1/2)", "r (nm)", |
431 | mat_nx, nbin, tickx, ticky, mat, 0, matmax, rlo, rhi, &nlev); |
432 | gmx_ffclose(fp); |
433 | } |
434 | |
435 | if (orfile) |
436 | { |
437 | fp = xvgropen(orfile, "Van Hove function", "r (nm)", "G (nm\\S-1\\N)", oenv); |
438 | fprintf(fp, "@ subtitle \"for particles in group %s\"\n", grpname); |
439 | snew(legend, nr)(legend) = save_calloc("legend", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_vanhove.c" , 439, (nr), sizeof(*(legend))); |
440 | for (fbin = 0; fbin < nr; fbin++) |
441 | { |
442 | sprintf(buf, "%g ps", (fbin + 1)*fshift*dt); |
443 | legend[fbin] = strdup(buf)(__extension__ (__builtin_constant_p (buf) && ((size_t )(const void *)((buf) + 1) - (size_t)(const void *)(buf) == 1 ) ? (((const char *) (buf))[0] == '\0' ? (char *) calloc ((size_t ) 1, (size_t) 1) : ({ size_t __len = strlen (buf) + 1; char * __retval = (char *) malloc (__len); if (__retval != ((void*)0 )) __retval = (char *) memcpy (__retval, buf, __len); __retval ; })) : __strdup (buf))); |
444 | } |
445 | xvgr_legend(fp, nr, (const char**)legend, oenv); |
446 | for (i = 0; i < nalloc; i++) |
447 | { |
448 | fprintf(fp, "%g", i*rbin); |
449 | for (fbin = 0; fbin < nr; fbin++) |
450 | { |
451 | fprintf(fp, " %g", |
452 | (real)pr[fbin][i]/(rcount[fbin]*isize*rbin*(i == 0 ? 0.5 : 1))); |
453 | } |
454 | fprintf(fp, "\n"); |
455 | } |
456 | gmx_ffclose(fp); |
457 | } |
458 | |
459 | if (otfile) |
460 | { |
461 | sprintf(buf, "Probability of moving less than %g nm", rint); |
462 | fp = xvgropen(otfile, buf, "t (ps)", "", oenv); |
463 | fprintf(fp, "@ subtitle \"for particles in group %s\"\n", grpname); |
464 | for (f = 0; f <= ftmax; f++) |
465 | { |
466 | fprintf(fp, "%g %g\n", f*dt, (real)pt[f]/(tcount[f]*isize)); |
467 | } |
468 | gmx_ffclose(fp); |
469 | } |
470 | |
471 | do_view(oenv, matfile, NULL((void*)0)); |
472 | do_view(oenv, orfile, NULL((void*)0)); |
473 | do_view(oenv, otfile, NULL((void*)0)); |
474 | |
475 | return 0; |
476 | } |