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++) {
54 for(f=0; f<nft; f++) {
55 if (idef->functype[i] == ft[f])
63 static void fill_ft_ind(int nft,int *ft,t_idef *idef,
64 int ft_ind[],char *grpnames[])
69 /* Loop over all the function types in the topology */
70 for(i=0; (i<idef->ntypes); i++) {
72 /* Check all the selected function types */
73 for(f=0; f<nft; f++) {
75 if (idef->functype[i] == ftype) {
79 sprintf(buf,"Theta=%.1f_%.2f",idef->iparams[i].harmonic.rA,
80 idef->iparams[i].harmonic.krA);
83 sprintf(buf,"Cos_th=%.1f_%.2f",idef->iparams[i].harmonic.rA,
84 idef->iparams[i].harmonic.krA);
87 sprintf(buf,"UB_th=%.1f_%.2f2f",idef->iparams[i].u_b.theta,
88 idef->iparams[i].u_b.ktheta);
90 case F_QUARTIC_ANGLES:
91 sprintf(buf,"Q_th=%.1f_%.2f_%.2f",idef->iparams[i].qangle.theta,
92 idef->iparams[i].qangle.c[0],idef->iparams[i].qangle.c[1]);
95 sprintf(buf,"Table=%d_%.2f",idef->iparams[i].tab.table,
96 idef->iparams[i].tab.kA);
99 sprintf(buf,"Phi=%.1f_%d_%.2f",idef->iparams[i].pdihs.phiA,
100 idef->iparams[i].pdihs.mult,idef->iparams[i].pdihs.cpA);
103 sprintf(buf,"Xi=%.1f_%.2f",idef->iparams[i].harmonic.rA,
104 idef->iparams[i].harmonic.krA);
107 sprintf(buf,"RB-A1=%.2f",idef->iparams[i].rbdihs.rbcA[1]);
110 gmx_fatal(FARGS,"Unsupported function type '%s' selected",
111 interaction_function[ftype].longname);
113 grpnames[ind]=strdup(buf);
120 static void fill_ang(int nft,int *ft,int fac,
121 int nr[],int *index[],int ft_ind[],t_topology *top,
124 int f,ftype,i,j,indg,nr_fac;
132 atom = top->atoms.atom;
134 for(f=0; f<nft; f++) {
136 ia = idef->il[ftype].iatoms;
137 for(i=0; (i<idef->il[ftype].nr); ) {
138 indg = ft_ind[ia[0]];
140 gmx_incons("Routine fill_ang");
143 for(j=0; j<fac; j++) {
144 if (atom[ia[1+j]].m < 1.5)
149 for(j=0; j<fac; j++) {
150 if (atom[ia[1+j]].m < 1.5 && fabs(atom[ia[1+j]].q) < hq)
155 if (nr[indg] % 1000 == 0) {
156 srenew(index[indg],fac*(nr[indg]+1000));
158 nr_fac = fac*nr[indg];
159 for(j=0; (j<fac); j++)
160 index[indg][nr_fac+j] = ia[j+1];
163 ia += interaction_function[ftype].nratoms+1;
164 i += interaction_function[ftype].nratoms+1;
169 static int *select_ftype(const char *opt,int *nft,int *mult)
175 for(ftype=0; ftype<F_NRE; ftype++) {
176 if (interaction_function[ftype].flags & IF_ATYPE ||
177 ftype == F_TABANGLES) {
205 int main(int argc,char *argv[])
207 static const char *desc[] = {
208 "mk_angndx makes an index file for calculation of",
209 "angle distributions etc. It uses a run input file ([TT].tpx[tt]) for the",
210 "definitions of the angles, dihedrals etc."
212 static const char *opt[] = { NULL, "angle", "dihedral", "improper", "ryckaert-bellemans", NULL };
216 { "-type", FALSE, etENUM, {opt},
218 { "-hyd", FALSE, etBOOL, {&bH},
219 "Include angles with atoms with mass < 1.5" },
220 { "-hq", FALSE, etREAL, {&hq},
221 "Ignore angles with atoms with mass < 1.5 and |q| < hq" }
228 int nft=0,*ft,mult=0;
234 { efTPX, NULL, NULL, ffREAD },
235 { efNDX, NULL, "angle", ffWRITE }
237 #define NFILE asize(fnm)
239 CopyRight(stderr,argv[0]);
240 parse_common_args(&argc,argv,0,NFILE,fnm,asize(pa),pa,
241 asize(desc),desc,0,NULL,&oenv);
244 ft = select_ftype(opt[0],&nft,&mult);
246 top = read_top(ftp2fn(efTPX,NFILE,fnm),NULL);
248 ntype = calc_ntype(nft,ft,&(top->idef));
249 snew(grpnames,ntype);
250 snew(ft_ind,top->idef.ntypes);
251 fill_ft_ind(nft,ft,&top->idef,ft_ind,grpnames);
255 fill_ang(nft,ft,mult,nr,index,ft_ind,top,!bH,hq);
257 out=ftp2FILE(efNDX,NFILE,fnm,"w");
258 for(i=0; (i<ntype); i++) {
260 fprintf(out,"[ %s ]\n",grpnames[i]);
261 for(j=0; (j<nr[i]*mult); j++) {
262 fprintf(out," %5d",index[i][j]+1);