GLYLIB
0.3.0b
|
00001 /* Function file read_amber8_mden.c begun on 20070911 by BLFoley 00002 00003 This function reads information contained in fileset F. The file 00004 pointer in the fileset must already be pointing to a file open with 00005 at least read capability. The file should be a mden (detailed 00006 energy) output file generated by sander. This function was written 00007 for the AMBER 8 mden file format and should work with any other 00008 version that uses the same format. 00009 00010 The other arguments are NENT, the number of entries to read, and 00011 a pointer to a set of NENT scharsets (ENT). These charsets contain 00012 the list of keywords for the type of data to be stored in the 00013 return array. NOTE!!! The entries in ENT can be in any order, 00014 but in the return array the information will be IN THE ORDER IN WHICH 00015 IT APPEARS IN THE MDEN FILE. The entries in ENT must exactly match 00016 the header keywords in an mden file. 00017 00018 The return array will be two-dimensional, of size NENT*NMAX. The 00019 entries will be ordered such that the array can be accessed by 00020 ARYP[column][row]. 00021 00022 If NENT is set to zero (and ENT is an empty char**), then the program 00023 will read the contents of the entire mden file. 00024 00025 Example: 00026 00027 You want to read in Nsteps, Temp, pres_Z and E_hbon from the mden file. 00028 The mden file contains 430 entries. 00029 00030 Call the function: 00031 00032 double *myarray; 00033 char **ent; 00034 fileset MyMDEN; 00035 00036 ent=(char*)calloc(4,sizeof(char*)); 00037 ent[0]=strdup("Nsteps"); 00038 ent[1]=strdup("Temp"); 00039 ent[2]=strdup("pres_Z"); 00040 ent[3]=strdup("E_hbon"); 00041 MyMDEN.N=strdup("My_sander_eout.en"); 00042 MyMDEN.F=myfopen(MyMDEN.N,"r"); 00043 00044 myarray=read_amber8_mden(4,ent,MyMDEN); 00045 00046 The contents of myarray will look like: 00047 00048 column contains --> [index] Nsteps Temp pres_Z E_hbon 00049 mden entry 00050 | [0] 1 320 0.98 42.6 00051 V [1] 300 365 1.54 54.3 00052 [...] ... ... ... ... 00053 [429] 129000 319 1.25 48.9 00054 00055 The mden file looks like: 00056 00057 L0 Nsteps time(ps) Etot EKinetic 00058 L1 Temp T_solute T_solv Pres_scal_solu 00059 L2 Pres_scal_solv BoxX BoxY BoxZ 00060 L3 volume pres_X pres_Y pres_Z 00061 L4 Pressure EKCoM_x EKCoM_y EKCoM_z 00062 L5 EKComTot VIRIAL_x VIRIAL_y VIRIAL_z 00063 L6 VIRIAL_tot E_pot E_vdw E_el 00064 L7 E_hbon E_bon E_angle E_dih 00065 L8 E_14vdw E_14el E_const E_pol 00066 L9 AV_permMoment AV_indMoment AV_totMoment Density dV/dlambda 00067 L0 1 0.2330001000E+04 -.6429480312E+05 0.1684971086E+05 00068 L1 0.2977702692E+03 0.2977702692E+03 -.6434448936E+11 0.1000000000E+01 00069 L2 0.1000000000E+01 0.6062081490E+02 0.7469235400E+02 0.6063588770E+02 00070 L3 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 00071 L4 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 00072 L5 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 00073 L6 0.0000000000E+00 -.8114451398E+05 0.9217360369E+04 -.1082986506E+06 00074 L7 0.0000000000E+00 0.8145159202E+03 0.2010802568E+04 0.9283273104E+03 00075 L8 0.1272671427E+04 0.1291045899E+05 0.0000000000E+00 0.0000000000E+00 00076 L9 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 0.0000000000E+00 -.3245291868E+02 00077 (etc) 00078 */ 00079 #include <AMBER/read_amber8_mden.h> 00080 //#include "../inc/read_amber8_mden.h" 00081 //#include "mylib.h" 00082 00083 double **read_amber8_mden(int NENT,char **ENT,int NDATA,fileset F){ 00084 int ra=0,rb=0,rc=0,NMAX=0,S[41]; 00085 double **ARYP,toss; // pointer to return array 00086 char line[501]; 00087 char *MDEN_ENT[41]; 00088 00089 //printf("Here, read mden TOP\n"); 00090 00091 // initialize some stuff 00092 MDEN_ENT[0]=strdup("Nsteps"); 00093 MDEN_ENT[1]=strdup("time(ps)"); 00094 MDEN_ENT[2]=strdup("Etot"); 00095 MDEN_ENT[3]=strdup("EKinetic"); 00096 MDEN_ENT[4]=strdup("Temp"); 00097 MDEN_ENT[5]=strdup("T_solute"); 00098 MDEN_ENT[6]=strdup("T_solv"); 00099 MDEN_ENT[7]=strdup("Pres_scal_solu"); 00100 MDEN_ENT[8]=strdup("Pres_scal_solv"); 00101 MDEN_ENT[9]=strdup("BoxX"); 00102 00103 MDEN_ENT[10]=strdup("BoxY"); 00104 MDEN_ENT[11]=strdup("BoxZ"); 00105 MDEN_ENT[12]=strdup("volume"); 00106 MDEN_ENT[13]=strdup("pres_X"); 00107 MDEN_ENT[14]=strdup("pres_Y"); 00108 MDEN_ENT[15]=strdup("pres_Z"); 00109 MDEN_ENT[16]=strdup("Pressure"); 00110 MDEN_ENT[17]=strdup("EKCoM_x"); 00111 MDEN_ENT[18]=strdup("EKCoM_y"); 00112 MDEN_ENT[19]=strdup("EKCoM_z"); 00113 00114 MDEN_ENT[20]=strdup("EKComTot"); 00115 MDEN_ENT[21]=strdup("VIRIAL_x"); 00116 MDEN_ENT[22]=strdup("VIRIAL_y"); 00117 MDEN_ENT[23]=strdup("VIRIAL_z"); 00118 MDEN_ENT[24]=strdup("VIRIAL_tot"); 00119 MDEN_ENT[25]=strdup("E_pot"); 00120 MDEN_ENT[26]=strdup("E_vdw"); 00121 MDEN_ENT[27]=strdup("E_el"); 00122 MDEN_ENT[28]=strdup("E_hbon"); 00123 MDEN_ENT[29]=strdup("E_bon"); 00124 00125 MDEN_ENT[30]=strdup("E_angle"); 00126 MDEN_ENT[31]=strdup("E_dih"); 00127 MDEN_ENT[32]=strdup("E_14vdw"); 00128 MDEN_ENT[33]=strdup("E_14el"); 00129 MDEN_ENT[34]=strdup("E_const"); 00130 MDEN_ENT[35]=strdup("E_pol"); 00131 MDEN_ENT[36]=strdup("AV_permMoment"); 00132 MDEN_ENT[37]=strdup("AV_indMoment"); 00133 MDEN_ENT[38]=strdup("AV_totMoment"); 00134 MDEN_ENT[39]=strdup("Density"); 00135 MDEN_ENT[40]=strdup("dV/dlambda"); 00136 00137 if(NENT==0){ 00138 NENT=41; 00139 ENT=MDEN_ENT; 00140 } 00141 ARYP=(double**)calloc(NENT,sizeof(double*)); 00142 //printf("Here, read mden 1\n"); 00143 // see which entries we need to grab 00144 rb=0; 00145 for(ra=0;ra<41;ra++){ 00146 if(strcmp(MDEN_ENT[ra],ENT[rb])==0){ 00147 S[ra]=0; 00148 rb++;} 00149 else{S[ra]=1;} 00150 } 00151 if(rb!=NENT) mywhine("Mismatch in NENT & rb in mden read"); 00152 //printf("Here, read mden 2\n"); 00153 00154 F.F=myfreopen(F.N,"r",F.F); 00155 // look through the file to see how much memory to allocate 00156 rewind(F.F); 00157 //fgets(line,500,F.F); // get zeroth line 00158 ra=0; 00159 while(fgets(line,500,F.F)!=NULL){ if(strncmp(line,"L9",2)==0){ ra++; } } 00160 NMAX=ra-1; 00161 //printf("NMAX is %d\n",NMAX); 00162 if(NMAX<NDATA){mywhine("read_amber8_mden mismatch in number of expected data points");} 00163 for(ra=0;ra<NENT;ra++){ ARYP[ra]=(double*)calloc(NMAX,sizeof(double)); } 00164 rewind(F.F); 00165 00166 //// grab first 51 file entries and toss 00167 //for(ra=0;ra<51;ra++){fscanf(F.F,"%s",line);} 00168 // grab first 51 file entries and toss 00169 for(ra=0;ra<10;ra++){fgets(line,500,F.F);} 00170 // step through file, recording entries as needed 00171 for(ra=0;ra<NMAX;ra++){ 00172 rc=0; 00173 fscanf(F.F,"%s",line); // toss "L0" 00174 if(strcmp(line,"L0")!=0){ 00175 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00176 exit(1); 00177 } 00178 for(rb=0;rb<4;rb++){ 00179 fscanf(F.F,"%lf",&toss); 00180 if(S[rb]==0){ 00181 ARYP[rc][ra]=toss; 00182 rc++; 00183 } 00184 } 00185 fscanf(F.F,"%s",line); // toss "L1" 00186 if(strcmp(line,"L1")!=0){ 00187 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00188 exit(1); 00189 } 00190 for(rb=4;rb<8;rb++){ 00191 fscanf(F.F,"%lf",&toss); 00192 if(S[rb]==0){ 00193 ARYP[rc][ra]=toss; 00194 rc++; 00195 } 00196 } 00197 fscanf(F.F,"%s",line); // toss "L2" 00198 if(strcmp(line,"L2")!=0){ 00199 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00200 exit(1); 00201 } 00202 for(rb=8;rb<12;rb++){ 00203 fscanf(F.F,"%lf",&toss); 00204 if(S[rb]==0){ 00205 ARYP[rc][ra]=toss; 00206 rc++; 00207 } 00208 } 00209 fscanf(F.F,"%s",line); // toss "L3" 00210 if(strcmp(line,"L3")!=0){ 00211 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00212 exit(1); 00213 } 00214 for(rb=12;rb<16;rb++){ 00215 fscanf(F.F,"%lf",&toss); 00216 if(S[rb]==0){ 00217 ARYP[rc][ra]=toss; 00218 rc++; 00219 } 00220 } 00221 fscanf(F.F,"%s",line); // toss "L4" 00222 if(strcmp(line,"L4")!=0){ 00223 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00224 exit(1); 00225 } 00226 for(rb=16;rb<20;rb++){ 00227 fscanf(F.F,"%lf",&toss); 00228 if(S[rb]==0){ 00229 ARYP[rc][ra]=toss; 00230 rc++; 00231 } 00232 } 00233 fscanf(F.F,"%s",line); // toss "L5" 00234 if(strcmp(line,"L5")!=0){ 00235 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00236 exit(1); 00237 } 00238 for(rb=20;rb<24;rb++){ 00239 fscanf(F.F,"%lf",&toss); 00240 if(S[rb]==0){ 00241 ARYP[rc][ra]=toss; 00242 rc++; 00243 } 00244 } 00245 fscanf(F.F,"%s",line); // toss "L6" 00246 if(strcmp(line,"L6")!=0){ 00247 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00248 exit(1); 00249 } 00250 for(rb=24;rb<28;rb++){ 00251 fscanf(F.F,"%lf",&toss); 00252 if(S[rb]==0){ 00253 ARYP[rc][ra]=toss; 00254 rc++; 00255 } 00256 } 00257 fscanf(F.F,"%s",line); // toss "L7" 00258 if(strcmp(line,"L7")!=0){ 00259 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00260 exit(1); 00261 } 00262 for(rb=28;rb<32;rb++){ 00263 fscanf(F.F,"%lf",&toss); 00264 if(S[rb]==0){ 00265 ARYP[rc][ra]=toss; 00266 rc++; 00267 } 00268 } 00269 fscanf(F.F,"%s",line); // toss "L8" 00270 if(strcmp(line,"L8")!=0){ 00271 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00272 exit(1); 00273 } 00274 for(rb=32;rb<36;rb++){ 00275 fscanf(F.F,"%lf",&toss); 00276 if(S[rb]==0){ 00277 ARYP[rc][ra]=toss; 00278 rc++; 00279 } 00280 } 00281 fscanf(F.F,"%s",line); // toss "L9" 00282 if(strcmp(line,"L9")!=0){ 00283 printf("In read_amber8_mden, expect >>L0<< found >>%s<<\n",line); 00284 exit(1); 00285 } 00286 for(rb=36;rb<41;rb++){ 00287 fscanf(F.F,"%lf",&toss); 00288 if(S[rb]==0){ 00289 ARYP[rc][ra]=toss; 00290 rc++; 00291 } 00292 } 00293 } 00294 00295 return ARYP; 00296 }