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.
48 #include "gromacs/fileio/confio.h"
49 #include "gromacs/utility/smalloc.h"
52 #include "gromacs/commandline/pargs.h"
56 #include "gromacs/fileio/pdbio.h"
59 #include "gromacs/fileio/tpxio.h"
60 #include "gromacs/fileio/trxio.h"
63 #include "gromacs/fileio/matio.h"
64 #include "mtop_util.h"
67 #include "gromacs/math/do_fit.h"
68 #include "gromacs/legacyheaders/gmx_fatal.h"
80 real sumv, averv, maxv;
81 real *aver1, *aver2, *aver_3, *aver_6;
84 static void init5(int n)
90 static void reset5(void)
94 for (i = 0; (i < ntop); i++)
101 int tpcomp(const void *a, const void *b)
109 return (1e7*(tpb->v-tpa->v));
112 static void add5(int ndr, real viol)
117 for (i = 1; (i < ntop); i++)
119 if (top[i].v < top[mini].v)
124 if (viol > top[mini].v)
131 static void print5(FILE *fp)
135 qsort(top, ntop, sizeof(top[0]), tpcomp);
136 fprintf(fp, "Index:");
137 for (i = 0; (i < ntop); i++)
139 fprintf(fp, " %6d", top[i].n);
141 fprintf(fp, "\nViol: ");
142 for (i = 0; (i < ntop); i++)
144 fprintf(fp, " %6.3f", top[i].v);
149 static void check_viol(FILE *log,
150 t_ilist *disres, t_iparams forceparams[],
152 t_pbc *pbc, t_graph *g, t_dr_result dr[],
153 int clust_id, int isize, atom_id index[], real vvindex[],
157 int i, j, nat, n, type, nviol, ndr, label;
158 real ener, rt, mviol, tviol, viol, lam, dvdl, drt;
160 static gmx_bool bFirst = TRUE;
172 forceatoms = disres->iatoms;
173 for (j = 0; (j < isize); j++)
177 nat = interaction_function[F_DISRES].nratoms+1;
178 for (i = 0; (i < disres->nr); )
180 type = forceatoms[i];
182 label = forceparams[type].disres.label;
185 fprintf(debug, "DISRE: ndr = %d, label = %d i=%d, n =%d\n",
190 gmx_fatal(FARGS, "tpr inconsistency. ndr = %d, label = %d\n", ndr, label);
196 while (((i+n) < disres->nr) &&
197 (forceparams[forceatoms[i+n]].disres.label == label));
199 calc_disres_R_6(n, &forceatoms[i], forceparams,
200 (const rvec*)x, pbc, fcd, NULL);
202 if (fcd->disres.Rt_6[0] <= 0)
204 gmx_fatal(FARGS, "ndr = %d, rt_6 = %f", ndr, fcd->disres.Rt_6[0]);
207 rt = pow(fcd->disres.Rt_6[0], -1.0/6.0);
208 dr[clust_id].aver1[ndr] += rt;
209 dr[clust_id].aver2[ndr] += sqr(rt);
211 dr[clust_id].aver_3[ndr] += drt;
212 dr[clust_id].aver_6[ndr] += fcd->disres.Rt_6[0];
214 snew(fshift, SHIFTS);
215 ener = interaction_function[F_DISRES].ifunc(n, &forceatoms[i], forceparams,
216 (const rvec*)x, f, fshift,
217 pbc, g, lam, &dvdl, NULL, fcd, NULL);
219 viol = fcd->disres.sumviol;
226 add5(forceparams[type].disres.label, viol);
233 for (j = 0; (j < isize); j++)
235 if (index[j] == forceparams[type].disres.label)
237 vvindex[j] = pow(fcd->disres.Rt_6[0], -1.0/6.0);
244 dr[clust_id].nv = nviol;
245 dr[clust_id].maxv = mviol;
246 dr[clust_id].sumv = tviol;
247 dr[clust_id].averv = tviol/ndr;
248 dr[clust_id].nframes++;
252 fprintf(stderr, "\nThere are %d restraints and %d pairs\n",
253 ndr, disres->nr/nat);
265 real up1, r, rT3, rT6, viol, violT3, violT6;
268 static int drs_comp(const void *a, const void *b)
272 da = (t_dr_stats *)a;
273 db = (t_dr_stats *)b;
275 if (da->viol > db->viol)
279 else if (da->viol < db->viol)
289 static void dump_dump(FILE *log, int ndr, t_dr_stats drs[])
291 static const char *core[] = { "All restraints", "Core restraints" };
292 static const char *tp[] = { "linear", "third power", "sixth power" };
293 real viol_tot, viol_max, viol = 0;
298 for (bCore = FALSE; (bCore <= TRUE); bCore++)
300 for (kkk = 0; (kkk < 3); kkk++)
306 for (i = 0; (i < ndr); i++)
308 if (!bCore || (bCore && drs[i].bCore))
316 viol = drs[i].violT3;
319 viol = drs[i].violT6;
322 gmx_incons("Dumping violations");
324 viol_max = max(viol_max, viol);
333 if ((nrestr > 0) || (bCore && (nrestr < ndr)))
336 fprintf(log, "+++++++ %s ++++++++\n", core[bCore]);
337 fprintf(log, "+++++++ Using %s averaging: ++++++++\n", tp[kkk]);
338 fprintf(log, "Sum of violations: %8.3f nm\n", viol_tot);
341 fprintf(log, "Average violation: %8.3f nm\n", viol_tot/nrestr);
343 fprintf(log, "Largest violation: %8.3f nm\n", viol_max);
344 fprintf(log, "Number of violated restraints: %d/%d\n", nviol, nrestr);
350 static void dump_viol(FILE *log, int ndr, t_dr_stats *drs, gmx_bool bLinear)
354 fprintf(log, "Restr. Core Up1 <r> <rT3> <rT6> <viol><violT3><violT6>\n");
355 for (i = 0; (i < ndr); i++)
357 if (bLinear && (drs[i].viol == 0))
361 fprintf(log, "%6d%5s%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f\n",
362 drs[i].index, yesno_names[drs[i].bCore],
363 drs[i].up1, drs[i].r, drs[i].rT3, drs[i].rT6,
364 drs[i].viol, drs[i].violT3, drs[i].violT6);
368 static gmx_bool is_core(int i, int isize, atom_id index[])
371 gmx_bool bIC = FALSE;
373 for (kk = 0; !bIC && (kk < isize); kk++)
375 bIC = (index[kk] == i);
381 static void dump_stats(FILE *log, int nsteps, int ndr, t_ilist *disres,
382 t_iparams ip[], t_dr_result *dr,
383 int isize, atom_id index[], t_atoms *atoms)
389 fprintf(log, "++++++++++++++ STATISTICS ++++++++++++++++++++++++\n");
392 nra = interaction_function[F_DISRES].nratoms+1;
393 for (i = 0; (i < ndr); i++)
395 /* Search for the appropriate restraint */
396 for (; (j < disres->nr); j += nra)
398 if (ip[disres->iatoms[j]].disres.label == i)
404 drs[i].bCore = is_core(i, isize, index);
405 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
406 drs[i].r = dr->aver1[i]/nsteps;
407 drs[i].rT3 = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
408 drs[i].rT6 = pow(dr->aver_6[i]/nsteps, -1.0/6.0);
409 drs[i].viol = max(0, drs[i].r-drs[i].up1);
410 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
411 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
414 int j1 = disres->iatoms[j+1];
415 int j2 = disres->iatoms[j+2];
416 atoms->pdbinfo[j1].bfac += drs[i].violT3*5;
417 atoms->pdbinfo[j2].bfac += drs[i].violT3*5;
420 dump_viol(log, ndr, drs, FALSE);
422 fprintf(log, "+++ Sorted by linear averaged violations: +++\n");
423 qsort(drs, ndr, sizeof(drs[0]), drs_comp);
424 dump_viol(log, ndr, drs, TRUE);
426 dump_dump(log, ndr, drs);
431 static void dump_clust_stats(FILE *fp, int ndr, t_ilist *disres,
432 t_iparams ip[], t_blocka *clust, t_dr_result dr[],
433 char *clust_name[], int isize, atom_id index[])
435 int i, j, k, nra, mmm = 0;
436 double sumV, maxV, sumVT3, sumVT6, maxVT3, maxVT6;
440 fprintf(fp, "++++++++++++++ STATISTICS ++++++++++++++++++++++\n");
441 fprintf(fp, "Cluster NFrames SumV MaxV SumVT MaxVT SumVS MaxVS\n");
445 for (k = 0; (k < clust->nr); k++)
447 if (dr[k].nframes == 0)
451 if (dr[k].nframes != (clust->index[k+1]-clust->index[k]))
453 gmx_fatal(FARGS, "Inconsistency in cluster %s.\n"
454 "Found %d frames in trajectory rather than the expected %d\n",
455 clust_name[k], dr[k].nframes,
456 clust->index[k+1]-clust->index[k]);
460 gmx_fatal(FARGS, "Inconsistency with cluster %d. Invalid name", k);
463 nra = interaction_function[F_DISRES].nratoms+1;
464 sumV = sumVT3 = sumVT6 = maxV = maxVT3 = maxVT6 = 0;
465 for (i = 0; (i < ndr); i++)
467 /* Search for the appropriate restraint */
468 for (; (j < disres->nr); j += nra)
470 if (ip[disres->iatoms[j]].disres.label == i)
476 drs[i].bCore = is_core(i, isize, index);
477 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
478 drs[i].r = dr[k].aver1[i]/dr[k].nframes;
479 if ((dr[k].aver_3[i] <= 0) || (dr[k].aver_3[i] != dr[k].aver_3[i]))
481 gmx_fatal(FARGS, "dr[%d].aver_3[%d] = %f", k, i, dr[k].aver_3[i]);
483 drs[i].rT3 = pow(dr[k].aver_3[i]/dr[k].nframes, -1.0/3.0);
484 drs[i].rT6 = pow(dr[k].aver_6[i]/dr[k].nframes, -1.0/6.0);
485 drs[i].viol = max(0, drs[i].r-drs[i].up1);
486 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
487 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
489 sumVT3 += drs[i].violT3;
490 sumVT6 += drs[i].violT6;
491 maxV = max(maxV, drs[i].viol);
492 maxVT3 = max(maxVT3, drs[i].violT3);
493 maxVT6 = max(maxVT6, drs[i].violT6);
495 if (strcmp(clust_name[k], "1000") == 0)
499 fprintf(fp, "%-10s%6d%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n",
501 dr[k].nframes, sumV, maxV, sumVT3, maxVT3, sumVT6, maxVT6);
508 static void init_dr_res(t_dr_result *dr, int ndr)
510 snew(dr->aver1, ndr+1);
511 snew(dr->aver2, ndr+1);
512 snew(dr->aver_3, ndr+1);
513 snew(dr->aver_6, ndr+1);
521 static void dump_disre_matrix(const char *fn, t_dr_result *dr, int ndr,
522 int nsteps, t_idef *idef, gmx_mtop_t *mtop,
523 real max_dr, int nlevels, gmx_bool bThird)
527 int n_res, a_offset, mb, mol, a;
529 int iii, i, j, nra, nratoms, tp, ri, rj, index, nlabel, label;
530 atom_id ai, aj, *ptr;
531 real **matrix, *t_res, hi, *w_dr, rav, rviol;
532 t_rgb rlo = { 1, 1, 1 };
533 t_rgb rhi = { 0, 0, 0 };
539 snew(resnr, mtop->natoms);
542 for (mb = 0; mb < mtop->nmolblock; mb++)
544 atoms = &mtop->moltype[mtop->molblock[mb].type].atoms;
545 for (mol = 0; mol < mtop->molblock[mb].nmol; mol++)
547 for (a = 0; a < atoms->nr; a++)
549 resnr[a_offset+a] = n_res + atoms->atom[a].resind;
551 n_res += atoms->nres;
552 a_offset += atoms->nr;
557 for (i = 0; (i < n_res); i++)
562 for (i = 0; (i < n_res); i++)
564 snew(matrix[i], n_res);
566 nratoms = interaction_function[F_DISRES].nratoms;
567 nra = (idef->il[F_DISRES].nr/(nratoms+1));
573 for (i = 0; (i < idef->il[F_DISRES].nr); i += nratoms+1)
575 tp = idef->il[F_DISRES].iatoms[i];
576 label = idef->iparams[tp].disres.label;
580 /* Set index pointer */
584 gmx_fatal(FARGS, "nlabel is %d, label = %d", nlabel, label);
588 gmx_fatal(FARGS, "ndr = %d, index = %d", ndr, index);
590 /* Update the weight */
591 w_dr[index] = 1.0/nlabel;
600 printf("nlabel = %d, index = %d, ndr = %d\n", nlabel, index, ndr);
602 for (i = 0; (i < ndr); i++)
604 for (j = ptr[i]; (j < ptr[i+1]); j += nratoms+1)
606 tp = idef->il[F_DISRES].iatoms[j];
607 ai = idef->il[F_DISRES].iatoms[j+1];
608 aj = idef->il[F_DISRES].iatoms[j+2];
614 rav = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
618 rav = dr->aver1[i]/nsteps;
622 fprintf(debug, "DR %d, atoms %d, %d, distance %g\n", i, ai, aj, rav);
624 rviol = max(0, rav-idef->iparams[tp].disres.up1);
625 matrix[ri][rj] += w_dr[i]*rviol;
626 matrix[rj][ri] += w_dr[i]*rviol;
627 hi = max(hi, matrix[ri][rj]);
628 hi = max(hi, matrix[rj][ri]);
638 printf("Warning: the maxdr that you have specified (%g) is smaller than\nthe largest value in your simulation (%g)\n", max_dr, hi);
642 printf("Highest level in the matrix will be %g\n", hi);
643 fp = gmx_ffopen(fn, "w");
644 write_xpm(fp, 0, "Distance Violations", "<V> (nm)", "Residue", "Residue",
645 n_res, n_res, t_res, t_res, matrix, 0, hi, rlo, rhi, &nlevels);
649 int gmx_disre(int argc, char *argv[])
651 const char *desc[] = {
652 "[THISMODULE] computes violations of distance restraints.",
653 "If necessary, all protons can be added to a protein molecule ",
654 "using the [gmx-protonate] program.[PAR]",
655 "The program always",
656 "computes the instantaneous violations rather than time-averaged,",
657 "because this analysis is done from a trajectory file afterwards",
658 "it does not make sense to use time averaging. However,",
659 "the time averaged values per restraint are given in the log file.[PAR]",
660 "An index file may be used to select specific restraints for",
662 "When the optional [TT]-q[tt] flag is given a [TT].pdb[tt] file coloured by the",
663 "amount of average violations.[PAR]",
664 "When the [TT]-c[tt] option is given, an index file will be read",
665 "containing the frames in your trajectory corresponding to the clusters",
666 "(defined in another manner) that you want to analyze. For these clusters",
667 "the program will compute average violations using the third power",
668 "averaging algorithm and print them in the log file."
671 static int nlevels = 20;
672 static real max_dr = 0;
673 static gmx_bool bThird = TRUE;
675 { "-ntop", FALSE, etINT, {&ntop},
676 "Number of large violations that are stored in the log file every step" },
677 { "-maxdr", FALSE, etREAL, {&max_dr},
678 "Maximum distance violation in matrix output. If less than or equal to 0 the maximum will be determined by the data." },
679 { "-nlevels", FALSE, etINT, {&nlevels},
680 "Number of levels in the matrix output" },
681 { "-third", FALSE, etBOOL, {&bThird},
682 "Use inverse third power averaging or linear for matrix output" }
685 FILE *out = NULL, *aver = NULL, *numv = NULL, *maxxv = NULL, *xvg = NULL;
691 t_atoms *atoms = NULL;
695 int ntopatoms, natoms, i, j, kkk;
698 rvec *x, *f, *xav = NULL;
702 atom_id *index = NULL, *ind_fit = NULL;
704 t_cluster_ndx *clust = NULL;
705 t_dr_result dr, *dr_clust = NULL;
707 real *vvindex = NULL, *w_rls = NULL;
709 t_pbc pbc, *pbc_null;
713 gmx_rmpbc_t gpbc = NULL;
716 { efTPX, NULL, NULL, ffREAD },
717 { efTRX, "-f", NULL, ffREAD },
718 { efXVG, "-ds", "drsum", ffWRITE },
719 { efXVG, "-da", "draver", ffWRITE },
720 { efXVG, "-dn", "drnum", ffWRITE },
721 { efXVG, "-dm", "drmax", ffWRITE },
722 { efXVG, "-dr", "restr", ffWRITE },
723 { efLOG, "-l", "disres", ffWRITE },
724 { efNDX, NULL, "viol", ffOPTRD },
725 { efPDB, "-q", "viol", ffOPTWR },
726 { efNDX, "-c", "clust", ffOPTRD },
727 { efXPM, "-x", "matrix", ffOPTWR }
729 #define NFILE asize(fnm)
731 if (!parse_common_args(&argc, argv, PCA_CAN_TIME | PCA_CAN_VIEW | PCA_BE_NICE,
732 NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &oenv))
737 fplog = ftp2FILE(efLOG, NFILE, fnm, "w");
744 read_tpxheader(ftp2fn(efTPX, NFILE, fnm), &header, FALSE, NULL, NULL);
745 snew(xtop, header.natoms);
746 read_tpx(ftp2fn(efTPX, NFILE, fnm), &ir, box, &ntopatoms, xtop, NULL, NULL, &mtop);
747 bPDB = opt2bSet("-q", NFILE, fnm);
750 snew(xav, ntopatoms);
751 snew(ind_fit, ntopatoms);
752 snew(w_rls, ntopatoms);
753 for (kkk = 0; (kkk < ntopatoms); kkk++)
760 *atoms = gmx_mtop_global_atoms(&mtop);
762 if (atoms->pdbinfo == NULL)
764 snew(atoms->pdbinfo, atoms->nr);
768 top = gmx_mtop_generate_local_top(&mtop, &ir);
772 if (ir.ePBC != epbcNONE)
774 if (ir.bPeriodicMols)
780 g = mk_graph(fplog, &top->idef, 0, mtop.natoms, FALSE, FALSE);
784 if (ftp2bSet(efNDX, NFILE, fnm))
786 rd_index(ftp2fn(efNDX, NFILE, fnm), 1, &isize, &index, &grpname);
787 xvg = xvgropen(opt2fn("-dr", NFILE, fnm), "Individual Restraints", "Time (ps)",
789 snew(vvindex, isize);
791 for (i = 0; (i < isize); i++)
795 sprintf(leg[i], "index %d", index[i]);
797 xvgr_legend(xvg, isize, (const char**)leg, oenv);
805 init_disres(fplog, &mtop, &ir, NULL, &fcd, NULL, FALSE);
807 natoms = read_first_x(oenv, &status, ftp2fn(efTRX, NFILE, fnm), &t, &x, box);
810 init_dr_res(&dr, fcd.disres.nres);
811 if (opt2bSet("-c", NFILE, fnm))
813 clust = cluster_index(fplog, opt2fn("-c", NFILE, fnm));
814 snew(dr_clust, clust->clust->nr+1);
815 for (i = 0; (i <= clust->clust->nr); i++)
817 init_dr_res(&dr_clust[i], fcd.disres.nres);
822 out = xvgropen(opt2fn("-ds", NFILE, fnm),
823 "Sum of Violations", "Time (ps)", "nm", oenv);
824 aver = xvgropen(opt2fn("-da", NFILE, fnm),
825 "Average Violation", "Time (ps)", "nm", oenv);
826 numv = xvgropen(opt2fn("-dn", NFILE, fnm),
827 "# Violations", "Time (ps)", "#", oenv);
828 maxxv = xvgropen(opt2fn("-dm", NFILE, fnm),
829 "Largest Violation", "Time (ps)", "nm", oenv);
832 mdatoms = init_mdatoms(fplog, &mtop, ir.efep != efepNO);
833 atoms2md(&mtop, &ir, 0, NULL, mtop.natoms, mdatoms);
834 update_mdatoms(mdatoms, ir.fepvals->init_lambda);
836 if (ir.ePBC != epbcNONE)
838 gpbc = gmx_rmpbc_init(&top->idef, ir.ePBC, natoms);
844 if (ir.ePBC != epbcNONE)
846 if (ir.bPeriodicMols)
848 set_pbc(&pbc, ir.ePBC, box);
852 gmx_rmpbc(gpbc, natoms, box, x);
858 if (j > clust->maxframe)
860 gmx_fatal(FARGS, "There are more frames in the trajectory than in the cluster index file. t = %8f\n", t);
862 my_clust = clust->inv_clust[j];
863 range_check(my_clust, 0, clust->clust->nr);
864 check_viol(fplog, &(top->idef.il[F_DISRES]),
866 x, f, pbc_null, g, dr_clust, my_clust, isize, index, vvindex, &fcd);
870 check_viol(fplog, &(top->idef.il[F_DISRES]),
872 x, f, pbc_null, g, &dr, 0, isize, index, vvindex, &fcd);
876 reset_x(atoms->nr, ind_fit, atoms->nr, NULL, x, w_rls);
877 do_fit(atoms->nr, w_rls, x, x);
880 /* Store the first frame of the trajectory as 'characteristic'
881 * for colouring with violations.
883 for (kkk = 0; (kkk < atoms->nr); kkk++)
885 copy_rvec(x[kkk], xav[kkk]);
893 fprintf(xvg, "%10g", t);
894 for (i = 0; (i < isize); i++)
896 fprintf(xvg, " %10g", vvindex[i]);
900 fprintf(out, "%10g %10g\n", t, dr.sumv);
901 fprintf(aver, "%10g %10g\n", t, dr.averv);
902 fprintf(maxxv, "%10g %10g\n", t, dr.maxv);
903 fprintf(numv, "%10g %10d\n", t, dr.nv);
907 while (read_next_x(oenv, status, &t, x, box));
909 if (ir.ePBC != epbcNONE)
911 gmx_rmpbc_done(gpbc);
916 dump_clust_stats(fplog, fcd.disres.nres, &(top->idef.il[F_DISRES]),
917 top->idef.iparams, clust->clust, dr_clust,
918 clust->grpname, isize, index);
922 dump_stats(fplog, j, fcd.disres.nres, &(top->idef.il[F_DISRES]),
923 top->idef.iparams, &dr, isize, index,
924 bPDB ? atoms : NULL);
927 write_sto_conf(opt2fn("-q", NFILE, fnm),
928 "Coloured by average violation in Angstrom",
929 atoms, xav, NULL, ir.ePBC, box);
931 dump_disre_matrix(opt2fn_null("-x", NFILE, fnm), &dr, fcd.disres.nres,
932 j, &top->idef, &mtop, max_dr, nlevels, bThird);
940 do_view(oenv, opt2fn("-dr", NFILE, fnm), "-nxy");
942 do_view(oenv, opt2fn("-dn", NFILE, fnm), "-nxy");
943 do_view(oenv, opt2fn("-da", NFILE, fnm), "-nxy");
944 do_view(oenv, opt2fn("-ds", NFILE, fnm), "-nxy");
945 do_view(oenv, opt2fn("-dm", NFILE, fnm), "-nxy");
948 gmx_log_close(fplog);