Refactor CUDA kernel use of tables
[gromacs/AngularHB.git] / src / contrib / do_shift.c
blobfce2869ebb926dfdd4b76b03f1098ab510307429
1 /*
2 *
3 * This source code is part of
4 *
5 * G R O M A C S
6 *
7 * GROningen MAchine for Chemical Simulations
8 *
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
32 * And Hey:
33 * Groningen Machine for Chemical Simulation
35 #ifdef HAVE_CONFIG_H
36 #include <config.h>
37 #endif
39 #include <stdlib.h>
40 #include "errno.h"
41 #include "typedefs.h"
42 #include "gromacs/utility/cstringutil.h"
43 #include "macros.h"
44 #include "gromacs/utility/smalloc.h"
45 #include "gromacs/commandline/pargs.h"
46 #include "copyrite.h"
47 #include "gromacs/fileio/confio.h"
48 #include "gromacs/utility/fatalerror.h"
49 #include "gromacs/fileio/xvgr.h"
50 #include "gstat.h"
51 #include "gromacs/topology/index.h"
52 #include "gromacs/fileio/pdbio.h"
54 void cat(FILE *out,char *fn,real t)
56 FILE *in;
57 char *ptr,buf[256];
58 int anr,rnr;
59 char anm[24],rnm[24];
60 double f1,f2,f3,f4,f5,f6;
62 in=gmx_ffopen(fn,"r");
63 while ((ptr=fgets2(buf,255,in)) != NULL) {
64 sscanf(buf,"%d%d%s%s%lf%lf%lf%lf%lf%lf",
65 &anr,&rnr,rnm,anm,&f1,&f2,&f3,&f4,&f5,&f6);
66 fprintf(out,"%8g %10g %10g %10g %10g %10g %10g %s%d-%s%d\n",
67 t,f6,f1,f2,f3,f4,f5,rnm,rnr,anm,anr);
69 /*if ((int)strlen(buf) > 0)
70 fprintf(out,"%s\n",buf);*/
71 fflush(out);
72 gmx_ffclose(in);
75 int main(int argc,char *argv[])
77 static char *desc[] = {
78 "[TT]do_shift[tt] reads a trajectory file and computes the chemical",
79 "shift for each time frame (or every [BB]dt[bb] ps) by",
80 "calling the 'total' program. If you do not have the total program,",
81 "get it. do_shift assumes that the total executable is in",
82 "[TT]/home/mdgroup/total/total[tt]. If that is not the case, then you should",
83 "set an environment variable [BB]GMX_TOTAL[bb] as in: [PAR]",
84 "[TT]setenv GMX_TOTAL /usr/local/bin/total[tt][PAR]",
85 "where the right hand side should point to the total executable.[PAR]",
86 "Output is printed in files [TT]shift.out[tt] where t is the time of the frame.[PAR]",
87 "The program also needs an input file called [BB]random.dat[bb] which",
88 "contains the random coil chemical shifts of all protons."
90 static real dt=0.0;
91 t_pargs pa[] = {
92 { "-dt", FALSE, etREAL, { &dt }, "Time interval between frames." }
94 static char *bugs[] = {
95 "The program is very slow"
97 static char *OXYGEN="O";
98 FILE *out,*tot,*fp;
99 t_topology *top;
100 t_atoms *atoms;
101 int status,nres;
102 real t,nt;
103 int i,natoms,nframe=0;
104 matrix box;
105 int gnx;
106 char *grpnm,*randf;
107 atom_id *index;
108 rvec *x,*x_s;
109 char pdbfile[32],tmpfile[32];
110 char total[256],*dptr;
111 t_filenm fnm[] = {
112 { efTRX, "-f", NULL, ffREAD },
113 { efTPX, NULL, NULL, ffREAD },
114 { efNDX, NULL, NULL, ffREAD },
115 { efOUT, "-o", "shift", ffWRITE },
116 { efDAT, "-d", "random", ffREAD }
118 char *leg[] = { "shift","ring","anisCO","anisCN","sigmaE","sum" };
119 #define NFILE asize(fnm)
121 CopyRight(stdout,argv[0]);
122 parse_common_args(&argc,argv,PCA_CAN_TIME,NFILE,fnm,
123 asize(pa),pa,asize(desc),desc,asize(bugs),bugs);
125 top=read_top(ftp2fn(efTPX,NFILE,fnm));
126 atoms=&(top->atoms);
127 nres=atoms->nres;
128 for(i=0; (i<atoms->nr); i++)
129 if ((strcmp(*atoms->atomname[i],"O1") == 0) ||
130 (strcmp(*atoms->atomname[i],"O2") == 0) ||
131 (strcmp(*atoms->atomname[i],"OXT") == 0) ||
132 (strcmp(*atoms->atomname[i],"OT") == 0))
133 atoms->atomname[i]=&OXYGEN;
134 rd_index(ftp2fn(efNDX,NFILE,fnm),1,&gnx,&index,&grpnm);
136 snew(x_s,atoms->nr);
138 strcpy(pdbfile,"dsXXXXXX");
139 gmx_tmpnam(pdbfile);
140 strcpy(tmpfile,"dsXXXXXX");
141 gmx_tmpnam(tmpfile);
142 fprintf(stderr,"pdbfile = %s\ntmpfile = %s\n",pdbfile,tmpfile);
144 if ((dptr=getenv("GMX_TOTAL")) == NULL)
145 dptr="/home/mdgroup/total/total";
146 sprintf(total,"%s > /dev/null",dptr);
147 fprintf(stderr,"total cmd='%s'\n",total);
148 randf=ftp2fn(efDAT,NFILE,fnm);
150 natoms=read_first_x(&status,ftp2fn(efTRX,NFILE,fnm),&t,&x,box);
151 if (natoms != atoms->nr)
152 gmx_fatal(FARGS,"Trajectory does not match topology!");
153 out=ftp2FILE(efOUT,NFILE,fnm,"w");
154 xvgr_legend(out,asize(leg),leg);
155 nt=t;
156 do {
157 if (t >= nt) {
158 rm_pbc(&(top->idef),top->atoms.nr,box,x,x_s);
159 fp=gmx_ffopen(pdbfile,"w");
160 write_pdbfile_indexed(fp,"Generated by do_shift",
161 atoms,x_s,box,0,-1,gnx,index);
162 gmx_ffclose(fp);
164 if ((tot=popen(total,"w")) == NULL)
165 perror("opening pipe to total");
166 fprintf(tot,"%s\n",pdbfile);
167 fprintf(tot,"%s\n",tmpfile);
168 fprintf(tot,"3\n");
169 fprintf(tot,"N\n");
170 fprintf(tot,"%s\n",randf);
171 fprintf(tot,"N\n");
172 fprintf(tot,"N\n");
173 if (pclose(tot) != 0)
174 perror("closing pipe to total");
175 cat(out,tmpfile,t);
176 remove(pdbfile);
177 remove(tmpfile);
178 nt+=dt;
179 nframe++;
181 } while(read_next_x(status,&t,natoms,x,box));
182 close_trj(status);
183 gmx_ffclose(out);
185 gmx_thanx(stderr);
187 return 0;