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
47 #include "gmx_fatal.h"
49 static int calc_ntype(int nft, int *ft, t_idef *idef)
53 for (i = 0; (i < idef->ntypes); i++)
55 for (f = 0; f < nft; f++)
57 if (idef->functype[i] == ft[f])
67 static void fill_ft_ind(int nft, int *ft, t_idef *idef,
68 int ft_ind[], char *grpnames[])
71 int i, f, ftype, ind = 0;
73 /* Loop over all the function types in the topology */
74 for (i = 0; (i < idef->ntypes); i++)
77 /* Check all the selected function types */
78 for (f = 0; f < nft; f++)
81 if (idef->functype[i] == ftype)
87 sprintf(buf, "Theta=%.1f_%.2f", idef->iparams[i].harmonic.rA,
88 idef->iparams[i].harmonic.krA);
91 sprintf(buf, "Cos_th=%.1f_%.2f", idef->iparams[i].harmonic.rA,
92 idef->iparams[i].harmonic.krA);
95 sprintf(buf, "UB_th=%.1f_%.2f2f", idef->iparams[i].u_b.thetaA,
96 idef->iparams[i].u_b.kthetaA);
98 case F_QUARTIC_ANGLES:
99 sprintf(buf, "Q_th=%.1f_%.2f_%.2f", idef->iparams[i].qangle.theta,
100 idef->iparams[i].qangle.c[0], idef->iparams[i].qangle.c[1]);
103 sprintf(buf, "Table=%d_%.2f", idef->iparams[i].tab.table,
104 idef->iparams[i].tab.kA);
107 sprintf(buf, "Phi=%.1f_%d_%.2f", idef->iparams[i].pdihs.phiA,
108 idef->iparams[i].pdihs.mult, idef->iparams[i].pdihs.cpA);
111 sprintf(buf, "Xi=%.1f_%.2f", idef->iparams[i].harmonic.rA,
112 idef->iparams[i].harmonic.krA);
115 sprintf(buf, "RB-A1=%.2f", idef->iparams[i].rbdihs.rbcA[1]);
118 gmx_fatal(FARGS, "Unsupported function type '%s' selected",
119 interaction_function[ftype].longname);
121 grpnames[ind] = strdup(buf);
128 static void fill_ang(int nft, int *ft, int fac,
129 int nr[], int *index[], int ft_ind[], t_topology *top,
130 gmx_bool bNoH, real hq)
132 int f, ftype, i, j, indg, nr_fac;
140 atom = top->atoms.atom;
142 for (f = 0; f < nft; f++)
145 ia = idef->il[ftype].iatoms;
146 for (i = 0; (i < idef->il[ftype].nr); )
148 indg = ft_ind[ia[0]];
151 gmx_incons("Routine fill_ang");
156 for (j = 0; j < fac; j++)
158 if (atom[ia[1+j]].m < 1.5)
166 for (j = 0; j < fac; j++)
168 if (atom[ia[1+j]].m < 1.5 && fabs(atom[ia[1+j]].q) < hq)
176 if (nr[indg] % 1000 == 0)
178 srenew(index[indg], fac*(nr[indg]+1000));
180 nr_fac = fac*nr[indg];
181 for (j = 0; (j < fac); j++)
183 index[indg][nr_fac+j] = ia[j+1];
187 ia += interaction_function[ftype].nratoms+1;
188 i += interaction_function[ftype].nratoms+1;
193 static int *select_ftype(const char *opt, int *nft, int *mult)
195 int *ft = NULL, ftype;
200 for (ftype = 0; ftype < F_NRE; ftype++)
202 if ((interaction_function[ftype].flags & IF_ATYPE) ||
203 ftype == F_TABANGLES)
235 int gmx_mk_angndx(int argc, char *argv[])
237 static const char *desc[] = {
238 "[TT]mk_angndx[tt] makes an index file for calculation of",
239 "angle distributions etc. It uses a run input file ([TT].tpx[tt]) for the",
240 "definitions of the angles, dihedrals etc."
242 static const char *opt[] = { NULL, "angle", "dihedral", "improper", "ryckaert-bellemans", NULL };
243 static gmx_bool bH = TRUE;
246 { "-type", FALSE, etENUM, {opt},
248 { "-hyd", FALSE, etBOOL, {&bH},
249 "Include angles with atoms with mass < 1.5" },
250 { "-hq", FALSE, etREAL, {&hq},
251 "Ignore angles with atoms with mass < 1.5 and magnitude of their charge less than this value" }
258 int nft = 0, *ft, mult = 0;
264 { efTPX, NULL, NULL, ffREAD },
265 { efNDX, NULL, "angle", ffWRITE }
267 #define NFILE asize(fnm)
269 parse_common_args(&argc, argv, 0, NFILE, fnm, asize(pa), pa,
270 asize(desc), desc, 0, NULL, &oenv);
273 ft = select_ftype(opt[0], &nft, &mult);
275 top = read_top(ftp2fn(efTPX, NFILE, fnm), NULL);
277 ntype = calc_ntype(nft, ft, &(top->idef));
278 snew(grpnames, ntype);
279 snew(ft_ind, top->idef.ntypes);
280 fill_ft_ind(nft, ft, &top->idef, ft_ind, grpnames);
284 fill_ang(nft, ft, mult, nr, index, ft_ind, top, !bH, hq);
286 out = ftp2FILE(efNDX, NFILE, fnm, "w");
287 for (i = 0; (i < ntype); i++)
291 fprintf(out, "[ %s ]\n", grpnames[i]);
292 for (j = 0; (j < nr[i]*mult); j++)
294 fprintf(out, " %5d", index[i][j]+1);