3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Green Red Orange Magenta Azure Cyan Skyblue
63 #include "mtop_util.h"
77 real sumv, averv, maxv;
78 real *aver1, *aver2, *aver_3, *aver_6;
81 static void init5(int n)
87 static void reset5(void)
91 for (i = 0; (i < ntop); i++)
98 int tpcomp(const void *a, const void *b)
106 return (1e7*(tpb->v-tpa->v));
109 static void add5(int ndr, real viol)
114 for (i = 1; (i < ntop); i++)
116 if (top[i].v < top[mini].v)
121 if (viol > top[mini].v)
128 static void print5(FILE *fp)
132 qsort(top, ntop, sizeof(top[0]), tpcomp);
133 fprintf(fp, "Index:");
134 for (i = 0; (i < ntop); i++)
136 fprintf(fp, " %6d", top[i].n);
138 fprintf(fp, "\nViol: ");
139 for (i = 0; (i < ntop); i++)
141 fprintf(fp, " %6.3f", top[i].v);
146 static void check_viol(FILE *log, t_commrec *cr,
147 t_ilist *disres, t_iparams forceparams[],
148 t_functype functype[], rvec x[], rvec f[],
149 t_forcerec *fr, t_pbc *pbc, t_graph *g, t_dr_result dr[],
150 int clust_id, int isize, atom_id index[], real vvindex[],
154 int i, j, nat, n, type, nviol, ndr, label;
155 real ener, rt, mviol, tviol, viol, lam, dvdl, drt;
156 static gmx_bool bFirst = TRUE;
168 forceatoms = disres->iatoms;
169 for (j = 0; (j < isize); j++)
173 nat = interaction_function[F_DISRES].nratoms+1;
174 for (i = 0; (i < disres->nr); )
176 type = forceatoms[i];
178 label = forceparams[type].disres.label;
181 fprintf(debug, "DISRE: ndr = %d, label = %d i=%d, n =%d\n",
186 gmx_fatal(FARGS, "tpr inconsistency. ndr = %d, label = %d\n", ndr, label);
192 while (((i+n) < disres->nr) &&
193 (forceparams[forceatoms[i+n]].disres.label == label));
195 calc_disres_R_6(cr->ms, n, &forceatoms[i], forceparams,
196 (const rvec*)x, pbc, fcd, NULL);
198 if (fcd->disres.Rt_6[0] <= 0)
200 gmx_fatal(FARGS, "ndr = %d, rt_6 = %f", ndr, fcd->disres.Rt_6[0]);
203 rt = pow(fcd->disres.Rt_6[0], -1.0/6.0);
204 dr[clust_id].aver1[ndr] += rt;
205 dr[clust_id].aver2[ndr] += sqr(rt);
207 dr[clust_id].aver_3[ndr] += drt;
208 dr[clust_id].aver_6[ndr] += fcd->disres.Rt_6[0];
210 ener = interaction_function[F_DISRES].ifunc(n, &forceatoms[i], forceparams,
211 (const rvec*)x, f, fr->fshift,
212 pbc, g, lam, &dvdl, NULL, fcd, NULL);
213 viol = fcd->disres.sumviol;
220 add5(forceparams[type].disres.label, viol);
227 for (j = 0; (j < isize); j++)
229 if (index[j] == forceparams[type].disres.label)
231 vvindex[j] = pow(fcd->disres.Rt_6[0], -1.0/6.0);
238 dr[clust_id].nv = nviol;
239 dr[clust_id].maxv = mviol;
240 dr[clust_id].sumv = tviol;
241 dr[clust_id].averv = tviol/ndr;
242 dr[clust_id].nframes++;
246 fprintf(stderr, "\nThere are %d restraints and %d pairs\n",
247 ndr, disres->nr/nat);
259 real up1, r, rT3, rT6, viol, violT3, violT6;
262 static int drs_comp(const void *a, const void *b)
266 da = (t_dr_stats *)a;
267 db = (t_dr_stats *)b;
269 if (da->viol > db->viol)
273 else if (da->viol < db->viol)
283 static void dump_dump(FILE *log, int ndr, t_dr_stats drs[])
285 static const char *core[] = { "All restraints", "Core restraints" };
286 static const char *tp[] = { "linear", "third power", "sixth power" };
287 real viol_tot, viol_max, viol = 0;
292 for (bCore = FALSE; (bCore <= TRUE); bCore++)
294 for (kkk = 0; (kkk < 3); kkk++)
300 for (i = 0; (i < ndr); i++)
302 if (!bCore || (bCore && drs[i].bCore))
310 viol = drs[i].violT3;
313 viol = drs[i].violT6;
316 gmx_incons("Dumping violations");
318 viol_max = max(viol_max, viol);
327 if ((nrestr > 0) || (bCore && (nrestr < ndr)))
330 fprintf(log, "+++++++ %s ++++++++\n", core[bCore]);
331 fprintf(log, "+++++++ Using %s averaging: ++++++++\n", tp[kkk]);
332 fprintf(log, "Sum of violations: %8.3f nm\n", viol_tot);
335 fprintf(log, "Average violation: %8.3f nm\n", viol_tot/nrestr);
337 fprintf(log, "Largest violation: %8.3f nm\n", viol_max);
338 fprintf(log, "Number of violated restraints: %d/%d\n", nviol, nrestr);
344 static void dump_viol(FILE *log, int ndr, t_dr_stats *drs, gmx_bool bLinear)
348 fprintf(log, "Restr. Core Up1 <r> <rT3> <rT6> <viol><violT3><violT6>\n");
349 for (i = 0; (i < ndr); i++)
351 if (bLinear && (drs[i].viol == 0))
355 fprintf(log, "%6d%5s%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f\n",
356 drs[i].index, yesno_names[drs[i].bCore],
357 drs[i].up1, drs[i].r, drs[i].rT3, drs[i].rT6,
358 drs[i].viol, drs[i].violT3, drs[i].violT6);
362 static gmx_bool is_core(int i, int isize, atom_id index[])
365 gmx_bool bIC = FALSE;
367 for (kk = 0; !bIC && (kk < isize); kk++)
369 bIC = (index[kk] == i);
375 static void dump_stats(FILE *log, int nsteps, int ndr, t_ilist *disres,
376 t_iparams ip[], t_dr_result *dr,
377 int isize, atom_id index[], t_atoms *atoms)
383 fprintf(log, "++++++++++++++ STATISTICS ++++++++++++++++++++++++\n");
386 nra = interaction_function[F_DISRES].nratoms+1;
387 for (i = 0; (i < ndr); i++)
389 /* Search for the appropriate restraint */
390 for (; (j < disres->nr); j += nra)
392 if (ip[disres->iatoms[j]].disres.label == i)
398 drs[i].bCore = is_core(i, isize, index);
399 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
400 drs[i].r = dr->aver1[i]/nsteps;
401 drs[i].rT3 = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
402 drs[i].rT6 = pow(dr->aver_6[i]/nsteps, -1.0/6.0);
403 drs[i].viol = max(0, drs[i].r-drs[i].up1);
404 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
405 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
408 int j1 = disres->iatoms[j+1];
409 int j2 = disres->iatoms[j+2];
410 atoms->pdbinfo[j1].bfac += drs[i].violT3*5;
411 atoms->pdbinfo[j2].bfac += drs[i].violT3*5;
414 dump_viol(log, ndr, drs, FALSE);
416 fprintf(log, "+++ Sorted by linear averaged violations: +++\n");
417 qsort(drs, ndr, sizeof(drs[0]), drs_comp);
418 dump_viol(log, ndr, drs, TRUE);
420 dump_dump(log, ndr, drs);
425 static void dump_clust_stats(FILE *fp, int ndr, t_ilist *disres,
426 t_iparams ip[], t_blocka *clust, t_dr_result dr[],
427 char *clust_name[], int isize, atom_id index[])
429 int i, j, k, nra, mmm = 0;
430 double sumV, maxV, sumVT3, sumVT6, maxVT3, maxVT6;
434 fprintf(fp, "++++++++++++++ STATISTICS ++++++++++++++++++++++\n");
435 fprintf(fp, "Cluster NFrames SumV MaxV SumVT MaxVT SumVS MaxVS\n");
439 for (k = 0; (k < clust->nr); k++)
441 if (dr[k].nframes == 0)
445 if (dr[k].nframes != (clust->index[k+1]-clust->index[k]))
447 gmx_fatal(FARGS, "Inconsistency in cluster %s.\n"
448 "Found %d frames in trajectory rather than the expected %d\n",
449 clust_name[k], dr[k].nframes,
450 clust->index[k+1]-clust->index[k]);
454 gmx_fatal(FARGS, "Inconsistency with cluster %d. Invalid name", k);
457 nra = interaction_function[F_DISRES].nratoms+1;
458 sumV = sumVT3 = sumVT6 = maxV = maxVT3 = maxVT6 = 0;
459 for (i = 0; (i < ndr); i++)
461 /* Search for the appropriate restraint */
462 for (; (j < disres->nr); j += nra)
464 if (ip[disres->iatoms[j]].disres.label == i)
470 drs[i].bCore = is_core(i, isize, index);
471 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
472 drs[i].r = dr[k].aver1[i]/dr[k].nframes;
473 if ((dr[k].aver_3[i] <= 0) || (dr[k].aver_3[i] != dr[k].aver_3[i]))
475 gmx_fatal(FARGS, "dr[%d].aver_3[%d] = %f", k, i, dr[k].aver_3[i]);
477 drs[i].rT3 = pow(dr[k].aver_3[i]/dr[k].nframes, -1.0/3.0);
478 drs[i].rT6 = pow(dr[k].aver_6[i]/dr[k].nframes, -1.0/6.0);
479 drs[i].viol = max(0, drs[i].r-drs[i].up1);
480 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
481 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
483 sumVT3 += drs[i].violT3;
484 sumVT6 += drs[i].violT6;
485 maxV = max(maxV, drs[i].viol);
486 maxVT3 = max(maxVT3, drs[i].violT3);
487 maxVT6 = max(maxVT6, drs[i].violT6);
489 if (strcmp(clust_name[k], "1000") == 0)
493 fprintf(fp, "%-10s%6d%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n",
495 dr[k].nframes, sumV, maxV, sumVT3, maxVT3, sumVT6, maxVT6);
502 static void init_dr_res(t_dr_result *dr, int ndr)
504 snew(dr->aver1, ndr+1);
505 snew(dr->aver2, ndr+1);
506 snew(dr->aver_3, ndr+1);
507 snew(dr->aver_6, ndr+1);
515 static void dump_disre_matrix(const char *fn, t_dr_result *dr, int ndr,
516 int nsteps, t_idef *idef, gmx_mtop_t *mtop,
517 real max_dr, int nlevels, gmx_bool bThird)
521 int n_res, a_offset, mb, mol, a;
523 int iii, i, j, nra, nratoms, tp, ri, rj, index, nlabel, label;
524 atom_id ai, aj, *ptr;
525 real **matrix, *t_res, hi, *w_dr, rav, rviol;
526 t_rgb rlo = { 1, 1, 1 };
527 t_rgb rhi = { 0, 0, 0 };
533 snew(resnr, mtop->natoms);
536 for (mb = 0; mb < mtop->nmolblock; mb++)
538 atoms = &mtop->moltype[mtop->molblock[mb].type].atoms;
539 for (mol = 0; mol < mtop->molblock[mb].nmol; mol++)
541 for (a = 0; a < atoms->nr; a++)
543 resnr[a_offset+a] = n_res + atoms->atom[a].resind;
545 n_res += atoms->nres;
546 a_offset += atoms->nr;
551 for (i = 0; (i < n_res); i++)
556 for (i = 0; (i < n_res); i++)
558 snew(matrix[i], n_res);
560 nratoms = interaction_function[F_DISRES].nratoms;
561 nra = (idef->il[F_DISRES].nr/(nratoms+1));
567 for (i = 0; (i < idef->il[F_DISRES].nr); i += nratoms+1)
569 tp = idef->il[F_DISRES].iatoms[i];
570 label = idef->iparams[tp].disres.label;
574 /* Set index pointer */
578 gmx_fatal(FARGS, "nlabel is %d, label = %d", nlabel, label);
582 gmx_fatal(FARGS, "ndr = %d, index = %d", ndr, index);
584 /* Update the weight */
585 w_dr[index] = 1.0/nlabel;
594 printf("nlabel = %d, index = %d, ndr = %d\n", nlabel, index, ndr);
596 for (i = 0; (i < ndr); i++)
598 for (j = ptr[i]; (j < ptr[i+1]); j += nratoms+1)
600 tp = idef->il[F_DISRES].iatoms[j];
601 ai = idef->il[F_DISRES].iatoms[j+1];
602 aj = idef->il[F_DISRES].iatoms[j+2];
608 rav = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
612 rav = dr->aver1[i]/nsteps;
616 fprintf(debug, "DR %d, atoms %d, %d, distance %g\n", i, ai, aj, rav);
618 rviol = max(0, rav-idef->iparams[tp].disres.up1);
619 matrix[ri][rj] += w_dr[i]*rviol;
620 matrix[rj][ri] += w_dr[i]*rviol;
621 hi = max(hi, matrix[ri][rj]);
622 hi = max(hi, matrix[rj][ri]);
632 printf("Warning: the maxdr that you have specified (%g) is smaller than\nthe largest value in your simulation (%g)\n", max_dr, hi);
636 printf("Highest level in the matrix will be %g\n", hi);
637 fp = ffopen(fn, "w");
638 write_xpm(fp, 0, "Distance Violations", "<V> (nm)", "Residue", "Residue",
639 n_res, n_res, t_res, t_res, matrix, 0, hi, rlo, rhi, &nlevels);
643 int gmx_disre(int argc, char *argv[])
645 const char *desc[] = {
646 "[TT]g_disre[tt] computes violations of distance restraints.",
647 "If necessary, all protons can be added to a protein molecule ",
648 "using the [TT]g_protonate[tt] program.[PAR]",
649 "The program always",
650 "computes the instantaneous violations rather than time-averaged,",
651 "because this analysis is done from a trajectory file afterwards",
652 "it does not make sense to use time averaging. However,",
653 "the time averaged values per restraint are given in the log file.[PAR]",
654 "An index file may be used to select specific restraints for",
656 "When the optional [TT]-q[tt] flag is given a [TT].pdb[tt] file coloured by the",
657 "amount of average violations.[PAR]",
658 "When the [TT]-c[tt] option is given, an index file will be read",
659 "containing the frames in your trajectory corresponding to the clusters",
660 "(defined in another manner) that you want to analyze. For these clusters",
661 "the program will compute average violations using the third power",
662 "averaging algorithm and print them in the log file."
665 static int nlevels = 20;
666 static real max_dr = 0;
667 static gmx_bool bThird = TRUE;
669 { "-ntop", FALSE, etINT, {&ntop},
670 "Number of large violations that are stored in the log file every step" },
671 { "-maxdr", FALSE, etREAL, {&max_dr},
672 "Maximum distance violation in matrix output. If less than or equal to 0 the maximum will be determined by the data." },
673 { "-nlevels", FALSE, etINT, {&nlevels},
674 "Number of levels in the matrix output" },
675 { "-third", FALSE, etBOOL, {&bThird},
676 "Use inverse third power averaging or linear for matrix output" }
679 FILE *out = NULL, *aver = NULL, *numv = NULL, *maxxv = NULL, *xvg = NULL;
685 t_atoms *atoms = NULL;
691 int ntopatoms, natoms, i, j, kkk;
694 rvec *x, *f, *xav = NULL;
698 atom_id *index = NULL, *ind_fit = NULL;
700 t_cluster_ndx *clust = NULL;
701 t_dr_result dr, *dr_clust = NULL;
703 real *vvindex = NULL, *w_rls = NULL;
705 t_pbc pbc, *pbc_null;
709 gmx_rmpbc_t gpbc = NULL;
712 { efTPX, NULL, NULL, ffREAD },
713 { efTRX, "-f", NULL, ffREAD },
714 { efXVG, "-ds", "drsum", ffWRITE },
715 { efXVG, "-da", "draver", ffWRITE },
716 { efXVG, "-dn", "drnum", ffWRITE },
717 { efXVG, "-dm", "drmax", ffWRITE },
718 { efXVG, "-dr", "restr", ffWRITE },
719 { efLOG, "-l", "disres", ffWRITE },
720 { efNDX, NULL, "viol", ffOPTRD },
721 { efPDB, "-q", "viol", ffOPTWR },
722 { efNDX, "-c", "clust", ffOPTRD },
723 { efXPM, "-x", "matrix", ffOPTWR }
725 #define NFILE asize(fnm)
727 cr = init_par(&argc, &argv);
728 parse_common_args(&argc, argv, PCA_CAN_TIME | PCA_CAN_VIEW | PCA_BE_NICE,
729 NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &oenv);
731 gmx_log_open(ftp2fn(efLOG, NFILE, fnm), cr, FALSE, 0, &fplog);
738 read_tpxheader(ftp2fn(efTPX, NFILE, fnm), &header, FALSE, NULL, NULL);
739 snew(xtop, header.natoms);
740 read_tpx(ftp2fn(efTPX, NFILE, fnm), &ir, box, &ntopatoms, xtop, NULL, NULL, &mtop);
741 bPDB = opt2bSet("-q", NFILE, fnm);
744 snew(xav, ntopatoms);
745 snew(ind_fit, ntopatoms);
746 snew(w_rls, ntopatoms);
747 for (kkk = 0; (kkk < ntopatoms); kkk++)
754 *atoms = gmx_mtop_global_atoms(&mtop);
756 if (atoms->pdbinfo == NULL)
758 snew(atoms->pdbinfo, atoms->nr);
762 top = gmx_mtop_generate_local_top(&mtop, &ir);
766 if (ir.ePBC != epbcNONE)
768 if (ir.bPeriodicMols)
774 g = mk_graph(fplog, &top->idef, 0, mtop.natoms, FALSE, FALSE);
778 if (ftp2bSet(efNDX, NFILE, fnm))
780 rd_index(ftp2fn(efNDX, NFILE, fnm), 1, &isize, &index, &grpname);
781 xvg = xvgropen(opt2fn("-dr", NFILE, fnm), "Individual Restraints", "Time (ps)",
783 snew(vvindex, isize);
785 for (i = 0; (i < isize); i++)
789 sprintf(leg[i], "index %d", index[i]);
791 xvgr_legend(xvg, isize, (const char**)leg, oenv);
799 init_disres(fplog, &mtop, &ir, NULL, FALSE, &fcd, NULL, FALSE);
801 natoms = read_first_x(oenv, &status, ftp2fn(efTRX, NFILE, fnm), &t, &x, box);
804 init_dr_res(&dr, fcd.disres.nres);
805 if (opt2bSet("-c", NFILE, fnm))
807 clust = cluster_index(fplog, opt2fn("-c", NFILE, fnm));
808 snew(dr_clust, clust->clust->nr+1);
809 for (i = 0; (i <= clust->clust->nr); i++)
811 init_dr_res(&dr_clust[i], fcd.disres.nres);
816 out = xvgropen(opt2fn("-ds", NFILE, fnm),
817 "Sum of Violations", "Time (ps)", "nm", oenv);
818 aver = xvgropen(opt2fn("-da", NFILE, fnm),
819 "Average Violation", "Time (ps)", "nm", oenv);
820 numv = xvgropen(opt2fn("-dn", NFILE, fnm),
821 "# Violations", "Time (ps)", "#", oenv);
822 maxxv = xvgropen(opt2fn("-dm", NFILE, fnm),
823 "Largest Violation", "Time (ps)", "nm", oenv);
826 mdatoms = init_mdatoms(fplog, &mtop, ir.efep != efepNO);
827 atoms2md(&mtop, &ir, 0, NULL, 0, mtop.natoms, mdatoms);
828 update_mdatoms(mdatoms, ir.fepvals->init_lambda);
830 fprintf(fplog, "Made forcerec\n");
831 init_forcerec(fplog, oenv, fr, NULL, &ir, &mtop, cr, box, FALSE,
832 NULL, NULL, NULL, NULL, NULL, FALSE, -1);
834 if (ir.ePBC != epbcNONE)
836 gpbc = gmx_rmpbc_init(&top->idef, ir.ePBC, natoms, box);
842 if (ir.ePBC != epbcNONE)
844 if (ir.bPeriodicMols)
846 set_pbc(&pbc, ir.ePBC, box);
850 gmx_rmpbc(gpbc, natoms, box, x);
856 if (j > clust->maxframe)
858 gmx_fatal(FARGS, "There are more frames in the trajectory than in the cluster index file. t = %8f\n", t);
860 my_clust = clust->inv_clust[j];
861 range_check(my_clust, 0, clust->clust->nr);
862 check_viol(fplog, cr, &(top->idef.il[F_DISRES]),
863 top->idef.iparams, top->idef.functype,
864 x, f, fr, pbc_null, g, dr_clust, my_clust, isize, index, vvindex, &fcd);
868 check_viol(fplog, cr, &(top->idef.il[F_DISRES]),
869 top->idef.iparams, top->idef.functype,
870 x, f, fr, pbc_null, g, &dr, 0, isize, index, vvindex, &fcd);
874 reset_x(atoms->nr, ind_fit, atoms->nr, NULL, x, w_rls);
875 do_fit(atoms->nr, w_rls, x, x);
878 /* Store the first frame of the trajectory as 'characteristic'
879 * for colouring with violations.
881 for (kkk = 0; (kkk < atoms->nr); kkk++)
883 copy_rvec(x[kkk], xav[kkk]);
891 fprintf(xvg, "%10g", t);
892 for (i = 0; (i < isize); i++)
894 fprintf(xvg, " %10g", vvindex[i]);
898 fprintf(out, "%10g %10g\n", t, dr.sumv);
899 fprintf(aver, "%10g %10g\n", t, dr.averv);
900 fprintf(maxxv, "%10g %10g\n", t, dr.maxv);
901 fprintf(numv, "%10g %10d\n", t, dr.nv);
905 while (read_next_x(oenv, status, &t, natoms, x, box));
907 if (ir.ePBC != epbcNONE)
909 gmx_rmpbc_done(gpbc);
914 dump_clust_stats(fplog, fcd.disres.nres, &(top->idef.il[F_DISRES]),
915 top->idef.iparams, clust->clust, dr_clust,
916 clust->grpname, isize, index);
920 dump_stats(fplog, j, fcd.disres.nres, &(top->idef.il[F_DISRES]),
921 top->idef.iparams, &dr, isize, index,
922 bPDB ? atoms : NULL);
925 write_sto_conf(opt2fn("-q", NFILE, fnm),
926 "Coloured by average violation in Angstrom",
927 atoms, xav, NULL, ir.ePBC, box);
929 dump_disre_matrix(opt2fn_null("-x", NFILE, fnm), &dr, fcd.disres.nres,
930 j, &top->idef, &mtop, max_dr, nlevels, bThird);
938 do_view(oenv, opt2fn("-dr", NFILE, fnm), "-nxy");
940 do_view(oenv, opt2fn("-dn", NFILE, fnm), "-nxy");
941 do_view(oenv, opt2fn("-da", NFILE, fnm), "-nxy");
942 do_view(oenv, opt2fn("-ds", NFILE, fnm), "-nxy");
943 do_view(oenv, opt2fn("-dm", NFILE, fnm), "-nxy");
949 gmx_log_close(fplog);