3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
9 * VERSION 3.3.99_development_20071104
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-2006, 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 * Groningen Machine for Chemical Simulation
51 #include "gmx_fatal.h"
57 void cat(FILE *out
,char *fn
,real t
)
63 double f1
,f2
,f3
,f4
,f5
,f6
;
66 while ((ptr
=fgets2(buf
,255,in
)) != NULL
) {
67 sscanf(buf
,"%d%d%s%s%lf%lf%lf%lf%lf%lf",
68 &anr
,&rnr
,rnm
,anm
,&f1
,&f2
,&f3
,&f4
,&f5
,&f6
);
69 fprintf(out
,"%8g %10g %10g %10g %10g %10g %10g %s%d-%s%d\n",
70 t
,f6
,f1
,f2
,f3
,f4
,f5
,rnm
,rnr
,anm
,anr
);
72 /*if ((int)strlen(buf) > 0)
73 fprintf(out,"%s\n",buf);*/
78 int main(int argc
,char *argv
[])
80 static char *desc
[] = {
81 "[TT]do_shift[tt] reads a trajectory file and computes the chemical",
82 "shift for each time frame (or every [BB]dt[bb] ps) by",
83 "calling the 'total' program. If you do not have the total program,",
84 "get it. do_shift assumes that the total executable is in",
85 "[TT]/home/mdgroup/total/total[tt]. If that is not the case, then you should",
86 "set an environment variable [BB]TOTAL[bb] as in: [PAR]",
87 "[TT]setenv TOTAL /usr/local/bin/total[tt][PAR]",
88 "where the right hand side should point to the total executable.[PAR]",
89 "Output is printed in files [TT]shift.out[tt] where t is the time of the frame.[PAR]",
90 "The program also needs an input file called [BB]random.dat[bb] which",
91 "contains the random coil chemical shifts of all protons."
95 { "-dt", FALSE
, etREAL
, { &dt
}, "Time interval between frames." }
97 static char *bugs
[] = {
98 "The program is very slow"
100 static char *OXYGEN
="O";
106 int i
,natoms
,nframe
=0;
112 char pdbfile
[32],tmpfile
[32];
113 char total
[256],*dptr
;
115 { efTRX
, "-f", NULL
, ffREAD
},
116 { efTPX
, NULL
, NULL
, ffREAD
},
117 { efNDX
, NULL
, NULL
, ffREAD
},
118 { efOUT
, "-o", "shift", ffWRITE
},
119 { efDAT
, "-d", "random", ffREAD
}
121 char *leg
[] = { "shift","ring","anisCO","anisCN","sigmaE","sum" };
122 #define NFILE asize(fnm)
124 CopyRight(stdout
,argv
[0]);
125 parse_common_args(&argc
,argv
,PCA_CAN_TIME
| PCA_BE_NICE
,NFILE
,fnm
,
126 asize(pa
),pa
,asize(desc
),desc
,asize(bugs
),bugs
);
128 top
=read_top(ftp2fn(efTPX
,NFILE
,fnm
));
131 for(i
=0; (i
<atoms
->nr
); i
++)
132 if ((strcmp(*atoms
->atomname
[i
],"O1") == 0) ||
133 (strcmp(*atoms
->atomname
[i
],"O2") == 0) ||
134 (strcmp(*atoms
->atomname
[i
],"OXT") == 0) ||
135 (strcmp(*atoms
->atomname
[i
],"OT") == 0))
136 atoms
->atomname
[i
]=&OXYGEN
;
137 rd_index(ftp2fn(efNDX
,NFILE
,fnm
),1,&gnx
,&index
,&grpnm
);
141 strcpy(pdbfile
,"dsXXXXXX");
143 strcpy(tmpfile
,"dsXXXXXX");
145 fprintf(stderr
,"pdbfile = %s\ntmpfile = %s\n",pdbfile
,tmpfile
);
147 if ((dptr
=getenv("TOTAL")) == NULL
)
148 dptr
="/home/mdgroup/total/total";
149 sprintf(total
,"%s > /dev/null",dptr
);
150 fprintf(stderr
,"total cmd='%s'\n",total
);
151 randf
=ftp2fn(efDAT
,NFILE
,fnm
);
153 natoms
=read_first_x(&status
,ftp2fn(efTRX
,NFILE
,fnm
),&t
,&x
,box
);
154 if (natoms
!= atoms
->nr
)
155 gmx_fatal(FARGS
,"Trajectory does not match topology!");
156 out
=ftp2FILE(efOUT
,NFILE
,fnm
,"w");
157 xvgr_legend(out
,asize(leg
),leg
);
161 rm_pbc(&(top
->idef
),top
->atoms
.nr
,box
,x
,x_s
);
162 fp
=ffopen(pdbfile
,"w");
163 write_pdbfile_indexed(fp
,"Generated by do_shift",
164 atoms
,x_s
,box
,0,-1,gnx
,index
);
167 if ((tot
=popen(total
,"w")) == NULL
)
168 perror("opening pipe to total");
169 fprintf(tot
,"%s\n",pdbfile
);
170 fprintf(tot
,"%s\n",tmpfile
);
173 fprintf(tot
,"%s\n",randf
);
176 if (pclose(tot
) != 0)
177 perror("closing pipe to total");
184 } while(read_next_x(status
,&t
,natoms
,x
,box
));