2 * This file is part of the GROMACS molecular simulation package.
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.
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.
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.
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.
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.
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.
49 #include "gromacs/fileio/confio.h"
53 #include "gromacs/commandline/pargs.h"
57 #include "gromacs/fileio/pdbio.h"
60 #include "gromacs/fileio/tpxio.h"
61 #include "gromacs/fileio/trxio.h"
64 #include "gromacs/fileio/matio.h"
65 #include "mtop_util.h"
79 real sumv, averv, maxv;
80 real *aver1, *aver2, *aver_3, *aver_6;
83 static void init5(int n)
89 static void reset5(void)
93 for (i = 0; (i < ntop); i++)
100 int tpcomp(const void *a, const void *b)
108 return (1e7*(tpb->v-tpa->v));
111 static void add5(int ndr, real viol)
116 for (i = 1; (i < ntop); i++)
118 if (top[i].v < top[mini].v)
123 if (viol > top[mini].v)
130 static void print5(FILE *fp)
134 qsort(top, ntop, sizeof(top[0]), tpcomp);
135 fprintf(fp, "Index:");
136 for (i = 0; (i < ntop); i++)
138 fprintf(fp, " %6d", top[i].n);
140 fprintf(fp, "\nViol: ");
141 for (i = 0; (i < ntop); i++)
143 fprintf(fp, " %6.3f", top[i].v);
148 static void check_viol(FILE *log,
149 t_ilist *disres, t_iparams forceparams[],
151 t_pbc *pbc, t_graph *g, t_dr_result dr[],
152 int clust_id, int isize, atom_id index[], real vvindex[],
156 int i, j, nat, n, type, nviol, ndr, label;
157 real ener, rt, mviol, tviol, viol, lam, dvdl, drt;
159 static gmx_bool bFirst = TRUE;
171 forceatoms = disres->iatoms;
172 for (j = 0; (j < isize); j++)
176 nat = interaction_function[F_DISRES].nratoms+1;
177 for (i = 0; (i < disres->nr); )
179 type = forceatoms[i];
181 label = forceparams[type].disres.label;
184 fprintf(debug, "DISRE: ndr = %d, label = %d i=%d, n =%d\n",
189 gmx_fatal(FARGS, "tpr inconsistency. ndr = %d, label = %d\n", ndr, label);
195 while (((i+n) < disres->nr) &&
196 (forceparams[forceatoms[i+n]].disres.label == label));
198 calc_disres_R_6(n, &forceatoms[i], forceparams,
199 (const rvec*)x, pbc, fcd, NULL);
201 if (fcd->disres.Rt_6[0] <= 0)
203 gmx_fatal(FARGS, "ndr = %d, rt_6 = %f", ndr, fcd->disres.Rt_6[0]);
206 rt = pow(fcd->disres.Rt_6[0], -1.0/6.0);
207 dr[clust_id].aver1[ndr] += rt;
208 dr[clust_id].aver2[ndr] += sqr(rt);
210 dr[clust_id].aver_3[ndr] += drt;
211 dr[clust_id].aver_6[ndr] += fcd->disres.Rt_6[0];
213 snew(fshift, SHIFTS);
214 ener = interaction_function[F_DISRES].ifunc(n, &forceatoms[i], forceparams,
215 (const rvec*)x, f, fshift,
216 pbc, g, lam, &dvdl, NULL, fcd, NULL);
218 viol = fcd->disres.sumviol;
225 add5(forceparams[type].disres.label, viol);
232 for (j = 0; (j < isize); j++)
234 if (index[j] == forceparams[type].disres.label)
236 vvindex[j] = pow(fcd->disres.Rt_6[0], -1.0/6.0);
243 dr[clust_id].nv = nviol;
244 dr[clust_id].maxv = mviol;
245 dr[clust_id].sumv = tviol;
246 dr[clust_id].averv = tviol/ndr;
247 dr[clust_id].nframes++;
251 fprintf(stderr, "\nThere are %d restraints and %d pairs\n",
252 ndr, disres->nr/nat);
264 real up1, r, rT3, rT6, viol, violT3, violT6;
267 static int drs_comp(const void *a, const void *b)
271 da = (t_dr_stats *)a;
272 db = (t_dr_stats *)b;
274 if (da->viol > db->viol)
278 else if (da->viol < db->viol)
288 static void dump_dump(FILE *log, int ndr, t_dr_stats drs[])
290 static const char *core[] = { "All restraints", "Core restraints" };
291 static const char *tp[] = { "linear", "third power", "sixth power" };
292 real viol_tot, viol_max, viol = 0;
297 for (bCore = FALSE; (bCore <= TRUE); bCore++)
299 for (kkk = 0; (kkk < 3); kkk++)
305 for (i = 0; (i < ndr); i++)
307 if (!bCore || (bCore && drs[i].bCore))
315 viol = drs[i].violT3;
318 viol = drs[i].violT6;
321 gmx_incons("Dumping violations");
323 viol_max = max(viol_max, viol);
332 if ((nrestr > 0) || (bCore && (nrestr < ndr)))
335 fprintf(log, "+++++++ %s ++++++++\n", core[bCore]);
336 fprintf(log, "+++++++ Using %s averaging: ++++++++\n", tp[kkk]);
337 fprintf(log, "Sum of violations: %8.3f nm\n", viol_tot);
340 fprintf(log, "Average violation: %8.3f nm\n", viol_tot/nrestr);
342 fprintf(log, "Largest violation: %8.3f nm\n", viol_max);
343 fprintf(log, "Number of violated restraints: %d/%d\n", nviol, nrestr);
349 static void dump_viol(FILE *log, int ndr, t_dr_stats *drs, gmx_bool bLinear)
353 fprintf(log, "Restr. Core Up1 <r> <rT3> <rT6> <viol><violT3><violT6>\n");
354 for (i = 0; (i < ndr); i++)
356 if (bLinear && (drs[i].viol == 0))
360 fprintf(log, "%6d%5s%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f\n",
361 drs[i].index, yesno_names[drs[i].bCore],
362 drs[i].up1, drs[i].r, drs[i].rT3, drs[i].rT6,
363 drs[i].viol, drs[i].violT3, drs[i].violT6);
367 static gmx_bool is_core(int i, int isize, atom_id index[])
370 gmx_bool bIC = FALSE;
372 for (kk = 0; !bIC && (kk < isize); kk++)
374 bIC = (index[kk] == i);
380 static void dump_stats(FILE *log, int nsteps, int ndr, t_ilist *disres,
381 t_iparams ip[], t_dr_result *dr,
382 int isize, atom_id index[], t_atoms *atoms)
388 fprintf(log, "++++++++++++++ STATISTICS ++++++++++++++++++++++++\n");
391 nra = interaction_function[F_DISRES].nratoms+1;
392 for (i = 0; (i < ndr); i++)
394 /* Search for the appropriate restraint */
395 for (; (j < disres->nr); j += nra)
397 if (ip[disres->iatoms[j]].disres.label == i)
403 drs[i].bCore = is_core(i, isize, index);
404 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
405 drs[i].r = dr->aver1[i]/nsteps;
406 drs[i].rT3 = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
407 drs[i].rT6 = pow(dr->aver_6[i]/nsteps, -1.0/6.0);
408 drs[i].viol = max(0, drs[i].r-drs[i].up1);
409 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
410 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
413 int j1 = disres->iatoms[j+1];
414 int j2 = disres->iatoms[j+2];
415 atoms->pdbinfo[j1].bfac += drs[i].violT3*5;
416 atoms->pdbinfo[j2].bfac += drs[i].violT3*5;
419 dump_viol(log, ndr, drs, FALSE);
421 fprintf(log, "+++ Sorted by linear averaged violations: +++\n");
422 qsort(drs, ndr, sizeof(drs[0]), drs_comp);
423 dump_viol(log, ndr, drs, TRUE);
425 dump_dump(log, ndr, drs);
430 static void dump_clust_stats(FILE *fp, int ndr, t_ilist *disres,
431 t_iparams ip[], t_blocka *clust, t_dr_result dr[],
432 char *clust_name[], int isize, atom_id index[])
434 int i, j, k, nra, mmm = 0;
435 double sumV, maxV, sumVT3, sumVT6, maxVT3, maxVT6;
439 fprintf(fp, "++++++++++++++ STATISTICS ++++++++++++++++++++++\n");
440 fprintf(fp, "Cluster NFrames SumV MaxV SumVT MaxVT SumVS MaxVS\n");
444 for (k = 0; (k < clust->nr); k++)
446 if (dr[k].nframes == 0)
450 if (dr[k].nframes != (clust->index[k+1]-clust->index[k]))
452 gmx_fatal(FARGS, "Inconsistency in cluster %s.\n"
453 "Found %d frames in trajectory rather than the expected %d\n",
454 clust_name[k], dr[k].nframes,
455 clust->index[k+1]-clust->index[k]);
459 gmx_fatal(FARGS, "Inconsistency with cluster %d. Invalid name", k);
462 nra = interaction_function[F_DISRES].nratoms+1;
463 sumV = sumVT3 = sumVT6 = maxV = maxVT3 = maxVT6 = 0;
464 for (i = 0; (i < ndr); i++)
466 /* Search for the appropriate restraint */
467 for (; (j < disres->nr); j += nra)
469 if (ip[disres->iatoms[j]].disres.label == i)
475 drs[i].bCore = is_core(i, isize, index);
476 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
477 drs[i].r = dr[k].aver1[i]/dr[k].nframes;
478 if ((dr[k].aver_3[i] <= 0) || (dr[k].aver_3[i] != dr[k].aver_3[i]))
480 gmx_fatal(FARGS, "dr[%d].aver_3[%d] = %f", k, i, dr[k].aver_3[i]);
482 drs[i].rT3 = pow(dr[k].aver_3[i]/dr[k].nframes, -1.0/3.0);
483 drs[i].rT6 = pow(dr[k].aver_6[i]/dr[k].nframes, -1.0/6.0);
484 drs[i].viol = max(0, drs[i].r-drs[i].up1);
485 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
486 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
488 sumVT3 += drs[i].violT3;
489 sumVT6 += drs[i].violT6;
490 maxV = max(maxV, drs[i].viol);
491 maxVT3 = max(maxVT3, drs[i].violT3);
492 maxVT6 = max(maxVT6, drs[i].violT6);
494 if (strcmp(clust_name[k], "1000") == 0)
498 fprintf(fp, "%-10s%6d%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n",
500 dr[k].nframes, sumV, maxV, sumVT3, maxVT3, sumVT6, maxVT6);
507 static void init_dr_res(t_dr_result *dr, int ndr)
509 snew(dr->aver1, ndr+1);
510 snew(dr->aver2, ndr+1);
511 snew(dr->aver_3, ndr+1);
512 snew(dr->aver_6, ndr+1);
520 static void dump_disre_matrix(const char *fn, t_dr_result *dr, int ndr,
521 int nsteps, t_idef *idef, gmx_mtop_t *mtop,
522 real max_dr, int nlevels, gmx_bool bThird)
526 int n_res, a_offset, mb, mol, a;
528 int iii, i, j, nra, nratoms, tp, ri, rj, index, nlabel, label;
529 atom_id ai, aj, *ptr;
530 real **matrix, *t_res, hi, *w_dr, rav, rviol;
531 t_rgb rlo = { 1, 1, 1 };
532 t_rgb rhi = { 0, 0, 0 };
538 snew(resnr, mtop->natoms);
541 for (mb = 0; mb < mtop->nmolblock; mb++)
543 atoms = &mtop->moltype[mtop->molblock[mb].type].atoms;
544 for (mol = 0; mol < mtop->molblock[mb].nmol; mol++)
546 for (a = 0; a < atoms->nr; a++)
548 resnr[a_offset+a] = n_res + atoms->atom[a].resind;
550 n_res += atoms->nres;
551 a_offset += atoms->nr;
556 for (i = 0; (i < n_res); i++)
561 for (i = 0; (i < n_res); i++)
563 snew(matrix[i], n_res);
565 nratoms = interaction_function[F_DISRES].nratoms;
566 nra = (idef->il[F_DISRES].nr/(nratoms+1));
572 for (i = 0; (i < idef->il[F_DISRES].nr); i += nratoms+1)
574 tp = idef->il[F_DISRES].iatoms[i];
575 label = idef->iparams[tp].disres.label;
579 /* Set index pointer */
583 gmx_fatal(FARGS, "nlabel is %d, label = %d", nlabel, label);
587 gmx_fatal(FARGS, "ndr = %d, index = %d", ndr, index);
589 /* Update the weight */
590 w_dr[index] = 1.0/nlabel;
599 printf("nlabel = %d, index = %d, ndr = %d\n", nlabel, index, ndr);
601 for (i = 0; (i < ndr); i++)
603 for (j = ptr[i]; (j < ptr[i+1]); j += nratoms+1)
605 tp = idef->il[F_DISRES].iatoms[j];
606 ai = idef->il[F_DISRES].iatoms[j+1];
607 aj = idef->il[F_DISRES].iatoms[j+2];
613 rav = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
617 rav = dr->aver1[i]/nsteps;
621 fprintf(debug, "DR %d, atoms %d, %d, distance %g\n", i, ai, aj, rav);
623 rviol = max(0, rav-idef->iparams[tp].disres.up1);
624 matrix[ri][rj] += w_dr[i]*rviol;
625 matrix[rj][ri] += w_dr[i]*rviol;
626 hi = max(hi, matrix[ri][rj]);
627 hi = max(hi, matrix[rj][ri]);
637 printf("Warning: the maxdr that you have specified (%g) is smaller than\nthe largest value in your simulation (%g)\n", max_dr, hi);
641 printf("Highest level in the matrix will be %g\n", hi);
642 fp = gmx_ffopen(fn, "w");
643 write_xpm(fp, 0, "Distance Violations", "<V> (nm)", "Residue", "Residue",
644 n_res, n_res, t_res, t_res, matrix, 0, hi, rlo, rhi, &nlevels);
648 int gmx_disre(int argc, char *argv[])
650 const char *desc[] = {
651 "[THISMODULE] computes violations of distance restraints.",
652 "If necessary, all protons can be added to a protein molecule ",
653 "using the [gmx-protonate] program.[PAR]",
654 "The program always",
655 "computes the instantaneous violations rather than time-averaged,",
656 "because this analysis is done from a trajectory file afterwards",
657 "it does not make sense to use time averaging. However,",
658 "the time averaged values per restraint are given in the log file.[PAR]",
659 "An index file may be used to select specific restraints for",
661 "When the optional [TT]-q[tt] flag is given a [TT].pdb[tt] file coloured by the",
662 "amount of average violations.[PAR]",
663 "When the [TT]-c[tt] option is given, an index file will be read",
664 "containing the frames in your trajectory corresponding to the clusters",
665 "(defined in another manner) that you want to analyze. For these clusters",
666 "the program will compute average violations using the third power",
667 "averaging algorithm and print them in the log file."
670 static int nlevels = 20;
671 static real max_dr = 0;
672 static gmx_bool bThird = TRUE;
674 { "-ntop", FALSE, etINT, {&ntop},
675 "Number of large violations that are stored in the log file every step" },
676 { "-maxdr", FALSE, etREAL, {&max_dr},
677 "Maximum distance violation in matrix output. If less than or equal to 0 the maximum will be determined by the data." },
678 { "-nlevels", FALSE, etINT, {&nlevels},
679 "Number of levels in the matrix output" },
680 { "-third", FALSE, etBOOL, {&bThird},
681 "Use inverse third power averaging or linear for matrix output" }
684 FILE *out = NULL, *aver = NULL, *numv = NULL, *maxxv = NULL, *xvg = NULL;
690 t_atoms *atoms = NULL;
694 int ntopatoms, natoms, i, j, kkk;
697 rvec *x, *f, *xav = NULL;
701 atom_id *index = NULL, *ind_fit = NULL;
703 t_cluster_ndx *clust = NULL;
704 t_dr_result dr, *dr_clust = NULL;
706 real *vvindex = NULL, *w_rls = NULL;
708 t_pbc pbc, *pbc_null;
712 gmx_rmpbc_t gpbc = NULL;
715 { efTPX, NULL, NULL, ffREAD },
716 { efTRX, "-f", NULL, ffREAD },
717 { efXVG, "-ds", "drsum", ffWRITE },
718 { efXVG, "-da", "draver", ffWRITE },
719 { efXVG, "-dn", "drnum", ffWRITE },
720 { efXVG, "-dm", "drmax", ffWRITE },
721 { efXVG, "-dr", "restr", ffWRITE },
722 { efLOG, "-l", "disres", ffWRITE },
723 { efNDX, NULL, "viol", ffOPTRD },
724 { efPDB, "-q", "viol", ffOPTWR },
725 { efNDX, "-c", "clust", ffOPTRD },
726 { efXPM, "-x", "matrix", ffOPTWR }
728 #define NFILE asize(fnm)
730 if (!parse_common_args(&argc, argv, PCA_CAN_TIME | PCA_CAN_VIEW | PCA_BE_NICE,
731 NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &oenv))
736 fplog = ftp2FILE(efLOG, NFILE, fnm, "w");
743 read_tpxheader(ftp2fn(efTPX, NFILE, fnm), &header, FALSE, NULL, NULL);
744 snew(xtop, header.natoms);
745 read_tpx(ftp2fn(efTPX, NFILE, fnm), &ir, box, &ntopatoms, xtop, NULL, NULL, &mtop);
746 bPDB = opt2bSet("-q", NFILE, fnm);
749 snew(xav, ntopatoms);
750 snew(ind_fit, ntopatoms);
751 snew(w_rls, ntopatoms);
752 for (kkk = 0; (kkk < ntopatoms); kkk++)
759 *atoms = gmx_mtop_global_atoms(&mtop);
761 if (atoms->pdbinfo == NULL)
763 snew(atoms->pdbinfo, atoms->nr);
767 top = gmx_mtop_generate_local_top(&mtop, &ir);
771 if (ir.ePBC != epbcNONE)
773 if (ir.bPeriodicMols)
779 g = mk_graph(fplog, &top->idef, 0, mtop.natoms, FALSE, FALSE);
783 if (ftp2bSet(efNDX, NFILE, fnm))
785 rd_index(ftp2fn(efNDX, NFILE, fnm), 1, &isize, &index, &grpname);
786 xvg = xvgropen(opt2fn("-dr", NFILE, fnm), "Individual Restraints", "Time (ps)",
788 snew(vvindex, isize);
790 for (i = 0; (i < isize); i++)
794 sprintf(leg[i], "index %d", index[i]);
796 xvgr_legend(xvg, isize, (const char**)leg, oenv);
804 init_disres(fplog, &mtop, &ir, NULL, &fcd, NULL, FALSE);
806 natoms = read_first_x(oenv, &status, ftp2fn(efTRX, NFILE, fnm), &t, &x, box);
809 init_dr_res(&dr, fcd.disres.nres);
810 if (opt2bSet("-c", NFILE, fnm))
812 clust = cluster_index(fplog, opt2fn("-c", NFILE, fnm));
813 snew(dr_clust, clust->clust->nr+1);
814 for (i = 0; (i <= clust->clust->nr); i++)
816 init_dr_res(&dr_clust[i], fcd.disres.nres);
821 out = xvgropen(opt2fn("-ds", NFILE, fnm),
822 "Sum of Violations", "Time (ps)", "nm", oenv);
823 aver = xvgropen(opt2fn("-da", NFILE, fnm),
824 "Average Violation", "Time (ps)", "nm", oenv);
825 numv = xvgropen(opt2fn("-dn", NFILE, fnm),
826 "# Violations", "Time (ps)", "#", oenv);
827 maxxv = xvgropen(opt2fn("-dm", NFILE, fnm),
828 "Largest Violation", "Time (ps)", "nm", oenv);
831 mdatoms = init_mdatoms(fplog, &mtop, ir.efep != efepNO);
832 atoms2md(&mtop, &ir, 0, NULL, mtop.natoms, mdatoms);
833 update_mdatoms(mdatoms, ir.fepvals->init_lambda);
835 if (ir.ePBC != epbcNONE)
837 gpbc = gmx_rmpbc_init(&top->idef, ir.ePBC, natoms);
843 if (ir.ePBC != epbcNONE)
845 if (ir.bPeriodicMols)
847 set_pbc(&pbc, ir.ePBC, box);
851 gmx_rmpbc(gpbc, natoms, box, x);
857 if (j > clust->maxframe)
859 gmx_fatal(FARGS, "There are more frames in the trajectory than in the cluster index file. t = %8f\n", t);
861 my_clust = clust->inv_clust[j];
862 range_check(my_clust, 0, clust->clust->nr);
863 check_viol(fplog, &(top->idef.il[F_DISRES]),
865 x, f, pbc_null, g, dr_clust, my_clust, isize, index, vvindex, &fcd);
869 check_viol(fplog, &(top->idef.il[F_DISRES]),
871 x, f, pbc_null, g, &dr, 0, isize, index, vvindex, &fcd);
875 reset_x(atoms->nr, ind_fit, atoms->nr, NULL, x, w_rls);
876 do_fit(atoms->nr, w_rls, x, x);
879 /* Store the first frame of the trajectory as 'characteristic'
880 * for colouring with violations.
882 for (kkk = 0; (kkk < atoms->nr); kkk++)
884 copy_rvec(x[kkk], xav[kkk]);
892 fprintf(xvg, "%10g", t);
893 for (i = 0; (i < isize); i++)
895 fprintf(xvg, " %10g", vvindex[i]);
899 fprintf(out, "%10g %10g\n", t, dr.sumv);
900 fprintf(aver, "%10g %10g\n", t, dr.averv);
901 fprintf(maxxv, "%10g %10g\n", t, dr.maxv);
902 fprintf(numv, "%10g %10d\n", t, dr.nv);
906 while (read_next_x(oenv, status, &t, x, box));
908 if (ir.ePBC != epbcNONE)
910 gmx_rmpbc_done(gpbc);
915 dump_clust_stats(fplog, fcd.disres.nres, &(top->idef.il[F_DISRES]),
916 top->idef.iparams, clust->clust, dr_clust,
917 clust->grpname, isize, index);
921 dump_stats(fplog, j, fcd.disres.nres, &(top->idef.il[F_DISRES]),
922 top->idef.iparams, &dr, isize, index,
923 bPDB ? atoms : NULL);
926 write_sto_conf(opt2fn("-q", NFILE, fnm),
927 "Coloured by average violation in Angstrom",
928 atoms, xav, NULL, ir.ePBC, box);
930 dump_disre_matrix(opt2fn_null("-x", NFILE, fnm), &dr, fcd.disres.nres,
931 j, &top->idef, &mtop, max_dr, nlevels, bThird);
939 do_view(oenv, opt2fn("-dr", NFILE, fnm), "-nxy");
941 do_view(oenv, opt2fn("-dn", NFILE, fnm), "-nxy");
942 do_view(oenv, opt2fn("-da", NFILE, fnm), "-nxy");
943 do_view(oenv, opt2fn("-ds", NFILE, fnm), "-nxy");
944 do_view(oenv, opt2fn("-dm", NFILE, fnm), "-nxy");
947 gmx_log_close(fplog);