********>Master Bugfix File - AmberTools bugfixes for version 1.2 Programs: AmberTools Description: This is the master bugfix file for AmberTools-1.2. It contains all of the bugfixes released against version 1.2 in a single patch script. Usage: Save this file in your $AMBERHOME directory and then apply this patch file to your distribution as follows: cd $AMBERHOME patch -p0 -N bugfix 1. Author: Scott Brozell Date: 31 July 2008 Programs: antechamber Severity: minor Description: The ala2 charmm test parameter save output file was not updated for the phase change from 0.0 to 2.0. The patch below eliminates an incorrect possible failure in antechamber/charmm/ala2.prm.dif Fix: Apply the following fix to amber10/test/antechamber/charmm/ala2.prm.save --------------------------------------------------------------------------- *** test/antechamber/charmm/ala2.prm.save 15 Apr 2008 23:31:13 -0000 10.0 --- test/antechamber/charmm/ala2.prm.save 17 Jul 2008 19:34:35 -0000 10.1 *************** *** 51,58 **** c3 c3 c n 0.070 2 0.0 IMPHI ! X X c o 10.500 0 180.0 ! X X n hn 1.100 0 180.0 NONBONDED NBXMOD 5 GROUP SWITCH CDIEL - CUTNB 14.0 CTOFNB 12.0 CTONNB 10.0 EPS 1.0 E14FAC 0.83333333 WMIN 1.4 --- 51,58 ---- c3 c3 c n 0.070 2 0.0 IMPHI ! X X c o 10.500 2 180.0 ! X X n hn 1.100 2 180.0 NONBONDED NBXMOD 5 GROUP SWITCH CDIEL - CUTNB 14.0 CTOFNB 12.0 CTONNB 10.0 EPS 1.0 E14FAC 0.83333333 WMIN 1.4 --------------------------------------------------------------------------- Workarounds: Verify the antechamber charmm ala2 test results appropriately. ********> bugfix 2. Author: Andreas Svrcek-Seiler, Niel Michael Henriksen, Scott Brozell, reported by Claudio Carra Date: 8 August 2008 Programs: nab Severity: major for frcmod.parmbsc0 Description: The symptom is that NAB terminates with a failed assertion: sff2.c:1349: ephi2: Assertion `fabs(3.14159 - Phase[atyp]) < 0.01' failed. The problem is outdated code that validates dihedral phases in prmtop files. The patch below allows dihedral phases that are not integral multiples of pi. Such dihedral phases are used in the "Barcelona" changes to ff99 for nucleic acids which is most commonly accessed via leaprc.ff99bsc0. Fix: Apply the following fix to amber10/src/nab/sff2.c --------------------------------------------------------------------------- *** src/nab/sff2.c 30 Jul 2008 22:15:13 -0000 10.2 --- src/nab/sff2.c 11 Aug 2008 10:42:59 -0000 10.3 *************** *** 1342,1353 **** exit(1); } ! if (Phase[atyp] == 0.0) { ! arg = Pk[atyp]; ! } else { ! arg = -Pk[atyp]; ! assert(fabs(3.14159 - Phase[atyp]) < 0.01); ! } e = Pk[atyp] + arg * e; df = df * arg; --- 1342,1351 ---- exit(1); } ! if (fabs(Phase[atyp] - 3.142) < 0.01) ! arg = -Pk[atyp]; ! else ! arg = Pk[atyp]; e = Pk[atyp] + arg * e; df = df * arg; --------------------------------------------------------------------------- Workarounds: none. ********> bugfix 3. Author: Dan Roe Date: 29 August 2008 Programs: ptraj Description: Fixes multiple problems. 1) Ptraj did not correctly open netcdf files. 2) Ptraj misidentified Amber REMD restarts as trajectories. 3) io.c contained duplicate code in openFile() that resulted in strange behavior for bzipped files. 4) gzipped amber remd trajectories are now handled correctly. 5) remd trajectories can now have numerical extensions of arbitrary width, not just %03i. 6) Better error checking and more informative error messages for remd traj code. Affected files: ptraj.c, io.c, trajectory.h Fix: Run this patch in the $AMBERHOME/src/ptraj directory, with: $ patch -p0 -N -r patch_rejects < ------------------------------------------------------------------------------ --- src/ptraj/ptraj.c 2008-08-22 14:20:39.000000000 -0400 +++ src/ptraj/ptraj.c 2008-08-22 14:20:29.000000000 -0400 @@ -2081,9 +2081,7 @@ NCInfo->ncid = ncid; } else #endif - - err = openFile(&fp, filename, "r"); - if ( err == 0 ) { + if ( openFile(&fp, filename, "r") == 0 ) { fprintf(stdout, "trajin %s ignored; could not open file (%s)\n", filename, filename); safe_free(trajInfo->filename); diff -Naur io.c io.c --- src/ptraj/io.c 2008-08-29 11:31:32.000000000 -0400 +++ src/ptraj/io.c 2008-08-29 11:31:03.000000000 -0400 @@ -658,6 +658,7 @@ /* * Open file specified by filename with "mode" to the file pointer * fpp + * Returns 0 if an error occured, otherwise returns 1 */ int @@ -721,26 +722,9 @@ sprintf(buffer, "gzip >> %s", filename); else is_compressed = NULL; - } - - if (is_compressed != NULL) { - - /* - * trajectory is compressed, use popen() - */ - - if (prnlev > 2) - fprintf(stdout, "Opening compressed file: %s\n", buffer); - if ( ( *fpp = safe_popen(buffer, mode) ) == NULL ) { - if (prnlev > 4) - fprintf(stdout, "\nCould not open compressed file (%s) with mode (%s)\n", - buffer, mode); - return 0; - } } else if ( (is_compressed = strstr(filename, ".bz2")) != NULL && - (is_compressed[4] == (char) 0 || isspace(is_compressed[4])) ) - + (is_compressed[4] == (char) 0 || isspace(is_compressed[4])) ) { if (mode[0] == 'r') sprintf(buffer, "bunzip2 -c %s", filename); else if (mode[0] == 'w') @@ -749,16 +733,17 @@ sprintf(buffer, "bzip2 >> %s", filename); else is_compressed = NULL; - + } + if (is_compressed != NULL) { - + /* * trajectory is compressed, use popen() */ - + if (prnlev > 2) fprintf(stdout, "Opening compressed file: %s\n", buffer); - + if ( ( *fpp = safe_popen(buffer, mode) ) == NULL ) { if (prnlev > 4) fprintf(stdout, "Could not open compressed file (%s) with mode (%s)\n", diff -Naur ptraj.c ptraj.c --- src/ptraj/ptraj.c 2008-08-29 11:31:15.000000000 -0400 +++ src/ptraj/ptraj.c 2008-08-29 11:31:03.000000000 -0400 @@ -2190,11 +2190,13 @@ } /* - * Can we scan in three numbers from the first line? If so this is + * DAN ROE - Modified so that ptraj will recognize restarts from + * REMD, which have 3 numbers in the first line. + * Can we scan in four numbers from the first line? If so this is * an AMBER trajectory otherwise assume it is an AMBER restrt file */ if ( type == COORD_UNKNOWN ) { - if ( (i = sscanf(buffer2, "%f%f%f", &junk1, &junk2, &junk3)) == 3 ) + if ( (i = sscanf(buffer2, "%f%f%f%f", &junk1, &junk2, &junk3, &junk4)) == 4 ) type = COORD_AMBER_TRAJECTORY; else type = COORD_AMBER_RESTART; @@ -3380,8 +3382,8 @@ case COORD_AMBER_REMD: fprintf(stdout, " File (%s) is an AMBER REMD (new format) trajectory%s", p->filename, (p->isBox ? " (with box info)" : "")); - fprintf(stdout,", %i files total (First index is %03i), ",p->numREMDTRAJ, - p->firstREMDTRAJ); + fprintf(stdout,", %i files total (First index is %0*i), ",p->numREMDTRAJ, + p->EXTwidth, p->firstREMDTRAJ); if (p->isREMDTRAJ>0) fprintf(stdout,"frames at %lf K will be used, ",p->remdtrajtemp); else @@ -3730,8 +3732,9 @@ */ fprintf(stdout, "\nProcessing AMBER REMD trajectory (new format) "); if (currentCoordinateInfo->isREMDTRAJ>0) { - fprintf(stdout, "files %s -> %s%03i\n", + fprintf(stdout, "files %s -> %s%0*i\n", currentCoordinateInfo->filename,currentCoordinateInfo->baseFilename, + currentCoordinateInfo->EXTwidth, currentCoordinateInfo->firstREMDTRAJ+(currentCoordinateInfo->numREMDTRAJ-1)); /* * Now open up the rest of the REMD traj files @@ -3744,9 +3747,11 @@ /*NOTE: Currently spot 0 is not used; inefficient */ for (i=1; inumREMDTRAJ; i++) { if (currentCoordinateInfo->isREMDTRAJ==1) - sprintf(repFilename,"%s%03i",currentCoordinateInfo->baseFilename,j++); + sprintf(repFilename,"%s%0*i",currentCoordinateInfo->baseFilename, + currentCoordinateInfo->EXTwidth,j++); else if (currentCoordinateInfo->isREMDTRAJ==2) - sprintf(repFilename,"%s%03i.gz",currentCoordinateInfo->baseFilename,j++); + sprintf(repFilename,"%s%0*i.gz",currentCoordinateInfo->baseFilename, + currentCoordinateInfo->EXTwidth,j++); fprintf (stdout," Opening %s\n",repFilename); if (openFile(currentCoordinateInfo->remdtrajfiles+i,repFilename,"r")==0) { @@ -4014,12 +4019,10 @@ currentRep=currentCoordinateInfo->remdtrajfiles[i]; /* * Read in REMD header line + * */ - fscanf(currentRep,"%*s"); /*REMD*/ - fscanf(currentRep,"%*s"); /*Replica #*/ - fscanf(currentRep,"%*s"); /*exchange #*/ - fscanf(currentRep,"%*s"); /*step #*/ - fscanf(currentRep,"%lf",&repTemp); /*temp0*/ + fscanf(currentRep,"%*s %*s %*s %*s %lf",&repTemp); + /*fprintf(stdout,"REMD DEBUG: Rep%03i temp= %lf\n",i,repTemp); */ /* @@ -4053,7 +4056,12 @@ safe_free(currentCoordinateInfo->baseFilename); safe_fclose(currentCoordinateInfo->file); currentCoordinateInfo->file = NULL; - fprintf(stdout,"DAN DEBUG: Final repTemp value read= %lf\n",repTemp); + fprintf(stdout,"\nREMDTRAJ: Final repTemp value read= %lf, set %i\n",repTemp,set); + fprintf(stdout,"Could not find target %lf in any of the replica trajectories.\n", + currentCoordinateInfo->remdtrajtemp); + fprintf(stdout,"Check that all replica trajectory files were found and that\n"); + fprintf(stdout, + "none of the trajectories are corrupted (e.g. missing a temperature).\n"); error(ROUTINE, "Target replica temperature not found in traj!"); return; } @@ -4071,7 +4079,7 @@ } *processCoordinates = 1; } - /* cleanup if necessary */ + /* REMDTRAJ: cleanup if necessary */ /* NOTES: Should this go somewhere else? */ if ( *readCoordinates == 0 ) { *processCoordinates = 0; @@ -5746,7 +5754,6 @@ info->isREMDTRAJ=argumentStackContains(&argumentStack,"remdtraj"); info->remdtrajtemp=argumentStackKeyToDouble(&argumentStack, "remdtrajtemp", 0.0); - /*fprintf(stdout," REMD DEBUG: Frames at %lf K will be processed.\n",info->remdtrajtemp);*/ /* DAN ROE: If the coordinate type is COORD_AMBER_REMD (new remd format) * and remdtraj was specified in the trajin, we need to look * for and open any additional traj files. @@ -5755,91 +5762,125 @@ fprintf(stdout," REMDTRAJ: Frames at %lf K will be processed.\n",info->remdtrajtemp); /* * Scan for additional REMD traj files. - * Assume the extension of given trajectory is in the format %03i - * and represents the lowest replica, and that the other files are in - * sequence (i.e. 000, 001, 002 etc). - * - * Get filename 3 char extension - this is assumed to be the first replica + * Assume the extension of given trajectory is the number + * of the lowest replica, and that the other files are in + * sequence (e.g. rem.000, rem.001, rem.002 etc). */ - i=strlen(filename); - EXT=(char*) safe_malloc(4*sizeof(char)); - EXT[0]=filename[i-3]; - EXT[1]=filename[i-2]; - EXT[2]=filename[i-1]; - EXT[3]='\0'; - + loop=strlen(filename); /* Was the file zipped? Currently only the .gz extension is looked for */ - if (strcmp(EXT,".gz")==0) { + if ( (filename[loop-3]=='.')&&(filename[loop-2]=='g')&&(filename[loop-1]=='z') ) { fprintf(stdout," REMDTRAJ: File is gzipped.\n"); info->isREMDTRAJ=2; - EXT[0]=filename[i-6]; - EXT[1]=filename[i-5]; - EXT[2]=filename[i-4]; - EXT[3]='\0'; + loop=loop-3; } - /* Check that extension is a number */ - for (j=0; j<2; j++) { + /* First, find location of last '.' and store it in i*/ + for (j=0; j0) { + fprintf(stdout," REMDDEBUG: Last . in %s located at %i\n",filename,i); + fprintf(stdout," REMDDEBUG: Allocating %i for extension\n",loop-i); + fprintf(stdout," REMDDEBUG: EXTwidth=%i\n",loop-i-1); + } + + /* Get filename extension */ + EXT=(char*) safe_malloc((loop-i)*sizeof(char)); + info->EXTwidth=loop-i-1; + k=0; + for (j=i+1; j0) fprintf(stdout, " REMDDEBUG: Replica extension is %s\n",EXT); + + /* Check that all digits in extension are numbers */ + for (j=0; j57) ) { - fprintf(stdout,"WARNING: REMDTRAJ: Extension is not a number! %s\n",EXT); + fprintf(stdout, + "REMDTRAJ: WARNING: Character #%i (%c) in extension %s is not a number!\n", + j,EXT[j],EXT); fprintf(stdout," Trajectory will be processed as replica traj.\n"); safe_free(EXT); info->isREMDTRAJ=0; } } - /* Look for the other replica files, assuming the name is basefilename.XXX or - * basefilename.XXX.gz + /* Look for the other replica files, assuming the name is basefilename.num or + * basefilename.num.gz */ if (info->isREMDTRAJ>0) { j=atoi(EXT); - if (prnlev>0) fprintf(stderr,"REMD DEBUG: index of first replica = %i\n",j); + if (prnlev>0) fprintf(stderr," REMDDEBUG: index of first replica = %i\n",j); safe_free(EXT); info->firstREMDTRAJ=j; + /* - * allocate memory for the filename + * Allocate memory for the replica filenames. */ - info->baseFilename=(char*) safe_malloc((i+1)*sizeof(char)); - repFilename=(char*) safe_malloc((i+1)*sizeof(char)); + k=strlen(filename); + info->baseFilename=(char*) safe_malloc((k+1)*sizeof(char)); + repFilename=(char*) safe_malloc((k+1)*sizeof(char)); + /* - * store base filename + * Store base filename. + * Variable 'i' still contains location of '.' before # extension + */ + strncpy(info->baseFilename,filename,i+1); + if (prnlev>0) fprintf(stderr," REMDDEBUG: base filename = %s\n",info->baseFilename); + + /* + * Search for a replica number lower than this. Correct functioning + * of the replica code requires the file specified by trajin be the + * lowest # replica. */ if (info->isREMDTRAJ==1) - strncpy(info->baseFilename,filename,i-3); + sprintf(repFilename,"%s%0*i",info->baseFilename,info->EXTwidth,j-1); else if (info->isREMDTRAJ==2) - strncpy(info->baseFilename,filename,i-6); + sprintf(repFilename,"%s%0*i.gz",info->baseFilename,info->EXTwidth,j-1); + if ( (REMDTRAJFILE=fopen(repFilename,"r"))!=NULL ) { + fprintf(stdout, + " REMDTRAJ: WARNING: Replica# found lower than file specified with trajin!\n"); + fprintf(stdout, " (Found %s)\n",repFilename); + fprintf(stdout, " trajin remdtraj requires lowest # replica!\n"); + fclose(REMDTRAJFILE); + } - if (prnlev>0) fprintf(stderr,"REMD DEBUG: base filename = %s\n",info->baseFilename); /* * Now count REMD files */ - i=1; /* file pointer index */ - loop=1; - j++; /* replica file index, start at next replica file */ + i=1; /* # of replica files counted. */ + loop=1; /* Continue loop */ + j++; /* Replica file index, start at next replica file */ fprintf(stderr," REMDTRAJ: Scanning for other REMD files.\n"); while ( loop == 1 ) { /* * filename to scan for */ if (info->isREMDTRAJ==1) - sprintf(repFilename,"%s%03i",info->baseFilename,j); + sprintf(repFilename,"%s%0*i",info->baseFilename,info->EXTwidth,j); else if (info->isREMDTRAJ==2) - sprintf(repFilename,"%s%03i.gz",info->baseFilename,j); + sprintf(repFilename,"%s%0*i.gz",info->baseFilename,info->EXTwidth,j); - if ( openFile(&REMDTRAJFILE,repFilename,"r") == 0 ) { + /* DAN ROE - For some reason popen() called by openFile does not + * return NULL if a gzipped file does not exist. Use + * fopen() here for now since we dont need to do any + * reading yet anyway. + */ + if ( (REMDTRAJFILE=fopen(repFilename,"r"))==NULL ) { loop=0; - if (prnlev>0) fprintf(stderr,"REMD DEBUG: %s not found, exiting loop.\n",repFilename); + if (prnlev>0) fprintf(stderr," REMDDEBUG: %s not found, exiting loop.\n",repFilename); } else { i++; j++; - if (prnlev>0) fprintf(stderr," REMD DEBUG: Found %s\n",repFilename); - safe_fclose(REMDTRAJFILE); + if (prnlev>0) fprintf(stderr," REMDDEBUG: Found %s\n",repFilename); + fclose(REMDTRAJFILE); } - } + } /* End loop to find replica trajs */ fprintf(stderr," REMDTRAJ: Found %i replica traj files.\n",i); info->numREMDTRAJ=i; safe_free(repFilename); } - } + } /* END AMBER REPLICA TRAJECTORY SETUP */ /* * set start, stop and offset if they were specified and relevent diff -Naur trajectory.h trajectory.h --- src/ptraj/trajectory.h 2008-08-29 11:31:24.000000000 -0400 +++ src/ptraj/trajectory.h 2008-08-29 11:31:03.000000000 -0400 @@ -107,11 +107,15 @@ int nlescopy; LesStatus les_status; /* - * REMD info + * REMD Trajectory info */ - int isREMDTRAJ; + int isREMDTRAJ; /* 0 - Process as normal trajectory, + * >0 - Process as replica trajectory, + * 2 - Files have .gz extension + */ int numREMDTRAJ; /* How many replica trajectories are present */ int firstREMDTRAJ; /* Index of first replica */ + int EXTwidth; /* Length of replica extension */ FILE** remdtrajfiles; /* To hold replica trajectory files */ char* baseFilename; /* To hold replica traj base filename */ int linesperset; /* for fast scanthrough of other REMD files */ --------------------------------------------------------------------------- Workarounds: none. ********> bugfix 4 Author: Christine Cezard - q4md-forcefieldtools.org Date: 13 October 2008 Programs: LEaP: xLEaP and tLEaP Usage: Save this file in your $AMBERHOME directory and then apply this patch file to your distribution as follows: cd $AMBERHOME patch -p0 -N < patch.txt recompile LEaP Description: The savemol2 command did not work on some Linux distributions (i.e. Ubuntu, 32 bits), leading to a crash of LEaP with a "stack smashing detected" error when writing the atoms description. The patch below reduces the 'buffer' size when writing the ATOM description. --------------------------------------------------------------------------- --- src/leap/src/leap/pdb_format.c 2007-10-31 21:43:23.000000000 +0100 +++ src/leap/src/leap/pdb_format.c 2008-10-13 18:08:56.000000000 +0200 @@ -143,8 +143,8 @@ "%9 %5d", "MODEL %5d" }, { /* 35 PDB_ENDMDL */ "", "ENDMDL" }, - { /* MOL2_ATOM */ + { /* MOL2_ATOM */ "%3d %6s %12.6f%12.6f%12.6f %3s %3d %3s %12.4f", - "%3d %-4s %11.6f %11.6f %11.6f %-4s %2d %-5s %8.4f ****" }, + "%3d %-4s %11.6f %11.6f %11.6f %-4s %2d %-5s %8.4f ****" }, }; --------------------------------------------------------------------------- Workarounds: none. ********> bugfix 5 Author: Dave Case Date: 11/19/08 Programs: gaff.dat Description: During revisions of gaff.dat in Feb. 2006, some duplicates that had earlier been removed (CVS revisions 8.1 and 8.2) were re-introduced. This fix removes them again. Use the following patch to update $AMBERHOME/dat/leap/parm/gaff.dat --------------------------------------------------------------------------- --- dat/leap/parm/gaff.dat 15 Apr 2008 23:02:36 -0000 10.0 +++ dat/leap/parm/gaff.dat 19 Nov 2008 14:25:18 -0000 @@ -164,7 +164,6 @@ c1-na 452.0 1.362 SOURCE3 8 0.0031 0.0034 c1-ne 741.8 1.220 SOURCE3 1 0.0000 0.0000 c1-nf 741.8 1.220 SOURCE3 1 same as c1-ne -c1-ng 465.7 1.353 SOURCE3 3 0.0025 0.0027 c1-nh 465.7 1.353 SOURCE3 3 0.0025 0.0027 c1-no 393.0 1.405 SOURCE3 3 0.0004 0.0005 c1-o 777.0 1.166 SOURCE2 9 0.0044 0.0052 @@ -334,7 +333,6 @@ cc-pd 318.2 1.733 SOURCE3 84 same as cd-pc cc-ss 279.3 1.737 SOURCE3 52 0.0162 0.0194 cd-cd 418.3 1.429 SOURCE1 740 0.0058 0.0069 -cd-cd 504.0 1.371 SOURCE3 523 same as cc-cd cd-ce 519.2 1.362 SOURCE3 1 same as cc-cf cd-cg 420.9 1.427 SOURCE1 560 0.0000 0.0000 cd-ch 420.9 1.427 SOURCE1 560 0.0000 0.0000 @@ -369,7 +367,6 @@ ce-sx 239.7 1.797 SOURCE3 5 0.0080 0.0082 ce-sy 248.9 1.782 SOURCE3 5 0.0105 0.0114 cf-cf 390.5 1.451 SOURCE1 66 0.0060 0.0060 -cf-cf 562.4 1.338 SOURCE1 543 same as ce-cf cf-cg 415.6 1.431 SOURCE1 22 0.0000 0.0000 cf-ch 415.6 1.431 SOURCE1 22 0.0000 0.0000 cf-ha 341.5 1.089 SOURCE3 55 0.0030 0.0056 @@ -390,7 +387,6 @@ cg-ne 450.5 1.363 SOURCE3 4 0.0189 0.0204 cg-pe 429.8 1.621 SOURCE3 11 0.1928 0.2008 ch-ch 494.2 1.377 SOURCE1 42 0.0000 0.0000 -ch-ch 949.5 1.191 SOURCE1 80 same as cg-ch ch-n1 994.7 1.143 SOURCE1 316 0.0011 0.0018 ch-nf 450.5 1.363 SOURCE3 4 same as cg-ne ch-pf 429.8 1.621 SOURCE3 11 same as cg-pe @@ -539,7 +535,6 @@ n -sh 273.6 1.728 SOURCE3 4 0.0108 0.0128 n -ss 281.6 1.717 SOURCE3 7 0.0120 0.0133 n1-n1 1221.7 1.124 SOURCE1 19 0.0000 0.0000 -n1-n1 1365.7 1.100 SOURCE3 n1-n2 857.4 1.216 SOURCE1 19 0.0000 0.0000 n1-n3 535.7 1.350 SOURCE0 1 n1-n4 518.2 1.360 SOURCE0 1 @@ -651,7 +646,6 @@ nc-ss 433.5 1.560 SOURCE1 74 0.0000 0.0000 nc-sy 439.8 1.555 SOURCE3 2 nd-nd 486.8 1.379 SOURCE3 9 0.0145 0.0164 -nd-nd 602.9 1.315 SOURCE3 9 same as nc-nd nd-os 414.9 1.395 SOURCE1 46 0.0163 0.0188 nd-ss 433.5 1.560 SOURCE1 74 0.0000 0.0000 nd-sy 439.8 1.555 SOURCE3 2 @@ -666,7 +660,6 @@ ne-sx 207.3 1.838 SOURCE3 7 0.0760 0.1060 ne-sy 257.1 1.752 SOURCE3 7 0.0640 0.0814 nf-nf 355.3 1.479 SOURCE3 19 same as ne-ne -nf-nf 738.6 1.257 SOURCE3 1 same as ne-nf nf-o 736.4 1.228 SOURCE3 40 same as ne- o nf-p2 493.9 1.563 SOURCE3 14 same as ne-p2 nf-pf 327.8 1.712 SOURCE3 28 same as ne-pe @@ -787,7 +780,6 @@ pe-sx 145.9 2.168 SOURCE3 9 0.1092 0.1743 pe-sy 133.0 2.213 SOURCE3 6 0.0106 0.0127 pf-pf 240.7 2.092 SOURCE3 7 same as pe-pe -pf-pf 260.8 2.055 SOURCE3 1 same as pe-pf pf-px 291.4 2.005 SOURCE3 12 same as pe-px pf-py 278.6 2.025 SOURCE3 12 same as pe-py pf-s 374.7 1.758 SOURCE3 31 same as pe- s @@ -4014,6 +4006,7 @@ ss-s6-ss 60.5 101.82 SOURCE3 1 0.0000 0.0000 br-sh-hs 43.5 95.64 SOURCE3 1 0.0000 0.0000 c -sh-hs 46.0 96.07 SOURCE3 6 0.9968 1.1164 +c1-sh-hs 29.898 95.99 calculated based on C#C-SH c2-sh-hs 46.2 95.94 SOURCE3 3 0.3373 0.3578 c3-sh-hs 44.9 96.60 SOURCE3 12 0.6315 0.8009 ca-sh-hs 45.7 96.61 SOURCE3 2 0.2867 0.2867 @@ -4236,7 +4229,6 @@ c1-n1-s6 61.8 175.92 HF/6-31G* c1-n1-sh 55.7 174.25 HF/6-31G* c1-n1-ss 55.4 176.06 HF/6-31G* -c1-ng-ng 51.9 179.97 HF/6-31G* c2-c1-n1 63.1 180.00 HF/6-31G* c2-c2-n1 71.7 122.98 HF/6-31G* c2-n1-n1 61.6 180.00 HF/6-31G* @@ -4274,7 +4266,6 @@ n1-n1-os 64.0 176.12 HF/6-31G* n1-n1-p2 69.4 174.71 HF/6-31G* n1-n1-p3 70.2 174.27 HF/6-31G* -n1-n1-s 54.6 179.96 HF/6-31G* n1-n1-sh 57.0 175.07 HF/6-31G* n1-n1-ss 56.9 175.61 HF/6-31G* n1-p4-o 49.7 114.59 HF/6-31G* @@ -4282,7 +4273,6 @@ n1-s2-o 73.1 108.46 HF/6-31G* n1-s4-o 70.2 110.09 HF/6-31G* n1-s6-o 78.8 107.52 HF/6-31G* -n1-c1-n1 70.2 180.00 HF/6-31G* n1-c2-n1 73.6 124.15 HF/6-31G* n1-c3-n1 80.0 105.07 HF/6-31G* n1-ca-n1 75.8 117.03 HF/6-31G* @@ -4325,8 +4315,6 @@ nh-c1-nh 59.0 180.00 Guess no-c1-no 56.8 180.00 Guess na-c1-na 58.6 180.00 Guess -na-c1-na 58.6 180.00 Guess -na-c1-na 58.6 180.00 Guess X -c -c -X 4 1.200 180.000 2.000 X -c -c1-X 2 0.000 180.000 2.000 @@ -4891,7 +4879,6 @@ X -p3-py-X 6 11.000 180.000 2.000 X -p3-sh-X 2 9.200 0.000 2.000 X -p3-ss-X 2 2.300 0.000 3.000 -X -p3-ss-X 2 2.300 0.000 3.000 X -p3-s4-X 4 15.400 0.000 2.000 X -p3-sx-X 4 15.400 0.000 2.000 X -p3-s6-X 6 1.600 0.000 3.000 @@ -5040,7 +5027,6 @@ c -c3-n -o 1.1 180. 2. Junmei et al.1999 c2-c2-na-c3 1.1 180. 2. c2-c -c2-c3 1.1 180. 2. -c3-o -c -oh 1.1 180. 2. c2-c3-c2-hc 1.1 180. 2. Junmei et al.1999 c2-c3-ca-hc 1.1 180. 2. Junmei et al.1999 c2-hc-c -o 1.1 180. 2. Junmei et al.1999 @@ -5063,14 +5049,11 @@ hc-o -c -oh 1.1 180. 2. Junmei et al.1999 hc-o -c -os 1.1 180. 2. n2-c2-ca-n2 1.1 180. 2. dac guess, 9/94 -n2-ca-ca-c3 1.1 180. 2. n2-ca-ca-n2 1.1 180. 2. dac guess, 9/94 -na-ca-ca-c3 1.1 180. 2. na-n2-ca-n2 1.1 180. 2. dac, 10/94 hw ow 0000. 0000. 4. flag for fast water -n na n2 MOD4 RE h1 1.3870 0.0157 Veenstra et al JCC,8,(1992),963 --------------------------------------------------------------------------- Workarounds: none. ********> bugfix 6 Author: Ross Walker Date: 01/14/09 Programs: Ptraj Description: Ptraj can fail to identify some newer CHARMM PSF files such as those labelled as PSF CMAP but without the CHEQ flag. This should correct this problem. Use the following patch to update $AMBERHOME/src/ptraj/ptraj.c --------------------------------------------------------------------------- --- src/ptraj/ptraj.c 2009-01-14 08:50:33.000000000 -0800 +++ src/ptraj/ptraj.c 2009-01-14 09:12:08.000000000 -0800 @@ -3036,11 +3036,18 @@ /* * check to see if the file has "PSF" or "PSF " as the header - * or "PSF CMAP CHEQ" as seen in newer charmm psfs. + * Examples of CHARMM file headers are: + * PSF + * PSF CMAP + * PSF CMAP CHEQ * indicating a CHARMM PSF file, load it and return... */ + /* Look for "PSF\n" or "PSF " - we do this separately since + * It can help with false positives if we just searched for + * PSF for example. + */ if ( strcmp(buffer, "PSF \n") == 0 || - strcmp(buffer, "PSF\n") == 0 || strcmp(buffer, "PSF CMAP CHEQ\n") == 0) { + strncmp(buffer, "PSF ", 4) == 0 ) { state = loadCharmmPSF(fp, 1); *statep = state; return; --------------------------------------------------------------------------- Workarounds: Edit the header of your PSF to read PSF CMAP CHEQ ********>bugfix 7 Author: Daniel R. Roe Date: 02/09/2009 Programs: ptraj Description: Adds ability to read REMD netcdf format files to ptraj. Modifies ptraj.c and trajectory.h Usage: Save this file in your $AMBERHOME directory and then apply this patch file to your amber 10 distribution as follows: cd $AMBERHOME patch -p0 -N -r patch_rejects < bugfix ------------------------------------------------------------------------------ --- src/ptraj/ptraj.c 2008-08-29 11:50:53.000000000 -0400 +++ src/ptraj/ptraj.c 2008-09-03 12:16:26.000000000 -0400 @@ -2391,6 +2439,18 @@ err = nc_get_att_double(NCInfo->ncid, NCInfo->velocityVID, "scale_factor", &NCInfo->velocityScale); } } + + /* + * Are replica temperatures present? + */ + err=nc_inq_varid(NCInfo->ncid,"temp0",&NCInfo->TempVarID); + if (err == NC_NOERR) + fprintf(stdout,"\n\nBintraj has replica temperatures.\n\n"); + else { + fprintf(stdout,"\n\nBintraj DOES NOT have replica temperatures.\n\n"); + NCInfo->TempVarID=-1; + } + #endif break; @@ -3319,6 +3379,7 @@ { coordinateInfo *p; charmmTrajectoryInfo *charmmTraj; + netcdfTrajectoryInfo *NCInfo; byte u; int i; @@ -3382,30 +3443,29 @@ case COORD_AMBER_REMD: fprintf(stdout, " File (%s) is an AMBER REMD (new format) trajectory%s", p->filename, (p->isBox ? " (with box info)" : "")); - fprintf(stdout,", %i files total (First index is %0*i), ",p->numREMDTRAJ, - p->EXTwidth, p->firstREMDTRAJ); - if (p->isREMDTRAJ>0) - fprintf(stdout,"frames at %lf K will be used, ",p->remdtrajtemp); - else - fprintf(stdout,"only frames from this file will be used (remdtraj not specified), "); if (p->stop > 0) { i = (p->stop - p->start)/p->offset+1; if (i != p->stop) - fprintf(stdout, " each file has %i sets (processing only %i)\n", p->stop, i); + fprintf(stdout, " with %i sets (processing only %i)\n", p->stop, i); else - fprintf(stdout, " each file has %i sets\n", i); + fprintf(stdout, " with %i sets\n", i); } else { fprintf(stdout, "\n"); } break; case COORD_AMBER_NETCDF: + + NCInfo=(netcdfTrajectoryInfo*) p->info; fprintf(stdout, " File (%s) is a NetCDF AMBER trajectory%s%s", p->filename, (p->isBox ? " with box info" : ""), (p->isVelocity ? " and velocities" : "")); + if (NCInfo->TempVarID!=-1) + fprintf(stdout," with replica temperatures"); + if (p->stop > 0) { i = (p->stop - p->start)/p->offset+1; if (i != p->stop) @@ -3479,6 +3539,15 @@ break; } + /* REMDTRAJ information */ + if (p->isREMDTRAJ>0) { + fprintf(stdout," Replica processing by temperature will occur.\n"); + fprintf(stdout," %i files total (First index is %0*i), ",p->numREMDTRAJ, + p->EXTwidth, p->firstREMDTRAJ); + fprintf(stdout,"frames at %lf K will be used.\n",p->remdtrajtemp); + } /*else + fprintf(stdout,"only frames from this file will be used (remdtraj not specified), ");*/ + if (prnlev > 4) { if ( p->file == NULL ) { fprintf(stdout, " [the FILE is currently closed]\n"); @@ -3507,6 +3576,7 @@ analyzeInformation *an; coordinateInfo *f; int *mask; + int i; /* * Clean up the INPUT files, transformFileStack @@ -3521,6 +3591,12 @@ safe_free(f->x); safe_free(f->y); safe_free(f->z); + /* REMDTRAJ Cleanup */ + safe_free(f->remdtrajfiles); + safe_free(f->remdncids); + safe_free(f->baseFilename); + safe_free(f->compressEXT); + INITIALIZE_coordinateInfo(f); safe_free(f); @@ -3727,50 +3803,61 @@ break; case COORD_AMBER_REMD: - /* DAN ROE: Modified so that opening of other REMD files only occurs - * if the user specifies remdtraj with trajin. - */ fprintf(stdout, "\nProcessing AMBER REMD trajectory (new format) "); - if (currentCoordinateInfo->isREMDTRAJ>0) { - fprintf(stdout, "files %s -> %s%0*i\n", + if (currentCoordinateInfo->isREMDTRAJ==0) + fprintf(stdout,"file %s",currentCoordinateInfo->filename); + fprintf(stdout,"\n"); + break; + + case COORD_CHARMM_TRAJECTORY: + fprintf(stdout, "\nProcessing CHARMM trajectory file %s\n", + currentCoordinateInfo->filename); + break; + } + + /* REMDTRAJ : If replica temperature processing was requested, + * need to open the other replica files. + */ + if (currentCoordinateInfo->isREMDTRAJ>0) { + fprintf(stdout, " REMDTRAJ: Opening files %s -> %s%0*i\n", currentCoordinateInfo->filename,currentCoordinateInfo->baseFilename, currentCoordinateInfo->EXTwidth, currentCoordinateInfo->firstREMDTRAJ+(currentCoordinateInfo->numREMDTRAJ-1)); - /* - * Now open up the rest of the REMD traj files - */ - fprintf(stdout," Opening other REMD files.\n"); - repFilename=(char*) safe_malloc((strlen(currentCoordinateInfo->filename)+1)*sizeof(char)); + + repFilename=(char*) safe_malloc((strlen(currentCoordinateInfo->filename)+1)*sizeof(char)); + currentCoordinateInfo->remdncids=NULL; + currentCoordinateInfo->remdtrajfiles=NULL; + if (currentCoordinateInfo->type==COORD_AMBER_NETCDF) { + currentCoordinateInfo->remdncids=(int*) safe_malloc( + currentCoordinateInfo->numREMDTRAJ*sizeof(int)); + currentCoordinateInfo->remdncids[0]=NCInfo->ncid; + } else currentCoordinateInfo->remdtrajfiles=(FILE**) safe_malloc( - currentCoordinateInfo->numREMDTRAJ*sizeof(FILE*)); - j=currentCoordinateInfo->firstREMDTRAJ+1; - /*NOTE: Currently spot 0 is not used; inefficient */ - for (i=1; inumREMDTRAJ; i++) { - if (currentCoordinateInfo->isREMDTRAJ==1) - sprintf(repFilename,"%s%0*i",currentCoordinateInfo->baseFilename, - currentCoordinateInfo->EXTwidth,j++); - else if (currentCoordinateInfo->isREMDTRAJ==2) - sprintf(repFilename,"%s%0*i.gz",currentCoordinateInfo->baseFilename, - currentCoordinateInfo->EXTwidth,j++); + currentCoordinateInfo->numREMDTRAJ*sizeof(FILE*)); - fprintf (stdout," Opening %s\n",repFilename); + j=currentCoordinateInfo->firstREMDTRAJ+1; + /*NOTE: Currently spot 0 is not used for not NETCDF; inefficient */ + for (i=1; inumREMDTRAJ; i++) { + sprintf(repFilename,"%s%0*i%s",currentCoordinateInfo->baseFilename, + currentCoordinateInfo->EXTwidth, j++,currentCoordinateInfo->compressEXT); + + fprintf (stdout," Opening %s\n",repFilename); + if (currentCoordinateInfo->type==COORD_AMBER_NETCDF) { + err = nc_open(repFilename, NC_NOWRITE, currentCoordinateInfo->remdncids+i); + if (err != NC_NOERR) error(ROUTINE, "%s", nc_strerror(err)); + /* DEBUG */ + /*dan_netcdf_debug(currentCoordinateInfo->remdncids[i]); */ + } else { if (openFile(currentCoordinateInfo->remdtrajfiles+i,repFilename,"r")==0) { fprintf(stdout,"Error opening REMD traj file %s\n",repFilename); return 1; } } - fprintf(stdout," Done.\n"); - safe_free(repFilename); - } else - fprintf(stdout,"file %s\n",currentCoordinateInfo->filename); - - break; - - case COORD_CHARMM_TRAJECTORY: - fprintf(stdout, "\nProcessing CHARMM trajectory file %s\n", - currentCoordinateInfo->filename); - break; - } + } + fprintf(stdout," Done.\n"); + safe_free(repFilename); + } + /* END REMDTRAJ */ switch ( currentCoordinateInfo->type ) { @@ -3790,16 +3877,16 @@ case COORD_AMBER_REMD: /* - * we need to read each title line to set up for processing - * by readAmberTrajectory() + * First, we need to read the title line from the original file */ if ( fgets(buffer, BUFFER_SIZE, currentCoordinateInfo->file) == NULL ) { fprintf(stdout, "WARNING: Error on processing the title from the AMBER\n"); fprintf(stdout, "REMD trajectory (%s)\n", currentCoordinateInfo->filename); return 1; } - /* DAN ROE: Modified so that the files are read only if user specified - * remdtraj with the trajin. + /* + * Now, if replica temperature trajectory processing was requested, read + * the title line from the other replica files. */ if (currentCoordinateInfo->isREMDTRAJ>0) { for (i=1; inumREMDTRAJ; i++) { @@ -4000,13 +4087,11 @@ if ( set > currentCoordinateInfo->stop ) { *readCoordinates = 0; } else { - /* - * scan through all replica files for target temp - */ + /* + * scan through all replica files for target temp + * if isREMDTRAJ==0 then we just want this replica + */ err=1; - /* - * if isREMDTRAJ==0 then we just want this replica - */ /*fprintf(stdout,"REMD DEBUG: remtrajtemp= %lf\n",currentCoordinateInfo->remdtrajtemp);*/ if (currentCoordinateInfo->isREMDTRAJ==0) j=1; @@ -4052,8 +4137,6 @@ *processCoordinates = 0; for (i=1; inumREMDTRAJ; i++) safe_fclose(currentCoordinateInfo->remdtrajfiles[i]); - safe_free(currentCoordinateInfo->remdtrajfiles); - safe_free(currentCoordinateInfo->baseFilename); safe_fclose(currentCoordinateInfo->file); currentCoordinateInfo->file = NULL; fprintf(stdout,"\nREMDTRAJ: Final repTemp value read= %lf, set %i\n",repTemp,set); @@ -4079,14 +4162,13 @@ } *processCoordinates = 1; } - /* REMDTRAJ: cleanup if necessary */ - /* NOTES: Should this go somewhere else? */ + /* REMDTRAJ: Close files */ if ( *readCoordinates == 0 ) { *processCoordinates = 0; - for (i=1; inumREMDTRAJ; i++) - safe_fclose(currentCoordinateInfo->remdtrajfiles[i]); - safe_free(currentCoordinateInfo->remdtrajfiles); - safe_free(currentCoordinateInfo->baseFilename); + if (currentCoordinateInfo->isREMDTRAJ!=0) { + for (i=1; inumREMDTRAJ; i++) + safe_fclose(currentCoordinateInfo->remdtrajfiles[i]); + } safe_fclose(currentCoordinateInfo->file); currentCoordinateInfo->file = NULL; } @@ -4137,6 +4219,51 @@ if ( NCInfo->R == NULL ) { NCInfo->R = (float *) safe_malloc( sizeof(float) * state->atoms * 3 ); } + + /* DAN REMDTRAJ - Which ncid has our target temperature? + * Note: Temperatures have been checked + * for in checkCoordinates and setupREMDTRAJ. + */ + if (currentCoordinateInfo->isREMDTRAJ>0) { + /*fprintf(stdout,"\nTEMP SEARCH\n"); *DEBUG */ + start[0]=set-1; + count[0]=1; + + j=-1; + for (i=0; inumREMDTRAJ; i++) { + k=currentCoordinateInfo->remdncids[i]; + + err=nc_get_vara_double(k, NCInfo->TempVarID,start,count,&repTemp); + if (err!=NC_NOERR) + fprintf(stderr,"Could not get replica temperature from file# %i\n",i); + /*fprintf(stdout," Replica %i mytemp=%lf ",i,repTemp); DEBUG */ + if (repTemp==currentCoordinateInfo->remdtrajtemp) { + if (prnlev>0) fprintf(stdout, + " REMD target temperature (%lf) found (%lf) at replica index %i\n", + currentCoordinateInfo->remdtrajtemp,repTemp,i); + NCInfo->ncid=k; /* This is the NCID that has our temperature */ + j=0; + break; + } + } + /*fprintf(stdout,"\n"); * DEBUG */ + /* If we didn't find out target temperature, bad news. */ + if (j==-1) { + *readCoordinates = 0; + *processCoordinates = 0; + for (i=0; inumREMDTRAJ; i++) + err=nc_close(currentCoordinateInfo->remdncids[i]); + currentCoordinateInfo->file = NULL; + fprintf(stdout,"\nREMDTRAJ: Final repTemp value read= %lf, set %i\n",repTemp,set); + fprintf(stdout,"Could not find target %lf in any of the replica trajectories.\n", + currentCoordinateInfo->remdtrajtemp); + fprintf(stdout,"Check that all replica trajectory files were found and that\n"); + fprintf(stdout, + "none of the trajectories are corrupted (e.g. missing a temperature).\n"); + error(ROUTINE, "Target replica temperature not found in traj!"); + return; + } + } /* * scan in coordinates for this frame @@ -4226,6 +4353,17 @@ if ( *readCoordinates == 0 ) { *processCoordinates = 0; + + /* REMDTRAJ Cleanup */ + if (currentCoordinateInfo->isREMDTRAJ!=0) { + /* Restore original NCID */ + NCInfo->ncid=currentCoordinateInfo->remdncids[0]; + for (i=1; inumREMDTRAJ; i++) { + err=nc_close(currentCoordinateInfo->remdncids[i]); + if (err) fprintf(stdout,"Error closing replica NetCDF file %i, error: %s\b", + i,nc_strerror(err)); + } + } err = nc_close(NCInfo->ncid); if (err) fprintf(stdout, "Error closing NetCDF file %s (nc_close), error: %s\b", @@ -5252,7 +5390,158 @@ } +/* setupREMDTRAJ(): This is called during a TRAJIN transform in ptrajSetupIO + * to set up this trajectory for processing as an REMD trajectory. All other + * replica trajectory files are searched for, and the base filename, starting + * replica, and number of replicas is recorded so the files may be opened + * later during processing. + */ + + int +setupREMDTRAJ(coordinateInfo *info, char *filename) { + FILE* REMDTRAJFILE; /* FILE* to check for replica trajectories */ + char* EXT; /* Store Replica traj filename extension */ + int loop; /* For looping over replica files */ + char* repFilename; /* Store replica base filename */ + netcdfTrajectoryInfo *NCInfo=NULL; + int i,j,k; + + /* DAN ROE: Check that if remdtraj was specified in trajin this format + * is suitable for REMD processing. + */ + if (info->type==COORD_AMBER_NETCDF) { + NCInfo=(netcdfTrajectoryInfo*) info->info; + if (NCInfo->TempVarID==-1) { + fprintf(stdout,"WARNING: setupREMDTRAJ(): remdtraj requested in trajin but this NETCDF trajectory does not contain temperatures!\n"); + return 1; + } + } else if (info->type!=COORD_AMBER_REMD) { + fprintf(stdout,"WARNING: setupREMDTRAJ(): remdtraj requested in trajin but this is not an amber REMD trajectory!\n"); + return 1; + } + + fprintf(stdout," REMDTRAJ: Frames at %lf K will be processed.\n",info->remdtrajtemp); + /* + * Scan for additional REMD traj files. + * Assume the extension of given trajectory is the number + * of the lowest replica, and that the other files are in + * sequence (e.g. rem.000, rem.001, rem.002 etc). + */ + loop=strlen(filename); + /* Was the file zipped? */ + info->compressEXT=(char*) malloc(5*sizeof(char)); + strcpy(info->compressEXT,""); + if ( (filename[loop-3]=='.')&&(filename[loop-2]=='g')&&(filename[loop-1]=='z') ) { + fprintf(stdout," REMDTRAJ: File is gzipped.\n"); + strcpy(info->compressEXT,".gz"); + loop=loop-3; + } else if ( (filename[loop-4]=='.')&&(filename[loop-3]=='b') + &&(filename[loop-2]=='z')&&(filename[loop-1]=='2') ) { + fprintf(stdout," REMDTRAJ: File is bzipped.\n"); + strcpy(info->compressEXT,".bz2"); + loop=loop-4; + } + + /* First, find location of last '.' and store it in i*/ + for (j=0; j0) { + fprintf(stdout," REMDDEBUG: Last . in %s located at %i\n",filename,i); + fprintf(stdout," REMDDEBUG: Allocating %i for extension\n",loop-i); + fprintf(stdout," REMDDEBUG: EXTwidth=%i\n",loop-i-1); + } + + /* Get filename extension */ + EXT=(char*) safe_malloc((loop-i)*sizeof(char)); + info->EXTwidth=loop-i-1; + k=0; + for (j=i+1; j0) fprintf(stdout, " REMDDEBUG: Replica extension is %s\n",EXT); + + /* Check that all digits in extension are numbers */ + for (j=0; j57) ) { + fprintf(stdout, + "REMDTRAJ: WARNING: Character #%i (%c) in extension %s is not a number!\n", + j,EXT[j],EXT); + safe_free(EXT); + return 1; + } + } + + /* Look for the other replica files, assuming the name is basefilename.num or + * basefilename.num.gz + */ + j=atoi(EXT); + if (prnlev>0) fprintf(stderr," REMDDEBUG: index of first replica = %i\n",j); + safe_free(EXT); + info->firstREMDTRAJ=j; + + /* + * Allocate memory for the replica filenames. + */ + k=strlen(filename); + info->baseFilename=(char*) safe_malloc((k+1)*sizeof(char)); + repFilename=(char*) safe_malloc((k+1)*sizeof(char)); + + /* + * Store base filename. + * Variable 'i' still contains location of '.' before # extension + */ + strncpy(info->baseFilename,filename,i+1); + if (prnlev>0) fprintf(stderr," REMDDEBUG: base filename = %s\n",info->baseFilename); + /* + * Search for a replica number lower than this. Correct functioning + * of the replica code requires the file specified by trajin be the + * lowest # replica. + */ + sprintf(repFilename,"%s%0*i%s",info->baseFilename,info->EXTwidth,j-1,info->compressEXT); + if ( (REMDTRAJFILE=fopen(repFilename,"r"))!=NULL ) { + fprintf(stdout, + " REMDTRAJ: WARNING: Replica# found lower than file specified with trajin!\n"); + fprintf(stdout, " (Found %s)\n",repFilename); + fprintf(stdout, " trajin remdtraj requires lowest # replica!\n"); + fclose(REMDTRAJFILE); + } + + /* + * Now count REMD files + */ + i=1; /* # of replica files counted. */ + loop=1; /* Continue loop */ + j++; /* Replica file index, start at next replica file */ + fprintf(stderr," REMDTRAJ: Scanning for other REMD files.\n"); + while ( loop == 1 ) { + /* + * filename to scan for + */ + sprintf(repFilename,"%s%0*i%s",info->baseFilename,info->EXTwidth,j,info->compressEXT); + + /* DAN ROE - For some reason popen() called by openFile does not + * return NULL if a gzipped file does not exist. Use + * fopen() here for now since we dont need to do any + * reading yet anyway. + */ + if ( (REMDTRAJFILE=fopen(repFilename,"r"))==NULL ) { + loop=0; + if (prnlev>0) fprintf(stderr," REMDDEBUG: %s not found, exiting loop.\n",repFilename); + } else { + i++; j++; + if (prnlev>0) fprintf(stderr," REMDDEBUG: Found %s\n",repFilename); + fclose(REMDTRAJFILE); + } + } /* End loop to find replica trajs */ + fprintf(stderr," REMDTRAJ: Found %i replica traj files.\n",i); + info->numREMDTRAJ=i; + safe_free(repFilename); + + return 0; +} /* * ptrajSetupIO(): This is called upon receipt of the trigger strings @@ -5281,11 +5570,6 @@ stackType *sp; coordinateInfo *infiles; double boxtmp; - /* DAN ROE: Variables for processing new REMD traj format. */ - FILE* REMDTRAJFILE; /* FILE* to check for replica trajectories */ - char* EXT; /* Store Replica traj filename extension */ - int loop; /* For looping over replica files */ - char* repFilename; /* Store replica base filename */ statep = ptrajCurrentState(); state = *statep; @@ -5754,134 +6038,14 @@ info->isREMDTRAJ=argumentStackContains(&argumentStack,"remdtraj"); info->remdtrajtemp=argumentStackKeyToDouble(&argumentStack, "remdtrajtemp", 0.0); - /* DAN ROE: If the coordinate type is COORD_AMBER_REMD (new remd format) - * and remdtraj was specified in the trajin, we need to look - * for and open any additional traj files. - */ - if ((info->type==COORD_AMBER_REMD) && (info->isREMDTRAJ==1)) { - fprintf(stdout," REMDTRAJ: Frames at %lf K will be processed.\n",info->remdtrajtemp); - /* - * Scan for additional REMD traj files. - * Assume the extension of given trajectory is the number - * of the lowest replica, and that the other files are in - * sequence (e.g. rem.000, rem.001, rem.002 etc). - */ - loop=strlen(filename); - /* Was the file zipped? Currently only the .gz extension is looked for */ - if ( (filename[loop-3]=='.')&&(filename[loop-2]=='g')&&(filename[loop-1]=='z') ) { - fprintf(stdout," REMDTRAJ: File is gzipped.\n"); - info->isREMDTRAJ=2; - loop=loop-3; - } - - /* First, find location of last '.' and store it in i*/ - for (j=0; j0) { - fprintf(stdout," REMDDEBUG: Last . in %s located at %i\n",filename,i); - fprintf(stdout," REMDDEBUG: Allocating %i for extension\n",loop-i); - fprintf(stdout," REMDDEBUG: EXTwidth=%i\n",loop-i-1); - } - - /* Get filename extension */ - EXT=(char*) safe_malloc((loop-i)*sizeof(char)); - info->EXTwidth=loop-i-1; - k=0; - for (j=i+1; j0) fprintf(stdout, " REMDDEBUG: Replica extension is %s\n",EXT); - - /* Check that all digits in extension are numbers */ - for (j=0; j57) ) { - fprintf(stdout, - "REMDTRAJ: WARNING: Character #%i (%c) in extension %s is not a number!\n", - j,EXT[j],EXT); - fprintf(stdout," Trajectory will be processed as replica traj.\n"); - safe_free(EXT); - info->isREMDTRAJ=0; - } - } - - /* Look for the other replica files, assuming the name is basefilename.num or - * basefilename.num.gz - */ - if (info->isREMDTRAJ>0) { - j=atoi(EXT); - if (prnlev>0) fprintf(stderr," REMDDEBUG: index of first replica = %i\n",j); - safe_free(EXT); - info->firstREMDTRAJ=j; - - /* - * Allocate memory for the replica filenames. - */ - k=strlen(filename); - info->baseFilename=(char*) safe_malloc((k+1)*sizeof(char)); - repFilename=(char*) safe_malloc((k+1)*sizeof(char)); - - /* - * Store base filename. - * Variable 'i' still contains location of '.' before # extension - */ - strncpy(info->baseFilename,filename,i+1); - if (prnlev>0) fprintf(stderr," REMDDEBUG: base filename = %s\n",info->baseFilename); - - /* - * Search for a replica number lower than this. Correct functioning - * of the replica code requires the file specified by trajin be the - * lowest # replica. - */ - if (info->isREMDTRAJ==1) - sprintf(repFilename,"%s%0*i",info->baseFilename,info->EXTwidth,j-1); - else if (info->isREMDTRAJ==2) - sprintf(repFilename,"%s%0*i.gz",info->baseFilename,info->EXTwidth,j-1); - if ( (REMDTRAJFILE=fopen(repFilename,"r"))!=NULL ) { - fprintf(stdout, - " REMDTRAJ: WARNING: Replica# found lower than file specified with trajin!\n"); - fprintf(stdout, " (Found %s)\n",repFilename); - fprintf(stdout, " trajin remdtraj requires lowest # replica!\n"); - fclose(REMDTRAJFILE); - } - - /* - * Now count REMD files - */ - i=1; /* # of replica files counted. */ - loop=1; /* Continue loop */ - j++; /* Replica file index, start at next replica file */ - fprintf(stderr," REMDTRAJ: Scanning for other REMD files.\n"); - while ( loop == 1 ) { - /* - * filename to scan for - */ - if (info->isREMDTRAJ==1) - sprintf(repFilename,"%s%0*i",info->baseFilename,info->EXTwidth,j); - else if (info->isREMDTRAJ==2) - sprintf(repFilename,"%s%0*i.gz",info->baseFilename,info->EXTwidth,j); - - /* DAN ROE - For some reason popen() called by openFile does not - * return NULL if a gzipped file does not exist. Use - * fopen() here for now since we dont need to do any - * reading yet anyway. - */ - if ( (REMDTRAJFILE=fopen(repFilename,"r"))==NULL ) { - loop=0; - if (prnlev>0) fprintf(stderr," REMDDEBUG: %s not found, exiting loop.\n",repFilename); - } else { - i++; j++; - if (prnlev>0) fprintf(stderr," REMDDEBUG: Found %s\n",repFilename); - fclose(REMDTRAJFILE); - } - } /* End loop to find replica trajs */ - fprintf(stderr," REMDTRAJ: Found %i replica traj files.\n",i); - info->numREMDTRAJ=i; - safe_free(repFilename); + /* DAN ROE: Moved REMD setup routine to its own place */ + if (info->isREMDTRAJ==1) { + if ( setupREMDTRAJ(info,filename)!=0 ) { + info->isREMDTRAJ==0; + fprintf(stdout," Trajectory will be processed as replica traj.\n"); } - } /* END AMBER REPLICA TRAJECTORY SETUP */ - + } + /* * set start, stop and offset if they were specified and relevent */ --- src/ptraj/trajectory.h 2008-08-29 11:50:53.000000000 -0400 +++ src/ptraj/trajectory.h 2008-09-03 10:59:25.000000000 -0400 @@ -109,14 +109,15 @@ /* * REMD Trajectory info */ - int isREMDTRAJ; /* 0 - Process as normal trajectory, - * >0 - Process as replica trajectory, - * 2 - Files have .gz extension + int isREMDTRAJ; /* 0 - Process as normal trajectory, + * 1 - Process as replica trajectory */ + char *compressEXT; /* If not NULL, contains compressed traj ext.*/ int numREMDTRAJ; /* How many replica trajectories are present */ int firstREMDTRAJ; /* Index of first replica */ int EXTwidth; /* Length of replica extension */ FILE** remdtrajfiles; /* To hold replica trajectory files */ + int* remdncids; /* To hold netcdf replica trajectory files */ char* baseFilename; /* To hold replica traj base filename */ int linesperset; /* for fast scanthrough of other REMD files */ double remdtrajtemp; /* Target temperature */ @@ -210,6 +211,7 @@ int cellLengthVID; int cellAngleVID; int velocityVID; + int TempVarID; double velocityScale; char *timeUnits; char *coordinateUnits; @@ -239,6 +241,7 @@ _p_->cellLengthVID = 0; \ _p_->cellAngleVID = 0; \ _p_->velocityVID = 0; \ + _p_->TempVarID = 0; \ _p_->timeUnits = NULL; \ _p_->coordinateUnits = NULL; \ _p_->cellLengthUnits = NULL; \ ********> bugfix.8 Author: Daniel R. Roe Date: 02/10/2009 Programs: ptraj Description: Ptraj would segfault when attempting to print info on output trajectory if output was in netcdf format. Modifies ptraj.c Usage: Save this file in your $AMBERHOME directory and then apply this patch file to your amber 10 distribution as follows: cd $AMBERHOME patch -p0 -N -r patch_rejects < ptraj.remdbintraj.fix.patch ------------------------------------------------------------------------------ --- src/ptraj/ptraj.c 2009-02-10 14:57:36.000000000 -0500 +++ src/ptraj/ptraj.c 2009-02-10 14:55:28.000000000 -0500 @@ -3415,15 +3415,16 @@ case COORD_AMBER_NETCDF: - NCInfo=(netcdfTrajectoryInfo*) p->info; - fprintf(stdout, " File (%s) is a NetCDF AMBER trajectory%s%s", p->filename, (p->isBox ? " with box info" : ""), (p->isVelocity ? " and velocities" : "")); - if (NCInfo->TempVarID!=-1) - fprintf(stdout," with replica temperatures"); + if (p->info!=NULL) { + NCInfo=(netcdfTrajectoryInfo*) p->info; + if (NCInfo->TempVarID!=-1) + fprintf(stdout," with replica temperatures"); + } if (p->stop > 0) { i = (p->stop - p->start)/p->offset+1; ********> bugfix.9 Author: Dave Case Date: 04/23/2009 Programs: nab Description: The nab dihedral codes are *not* equipped to handle arbitrary phase angles (as in parmbsc0). Check for that. Note that this effectively undoes bugfix.2 Usage: Use this patch to update amber10/src/nab/prm.c ------------------------------------------------------------------------------ *** src/nab/prm.c 2008-04-08 17:10:11.000000000 -0400 --- src/nab/prm.c 2009-04-23 10:54:51.000000000 -0400 *************** *** 27,32 **** --- 27,33 ---- #include #include #include + #include #include "nab.h" #include "errormsg.h" *************** *** 845,850 **** --- 846,854 ---- #else fscanf(file, " %f", &prm->Phase[i]); #endif + /* the NAB dihedral codes are not (yet) set up for arbitrary phases */ + assert( prm->Phase[i] == 0.0 || fabs(3.14159 - prm->Phase[i]) < 0.01 ); + } } ********> bugfix.10 Author: Dave Case Date: 07/17/2009 Programs: tleap Description: fix up the residue numbers in loadMol2 if the input mol2 file has residue numbers outside the range 1 to number of substructures phase angles (as in parmbsc0). Check for that. Note that this effectively undoes bugfix.2 Usage: Use this patch to update amber10/src/leap/src/leap/tripos.c ------------------------------------------------------------------------------ --- src/leap/src/leap/tripos.c 15 Apr 2008 23:22:21 -0000 10.0 +++ src/leap/src/leap/tripos.c 7 Jul 2009 20:04:46 -0000 @@ -182,6 +182,11 @@ if (iRet != 9) break; + /* don't allow the residue number to be greater than no. of + substructures, nor less than 1: */ + iResidue = iResidue > iFileSubstructures ? iFileSubstructures : iResidue; + iResidue = iResidue < 1 ? 1 : iResidue; + MESSAGE((" Atom: %s\n", sName)); aAtom = (ATOM) oCreate(ATOMid); ContainerSetName(aAtom, sName); ------------------------------------------------------------------------------ Workarounds: none ********> bugfix.11 Author: Junmei Wang Date: 08/21/2009 Programs: antechamber Description: fix the assignement of am1-bcc bond types in SO2 and NO2 groups Usage: Use this patch to update amber10/src/antechamber/bondtype.c: ------------------------------------------------------------------------------ --- src/antechamber/bondtype.c.orig 2008-05-22 00:22:01.000000000 -0400 +++ src/antechamber/bondtype.c 2009-08-21 08:45:57.000000000 -0400 @@ -314,7 +314,9 @@ conjatom[i] = 0; if (index1 >= 2 || index2 >= 2) conjatom[i] = 1; - if (index3 >= 2 || index4 >= 3) + /* index3 >= 3 implies that S=O bonds in SO2 will not be delocalized; + index6 >= 2 makes NO2 bonds delocalized */ + if (index3 >= 3 || index6 >= 2) conjatom[i] = 1; for (;;) { if (fgets(line, MAXCHAR, fp) == NULL) ------------------------------------------------------------------------------ Workarounds: none ********> bugfix.12 Author: Dave Case Date: 08/21/2009 Programs: LEaP, force fields Description: update the ff99SB and ff99bsc0 leaprc files to handle PDB format 3 atom and residue names Usage: Use this patch to update amber10/dat/leap/cmd/leaprc.ff99{SB,bsc0}: ------------------------------------------------------------------------------ --- dat/leap/cmd/leaprc.ff99SB 15 Apr 2008 23:02:36 -0000 10.0 +++ dat/leap/cmd/leaprc.ff99SB 17 Jun 2009 11:55:47 -0000 10.2 @@ -1,6 +1,7 @@ logFile leap.log # # ----- leaprc for loading the ff99SB (Hornak & Simmerling) force field +# ----- this file is updated for PDB format 3 # # load atom type hybridizations # @@ -38,7 +39,6 @@ { "CD" "C" "sp2" } { "CE" "C" "sp2" } { "CF" "C" "sp2" } - { "CG" "C" "sp2" } { "CP" "C" "sp2" } { "CI" "C" "sp2" } { "CJ" "C" "sp2" } @@ -131,18 +131,14 @@ { 0 "ADE" "DA5" } { 1 "ADE" "DA3" } { "ADE" "DA" } { 0 "CYT" "DC5" } { 1 "CYT" "DC3" } { "CYT" "DC" } { 0 "THY" "DT5" } { 1 "THY" "DT3" } { "THY" "DT" } - { 0 "G" "DG5" } { 1 "G" "DG3" } { "G" "DG" } { "GN" "DGN" } - { 0 "A" "DA5" } { 1 "A" "DA3" } { "A" "DA" } { "AN" "DAN" } - { 0 "C" "DC5" } { 1 "C" "DC3" } { "C" "DC" } { "CN" "DCN" } - { 0 "T" "DT5" } { 1 "T" "DT3" } { "T" "DT" } { "TN" "DTN" } - { 0 "C5" "DC5" } - { 0 "G5" "DG5" } - { 0 "A5" "DA5" } - { 0 "T5" "DT5" } - { 1 "C3" "DC3" } - { 1 "G3" "DG3" } - { 1 "A3" "DA3" } - { 1 "T3" "DT3" } + { 0 "G" "RG5" } { 1 "G" "RG3" } { "G" "RG" } { "GN" "RGN" } + { 0 "A" "RA5" } { 1 "A" "RA3" } { "A" "RA" } { "AN" "RAN" } + { 0 "C" "RC5" } { 1 "C" "RC3" } { "C" "RC" } { "CN" "RCN" } + { 0 "U" "RU5" } { 1 "U" "RU3" } { "U" "RU" } { "UN" "RUN" } + { 0 "DG" "DG5" } { 1 "DG" "DG3" } + { 0 "DA" "DA5" } { 1 "DA" "DA3" } + { 0 "DC" "DC5" } { 1 "DC" "DC3" } + { 0 "DT" "DT5" } { 1 "DT" "DT3" } } @@ -156,17 +152,26 @@ { "C2*" "C2'" } { "C1*" "C1'" } { "C5M" "C7" } + { "O2*" "O2'" } { "H1*" "H1'" } { "H2*1" "H2'1" } { "H2*2" "H2'2" } + { "H2'" "H2'1" } + { "H2''" "H2'2" } { "H3*" "H3'" } { "H4*" "H4'" } { "H5*1" "H5'1" } { "H5*2" "H5'2" } -# old ff atom names -> new + { "H5'" "H5'1" } + { "H5''" "H5'2" } + { "HO2'" "HO'2" } + { "HO5'" "H5T" } + { "HO3'" "H3T" } { "O1'" "O4'" } { "OA" "O1P" } { "OB" "O2P" } + { "OP1" "O1P" } + { "OP2" "O2P" } } --- dat/leap/cmd/leaprc.ff99bsc0 30 May 2008 16:10:02 -0000 10.1 +++ dat/leap/cmd/leaprc.ff99bsc0 17 Jun 2009 11:55:47 -0000 10.3 @@ -1,7 +1,8 @@ logFile leap.log # # ----- leaprc for loading the ff99SB (Hornak & Simmerling) force field -# using parmbsc0 for nucleic acids (assumed to be DNA) +# using parmbsc0 for nucleic acid. +# ----- NOTE: this is designed for PDB format 3! # # load atom type hybridizations # @@ -39,7 +40,6 @@ { "CD" "C" "sp2" } { "CE" "C" "sp2" } { "CF" "C" "sp2" } - { "CG" "C" "sp2" } { "CP" "C" "sp2" } { "CI" "C" "sp2" } { "CJ" "C" "sp2" } @@ -134,18 +134,14 @@ { 0 "ADE" "DA5" } { 1 "ADE" "DA3" } { "ADE" "DA" } { 0 "CYT" "DC5" } { 1 "CYT" "DC3" } { "CYT" "DC" } { 0 "THY" "DT5" } { 1 "THY" "DT3" } { "THY" "DT" } - { 0 "G" "DG5" } { 1 "G" "DG3" } { "G" "DG" } { "GN" "DGN" } - { 0 "A" "DA5" } { 1 "A" "DA3" } { "A" "DA" } { "AN" "DAN" } - { 0 "C" "DC5" } { 1 "C" "DC3" } { "C" "DC" } { "CN" "DCN" } - { 0 "T" "DT5" } { 1 "T" "DT3" } { "T" "DT" } { "TN" "DTN" } - { 0 "DC" "DC5" } - { 0 "DG" "DG5" } - { 0 "DA" "DA5" } - { 0 "DT" "DT5" } - { 1 "DC" "DC3" } - { 1 "DG" "DG3" } - { 1 "DA" "DA3" } - { 1 "DT" "DT3" } + { 0 "G" "RG5" } { 1 "G" "RG3" } { "G" "RG" } { "GN" "RGN" } + { 0 "A" "RA5" } { 1 "A" "RA3" } { "A" "RA" } { "AN" "RAN" } + { 0 "C" "RC5" } { 1 "C" "RC3" } { "C" "RC" } { "CN" "RCN" } + { 0 "U" "RU5" } { 1 "U" "RU3" } { "U" "RU" } { "UN" "RUN" } + { 0 "DG" "DG5" } { 1 "DG" "DG3" } + { 0 "DA" "DA5" } { 1 "DA" "DA3" } + { 0 "DC" "DC5" } { 1 "DC" "DC3" } + { 0 "DT" "DT5" } { 1 "DT" "DT3" } } @@ -157,23 +153,29 @@ { "C3*" "C3'" } { "O3*" "O3'" } { "C2*" "C2'" } + { "O2*" "O2'" } { "C1*" "C1'" } { "C5M" "C7" } { "H1*" "H1'" } { "H2*1" "H2'1" } { "H2*2" "H2'2" } - { "H2'" "H2'1" } + { "H2'" "H2'1" } { "H2''" "H2'2" } { "H3*" "H3'" } { "H4*" "H4'" } { "H5*1" "H5'1" } { "H5*2" "H5'2" } - { "H5'" "H5'1" } + { "H5'" "H5'1" } { "H5''" "H5'2" } - { "OP1" "O1P" } - { "OP2" "O2P" } + { "HO2'" "HO'2" } { "HO5'" "H5T" } { "HO3'" "H3T" } +# old ff atom names -> new + { "O1'" "O4'" } + { "OA" "O1P" } + { "OB" "O2P" } + { "OP1" "O1P" } + { "OP2" "O2P" } } ------------------------------------------------------------------------------ Workarounds: none