24 #ifdef MEM_OPT_VERSION
26 PDB::PDB(
const char *pdbfilename,
int expectedNumAtoms){
31 atomArray =
new PDBCoreData[expectedNumAtoms];
34 atomListHead = atomListTail = NULL;
35 infile =
Fopen(pdbfilename,
"r");
38 sprintf(s,
"Cannot open file '%s' for input in PDB::PDB.", pdbfilename);
43 while ( fgets(buf, 150, infile) ) {
46 for (s=buf; *s && *s!=
'\n'; s++)
49 if ( s == (buf + 149) ) {
51 sprintf( s,
"Line too long in pdbfile %s:\n%s\n", pdbfilename, buf);
59 NAMD_die(
"Could not allocate PDBData.\n");
69 atomArray[atomCount].coor[0] = coor[0];
70 atomArray[atomCount].coor[1] = coor[1];
71 atomArray[atomCount].coor[2] = coor[2];
72 atomArray[atomCount].myoccupancy = thisAtom->
occupancy();
76 if(atomCount > expectedNumAtoms)
NAMD_die(
"Number of pdb and psf atoms are not the same!");
83 PDB::PDB(molfile_plugin_t *pIOHdl,
void *pIOFileHdl,
int numAtoms,
const float *occupancy,
const float *bfactor){
84 #ifdef MEM_OPT_VERSION
85 NAMD_die(
"Sorry, plugin IO is not supported in the memory optimized version.");
87 molfile_timestep_t ts;
89 memset(&ts, 0,
sizeof(molfile_timestep_t));
92 ts.A = ts.B = ts.C = 0.0f;
93 ts.alpha = ts.beta = ts.gamma = 90.0f;
97 atomcoords = (
float *) malloc(3*numAtoms*
sizeof(
float));
98 memset(atomcoords, 0, 3*numAtoms*
sizeof(
float));
99 ts.coords = atomcoords;
101 if (pIOHdl->read_next_timestep(pIOFileHdl, numAtoms, &ts)) {
103 pIOHdl->close_file_read(pIOFileHdl);
104 NAMD_die(
"ERROR: failed reading atom coordinates");
112 for(
int i=0; i<numAtoms; i++) {
113 PDBAtom *thisAtom = tmpAtoms+i;
114 atomArray[i] = thisAtom;
116 tmpCoords[0] = atomcoords[0];
117 tmpCoords[1] = atomcoords[1];
118 tmpCoords[2] = atomcoords[2];
125 for(
int i=0; i<numAtoms; i++) {
126 PDBAtom *thisAtom = tmpAtoms+i;
131 for(
int i=0; i<numAtoms; i++) {
132 PDBAtom *thisAtom = tmpAtoms+i;
146 atomListHead = atomListTail = NULL;
147 infile =
Fopen(pdbfilename,
"r");
150 sprintf(s,
"Cannot open file '%s' for input in PDB::PDB.", pdbfilename);
155 while ( fgets(buf, 150, infile) ) {
158 for (s=buf; *s && *s!=
'\n'; s++)
161 if ( s == (buf + 149) ) {
163 sprintf( s,
"Line too long in pdbfile %s:\n%s\n", pdbfilename, buf);
171 NAMD_die(
"Could not allocate PDBData.\n");
179 add_atom_element( (
PDBAtom *) newelement);
189 #ifdef MEM_OPT_VERSION
191 altlocArray =
new char[atomCount];
192 atomArray =
new PDBCoreData[atomCount];
193 if(atomArray == NULL)
194 NAMD_die(
"memory allocation failed in PDB::PDB");
197 for(
int i=0; tmp!=NULL; tmp=tmp->
next, i++){
199 atomArray[i].coor[0] = coor[0];
200 atomArray[i].coor[1] = coor[1];
201 atomArray[i].coor[2] = coor[2];
215 atomListHead = atomListTail = NULL;
219 if ( atomArray == NULL )
221 NAMD_die(
"memory allocation failed in PDB::PDB");
225 for (i=0, tmp = atomListHead; tmp != NULL; tmp = tmp ->
next, i++) {
226 atomArray[i] = tmp -> data;
230 for (tmp2 = tmp = atomListHead; tmp != NULL; tmp = tmp2) {
234 atomListHead = atomListTail = NULL;
244 #ifndef MEM_OPT_VERSION
249 for (i=atomCount-1; i>=0; i--)
253 delete [] altlocArray;
261 void PDB::write(
const char *outfilename,
const char *commentline)
266 if ((outfile = fopen(outfilename,
"w")) == NULL) {
267 sprintf(s,
"Cannot open file '%s' in PDB::write.", outfilename);
271 if (commentline != NULL)
273 sprintf(s,
"REMARK %s\n", commentline);
274 if (fputs(s, outfile) == EOF)
276 NAMD_err(
"EOF in PDB::write writing the comment line - file system full?");
280 for (i=0; i<atomCount; i++){
281 #ifdef MEM_OPT_VERSION
286 if ( (fputs(s, outfile) == EOF) ||
287 (fputs(
"\n", outfile) == EOF) ) {
288 sprintf(s,
"EOF in PDB::write line %d - file system full?", i);
292 if (fputs(
"END\n", outfile) == EOF) {
293 NAMD_err(
"EOF in PDB::write while printing 'END' -- file system full?");
295 if (fclose(outfile) == EOF) {
296 NAMD_err(
"EOF in PDB::write while closing -- file system full?");
302 void PDB::add_atom_element(
PDBAtom *newAtom)
307 NAMD_die(
"memory allocation failed in PDB::add_atom_element");
309 tmp -> data = newAtom;
312 if (atomListHead == NULL) {
313 atomListHead = atomListTail = tmp;
315 atomListTail ->
next = tmp;
333 #ifdef MEM_OPT_VERSION
334 for(
int i=0; i<atomCount; i++){
335 atomArray[i].coor[0] = pos[i].
x;
336 atomArray[i].coor[1] = pos[i].
y;
337 atomArray[i].coor[2] = pos[i].
z;
343 for (i=0, atomptr=atomArray; i<atomCount; atomptr++, i++)
345 (*atomptr)->
xcoor(pos[i].
x);
346 (*atomptr)->ycoor(pos[i].
y);
347 (*atomptr)->zcoor(pos[i].
z);
353 #ifdef MEM_OPT_VERSION
354 pos->
x = atomArray[aid].coor[0];
355 pos->
y = atomArray[aid].coor[1];
356 pos->
z = atomArray[aid].coor[2];
359 pos->
x = (*atomptr)->xcoor();
360 pos->
y = (*atomptr)->ycoor();
361 pos->
z = (*atomptr)->zcoor();
368 #ifdef MEM_OPT_VERSION
369 for(
int i=0; i<atomCount; i++){
370 pos[i].
x = atomArray[i].coor[0];
371 pos[i].
y = atomArray[i].coor[1];
372 pos[i].
z = atomArray[i].coor[2];
378 for (i=0, atomptr=atomArray; i<atomCount; atomptr++, i++)
380 pos[i].
x = (*atomptr)->xcoor();
381 pos[i].
y = (*atomptr)->ycoor();
382 pos[i].
z = (*atomptr)->zcoor();
388 #ifdef MEM_OPT_VERSION
390 if(place<0 || place>=atomCount)
return NULL;
391 return &atomArray[place];
396 if (place <0 || place >= atomCount)
398 return atomArray[place];
405 void PDB::find_extremes_helper(
411 #ifdef MEM_OPT_VERSION
412 PDBCoreData *atomptr = atomArray;
413 for(
int i=0; i<atomCount; i++, atomptr++){
414 Vector pos(atomptr->xcoor(),atomptr->ycoor(),atomptr->zcoor());
419 for (
int i=0; i<atomCount; ++i, ++atomptr) {
426 int ilow = (int)((1.0 - frac) * atomCount);
427 if ( ilow < 0 ) ilow = 0;
428 if ( ilow > atomCount/2 ) ilow = atomCount/2;
429 int ihigh = atomCount - ilow - 1;
430 BigReal span = coor[ihigh] - coor[ilow];
431 BigReal extension = (1.0 - frac) * span / (2.0 * frac - 1.0);
432 max = coor[ihigh] + extension;
433 min = coor[ilow] - extension;
439 if (atomCount == 0) {
442 else if (frac < 1.0) {
447 find_extremes_helper(coor, smin.
x, smax.
x, lattice.
a_r(), frac);
448 find_extremes_helper(coor, smin.
y, smax.
y, lattice.
b_r(), frac);
449 find_extremes_helper(coor, smin.
z, smax.
z, lattice.
c_r(), frac);
453 #ifdef MEM_OPT_VERSION
454 PDBCoreData *atomptr = atomArray;
455 Vector p(atomptr->xcoor(),atomptr->ycoor(),atomptr->zcoor());
464 for(
int i=1; i<atomCount; i++, atomptr++){
465 #ifdef MEM_OPT_VERSION
466 p =
Vector(atomptr->xcoor(),atomptr->ycoor(),atomptr->zcoor());
469 p =
Vector (atom->xcoor(),atom->ycoor(),atom->zcoor());
472 if (smin.x > s.
x) smin.
x = s.
x;
473 else if (smax.
x < s.
x) smax.
x = s.
x;
474 if (smin.y > s.
y) smin.y = s.
y;
475 else if (smax.
y < s.
y) smax.
y = s.
y;
476 if (smin.z > s.
z) smin.z = s.
z;
477 else if (smax.
z < s.
z) smax.
z = s.
z;
482 origin_shift = lattice.
a_r() * lattice.
origin();
483 smin.
x -= origin_shift;
484 smax.
x -= origin_shift;
485 origin_shift = lattice.
b_r() * lattice.
origin();
486 smin.
y -= origin_shift;
487 smax.
y -= origin_shift;
488 origin_shift = lattice.
c_r() * lattice.
origin();
489 smin.
z -= origin_shift;
490 smax.
z -= origin_shift;
494 #ifdef TEST_PDB_CLASS
498 PDB *pdb =
new PDB(
"pti.pdb");
499 if ( atomArray == NULL )
501 NAMD_die(
"memory allocation failed in main of test PDB class");
503 ilist = pdb->find_atom_name(
"CA");
505 printf(
"None found.\n");
509 for (i=0; i<ilist->num(); i++) {
510 pdb->
atom((*ilist)[i]) -> sprint(s);
516 printf(
"Now, search through space.\n");
518 ilist = pdb->find_atoms_in_region(4.38, 19.5, 3.0,
521 printf(
"None found.\n");
525 printf(
"%d found\n", ilist -> num());
526 for (i=0; i<ilist->num(); i++) {
527 pdb->
atom((*ilist)[i]) -> sprint(s);
533 #endif // TEST_PDB_CLASS
543 while((c = getc(f)) !=
'\n') {
557 char buf[13],resname[5],atomname[5];
561 if ((infile=
Fopen(filename,
"r")) == NULL)
562 NAMD_err(
"Can't open AMBER coordinate file!");
566 fscanf(infile,
"%d",&atomCount);
567 if (atomCount != amber_data->
Natom)
568 NAMD_die(
"Num of atoms in coordinate file is different from that in parm file!");
571 #ifdef MEM_OPT_VERSION
573 atomArray =
new PDBCoreData[atomCount];
579 if ( atomArray == NULL )
581 NAMD_die(
"memory allocation failed in PDB::PDB");
586 for (i=0; i<atomCount; ++i)
589 {
for (k=0; k<12; ++k)
590 { buf[k]=getc(infile);
591 if (buf[k]==
'\n' || buf[k]==
'\0' || buf[k]==EOF)
592 NAMD_die(
"Error reading AMBER coordinate file!");
599 #ifdef MEM_OPT_VERSION
600 atomArray[i].coor[0] = coor[0];
601 atomArray[i].coor[1] = coor[1];
602 atomArray[i].coor[2] = coor[2];
609 atomname[j] = amber_data->
AtomNames[i*4+j];
611 resname[4] = atomname[4] =
'\0';
628 char buf[13],resname[5],atomname[5];
632 if ((infile=
Fopen(filename,
"r")) == NULL)
633 NAMD_err(
"Can't open AMBER coordinate file!");
637 fscanf(infile,
"%d",&atomCount);
638 if (atomCount != amber_data->
Natom)
639 NAMD_die(
"Num of atoms in coordinate file is different from that in parm file!");
642 #ifdef MEM_OPT_VERSION
644 atomArray =
new PDBCoreData[atomCount];
650 if ( atomArray == NULL )
652 NAMD_die(
"memory allocation failed in PDB::PDB");
657 for (i=0; i<atomCount; ++i)
660 {
for (k=0; k<12; ++k)
661 { buf[k]=getc(infile);
662 if (buf[k]==
'\n' || buf[k]==
'\0' || buf[k]==EOF)
663 NAMD_die(
"Error reading AMBER coordinate file!");
670 #ifdef MEM_OPT_VERSION
671 atomArray[i].coor[0] = coor[0];
672 atomArray[i].coor[1] = coor[1];
673 atomArray[i].coor[2] = coor[2];
678 strncpy(resname, amber_data->
ResNames[amber_data->
AtomRes[i]].c_str(), 5);
679 strncpy(atomname, amber_data->
AtomNames[i].c_str(), 5);
683 pdb->residuename(resname);
684 pdb->serialnumber(i+1);
685 pdb->residueseq(amber_data->
AtomRes[i]+1);
686 pdb->coordinates(coor);
704 infile=
Fopen(filename,
"r");
706 NAMD_err(
"Can't open GROMACS coordinate file!");
713 sscanf(buf,
"%d",&atomCount);
715 NAMD_die(
"Num of atoms in coordinate file is different from that in topology file!");
718 #ifdef MEM_OPT_VERSION
720 atomArray =
new PDBCoreData[atomCount];
725 if ( atomArray == NULL )
726 NAMD_die(
"memory allocation failed in PDB::PDB");
728 #ifdef MEM_OPT_VERSION
729 for(i=0; i<atomCount; i++){
733 if(3 != sscanf(buf2,
"%lf%lf%lf", &coor[0],&coor[1],&coor[2]))
734 NAMD_die(
"Couldn't get three coordinates from file.");
740 atomArray[i].coor[0] = coor[0];
741 atomArray[i].coor[1] = coor[1];
742 atomArray[i].coor[2] = coor[2];
747 for (i=0;i<atomCount;i++) {
748 char *buf2, resname[11], atomname[11], atmtype[11];
756 if(3 != sscanf(buf2,
"%lf%lf%lf",
757 &coor[0],&coor[1],&coor[2]))
758 NAMD_die(
"Couldn't get three coordinates from file.");
759 topology->
getAtom(i,&resnum,resname,
760 atomname,atmtype,&typenum,&charge,&mass);
776 #ifdef MEM_OPT_VERSION
780 for(
int i=0; i<79; i++) outstr[i] = 32;
791 sprintf(tmpstr,
" #");
794 strcat(outstr, tmpstr);
796 sprintf(tmpstr,
" #");
799 strcat(outstr, tmpstr);
801 sprintf(tmpstr,
" #");
804 strcat(outstr, tmpstr);
807 strcat(outstr, tmpstr);
810 strcat(outstr, tmpstr);
static const BigReal default_coor
void NAMD_err(const char *err_msg)
static const BigReal default_temperaturefactor
void get_position_for_atom(Vector *, int)
vector< string > ResNames
BigReal temperaturefactor(void)
void write(const char *outfilename, const char *commentline=NULL)
void set_all_positions(Vector *)
static int readtoeoln(FILE *f)
const char * residuename(void)
static Units next(Units u)
static void sprintcol(char *s, int start, int len, const char *val)
PDBAtom * atom(int place)
const BigReal * coordinates(void)
FILE * Fopen(const char *filename, const char *mode)
const char * alternatelocation(void)
void NAMD_die(const char *err_msg)
PDB(const char *pdbfilename)
vector< string > AtomNames
void find_extremes(const Lattice &, BigReal frac=1.0)
int main(int argc, char *argv[])
k< npairi;++k){TABENERGY(const int numtypes=simParams->tableNumTypes;const float table_spacing=simParams->tableSpacing;const int npertype=(int)(namdnearbyint(simParams->tableMaxDist/simParams->tableSpacing)+1);) int table_i=(r2iilist[2 *k] >> 14)+r2_delta_expc;const int j=pairlisti[k];#define p_j BigReal diffa=r2list[k]-r2_table[table_i];#define table_four_i TABENERGY(register const int tabtype=-1-(lj_pars->A< 0?lj_pars->A:0);) BigReal kqq=kq_i *p_j-> charge
PDBData * new_PDBData(const char *data)
void get_all_positions(Vector *)
void sprint(char *s, PDBFormatStyle usestyle=COLUMNS)
void getAtom(int num, int *residue_number, char *residue_name, char *atom_name, char *atom_type, int *atom_typenum, Real *charge, Real *mass) const
static const BigReal default_occupancy