PDBData.C

Go to the documentation of this file.
00001 
00007 /*
00008    The code to implement the various PDBData classes.  These are:
00009    PDBData (base class), PDBUnknown (for unimplemented records),
00010    and PDBAtom, PDBAtomRecord, and PDBHetAtm, for the various types
00011    of atom records.
00012 */
00013 
00014 // Here are the routines to manupulate a PDB ATOM record
00015 // It can be created by hand or with a string in the PDB format
00016 
00017 #include <stdio.h>    // sprintf and sscanf
00018 #ifndef WIN32
00019 #include <strings.h>  // strncpy
00020 #endif
00021 #include <stdlib.h>   // atoi and atof
00022 #include "InfoStream.h"
00023 #include "PDBData.h"
00024 
00025 //  Define some constants for the class
00026 const int PDBAtom::default_serial = -1;
00027 const int PDBAtom::default_residueseq = -1;
00028 const BigReal PDBAtom::default_coor = 9999.000; 
00029 const BigReal PDBAtom::default_occupancy = 1.00;
00030 const BigReal PDBAtom::default_temperaturefactor = 0.00;
00031 const int PDBAtom::no_footnote = 0;
00032 
00033 // write down the names so I won't have to do so again
00034 const char *PDBData::PDBNames[PDBData::UNKNOWN+1] = {
00035        "HEADER", "OBSLTE", "COMPND",
00036        "SOURCE", "EXPDTA", "AUTHOR", "REVDAT", "SPRSDE", "JRNL",
00037        "REMARK", "SEQRES", "FTNOTE", "HET", "FORMUL", "HELIX",
00038        "SHEET", "TURN", "SSBOND", "SITE", "CRYST1", "ORIGX",
00039        "SCALE", "MTRIX", "TVECT", "MODEL", "ATOM", "HETATM",
00040        "SIGATM", "ANISOU", "SIGUIJ", "TER", "ENDMDL", "CONECT",
00041        "MASTER", "END", "UNKNOWN"};
00042 
00043 // Parse the input, a char *, for an integer.
00044 //   The input string is "data" and is of length 'length'.
00045 //   The integer starts at position "start" (ASSUMES 1st character
00046 // is at location 1 !*!*!) and is at most "size" characters long.
00047 //   If the string is not long enough, return the default integer value.
00048 //Example: scan("test 12345", 10, 5, 2, &tempint, 99);
00049 //  would set tempint to 12.
00050 //NOTE:  a blank string ("     ") will return the value 0.  If
00051 //   there is more than one integer in the string, it will return
00052 //   the first one seen and not complain.
00053 void PDBData::scan( const char *data, int length, int start, int size, 
00054                        int *ans, int defalt)
00055 {
00056   char tempbuffer[200];       // temp. string buffer
00057   if (length < start) {       // check if the string is long enough
00058     *ans = defalt;                   // return the default
00059     return;
00060   }
00061   if (size>199)                      // make sure I won't overflow my array
00062     size=199;
00063   strncpy(tempbuffer, data + start-1, size);  // convert the string to an int
00064   tempbuffer[size]= 0;
00065   int flg=0;
00066   for (int i=strlen(tempbuffer)-1; i>=0; i--) { // see if this is a blank string
00067     if (tempbuffer[i]>' ') {
00068        flg = 1;
00069        break;
00070     }
00071   }
00072   if (flg != 1) {  // then it was a blank string
00073     *ans = defalt;
00074   } else {
00075   *ans = atoi(tempbuffer);
00076   }
00077 }
00078 
00079 //  Parse the input for a string; as above, but there is no default string.
00080 // I assume that "*ans" has at least 'size+1' characters.
00081 void PDBData::scan( const char *data, int length, int start, int size, char *ans)
00082 {
00083   if (length < start) {             // check if the string is long enough
00084    ans[0] = 0;
00085    return;
00086   }
00087   strncpy(ans, data + start - 1, size);
00088   ans[size]=0;
00089                       //  fix up leading and trailing spaces
00090   int i,j;
00091   for (i=0,j=0; ans[i]; i++)  //   all this does is strip out _all_ the
00092     if (ans[i]!=' ' && ans[i]!='\t')   //   spaces -- this is important because I
00093       ans[j++]=ans[i];  // check that a string is empty by looking to see if 
00094   ans[j]=0;   // [0] == 0 instead of checking to see if all the elements are spaces
00095 }
00096 //  Parse the input for a BigReal
00097 void PDBData::scan( const char *data, int length, int start, int size,
00098                         BigReal *ans, BigReal defalt)
00099 {
00100   char tempbuffer[200];
00101   if (length < start) {               // check if the string is long enough
00102    *ans = defalt;                     // return the default
00103    return;
00104   }
00105   if (size>199)                       // make sure I won't overflow my array
00106     size=199;
00107   strncpy(tempbuffer, data + start - 1, size);// convert the string to a BigReal
00108   tempbuffer[size]= 0;
00109   int flg=0;
00110   for (int i=strlen(tempbuffer)-1; i>=0; i--) { // see if this is a blank string
00111     if (tempbuffer[i]>' ') {
00112        flg = 1;
00113        break;
00114     }
00115   }
00116   if (flg != 1) {  // then it was a blank string
00117     *ans = defalt;
00118   } else {
00119     *ans = atof(tempbuffer);  // WARNING : ASSUMES BigReal <= double!!!
00120   }
00121 }
00122 
00124 // Parse the input data looking for field number 'N'
00125 // where a field is defined as a collection of
00126 // non-whitespace (ws == {space, tab}) characters
00127 // and each field is seperated by one or more whitespace characters
00128 // The result is either
00129 //  1)  field 'N' if it exists
00130 //  2)  some string starting with '#' if it does not
00131 // Note: this means that you can use # to mark fields which don't exist
00132 // Also note that this ASSUMES the first field is field 1 !*!*!*
00133 void PDBData::field(const char *data, int fld, char *result)
00134 {
00135   int i;
00136 
00137   Bool onword = FALSE;
00138   if (fld<=0)  {     // ask a stupid question, get a stupid answer
00139     result[0]='#';
00140     result[1]=0;
00141     return;
00142   }
00143   for (i=0; data[i]; i++)
00144     if (!onword && data[i] != ' ' && data[i] != '\t')  {  // if I found a field
00145        onword = TRUE;     // mark that I'm on it
00146        if (--fld <= 0)    // am I done?
00147          break;
00148     } else {
00149       if (onword && (data[i] == ' ' || data[i] == '\t')) { // left a field
00150         onword = FALSE;   // mark that I left
00151       }
00152     }
00153   if (fld>0) {  // oh no, didn't find the field!
00154      result[0] = '#';
00155      result[1] = 0;
00156      return;
00157   }
00158   
00159   int cpy=0;  // copy the field to the output
00160   while (data[i] != ' ' && data[i] != '\t' && data[i])
00161      result[cpy++] = data[i++];
00162   result[cpy] = 0;  // terminate and I'm done
00163 }
00164 
00165 
00167 // Note that this ASSUMES the first character is column 1 !*!*!*
00168 // print an integer
00169 void PDBData::sprintcol( char *s, int start, int len, int val)
00170 {
00171  char temps[100];
00172  sprintf(temps, "%*d", len, val);  // convert the int to a string
00173  sprintcol( s, start, len, temps); // copy to the output string
00174 }
00175 // print a BigReal
00176 void PDBData::sprintcol( char *s, int start, int len, int prec, BigReal val)
00177 {
00178  char temps[100];
00179  sprintf(temps, "%*.*f", len, prec, val);
00180  sprintcol( s, start, len, temps); // copy to the output string
00181 }
00182 // print a string
00183 void PDBData::sprintcol( char *s, int start, int len, const char *val)
00184 {
00185  s+=start-1;
00186  while (len-- >0 && *val)  // copy string up to end of string or len
00187   *s++ = *val++;
00188 }
00189 
00190 /*********************************************************/
00191 //   base class for both PDB ATOM and HETATM      //
00192 /*********************************************************/
00193 void PDBAtom::parse(const char *data)
00194 {
00195  char tempstr[100];
00196  field(data, 1, tempstr);  // get info about field #1 (the first one)
00197  if (tempstr[0] == '#')  {
00198    parse_field_data(data);
00199  } else {
00200    parse_column_data(data);
00201  } 
00202 }
00203 PDBAtom::PDBAtom( const char *data, PDBPossibleAtoms atomclass)
00204         : PDBData( atomclass==USE_HETATM ? PDBData::HETATM : PDBData::ATOM)
00205 {
00206  parse(data);
00207 };
00208 
00209 //  This constructor does nothing except default the record
00210 PDBAtom::PDBAtom( void) : PDBData( PDBData:: ATOM)
00211 {
00212   serialnumber(default_serial);
00213   name("");
00214   alternatelocation("");
00215   residuename("");
00216   chain("");
00217   residueseq(default_residueseq);
00218   insertioncode("");
00219   xcoor(default_coor);
00220   ycoor(default_coor);
00221   zcoor(default_coor);
00222   occupancy(default_occupancy);
00223   temperaturefactor(default_temperaturefactor);
00224   footnote(no_footnote);
00225   segmentname("");
00226   element("");
00227 }
00228 
00229 
00230 PDBAtom::~PDBAtom( void) {
00231 }
00232 
00233 
00234 // Create an atom or hetatm record given that it is column based data
00235 void PDBAtom::parse_column_data( const char *data)
00236 {
00237  int len = strlen(data);  // to check that there is info
00238  char tempstr[100];   
00239  int tempint;
00240  BigReal tempBigReal;
00241 
00242      // set the serial number
00243  scan(data, len, SSERIAL, LSERIAL, &tempint, default_serial);
00244  serialnumber( tempint );
00245  
00246      // set the name
00247  scan(data, len, SNAME, LNAME, tempstr);
00248  name( tempstr);
00249  
00250      // set the alternate location
00251  scan(data, len, SALT, LALT, tempstr);
00252  alternatelocation( tempstr);
00253 
00254      // set the residue name
00255  scan(data, len, SRESNAME, LRESNAME, tempstr);
00256  residuename( tempstr);
00257  
00258      // set the chain 
00259  scan(data, len, SCHAIN, LCHAIN, tempstr);
00260  chain( tempstr);
00261  
00262      // set the residue sequence
00263  {
00264     // If there are more than 9999 residues, X-Plor uses A000, A001,
00265     // ..., A999, B000.  Since we try to be X-Plor compatible,
00266     // we doo the following
00267     char s[10];
00268     scan(data, len, SRESSEQ, LRESSEQ, s);
00269     if (s[0] < '0' || s[0] > '9') {
00270        static int elvis_count = 0;
00271        int num = (s[0] - 'A') * 1000 + 10000;
00272        num += atoi(s+1);
00273        if (s[0] < 'A' || s[0] > 'Z') {
00274           if (elvis_count == 0) {
00275      iout << iWARN << "Man, tiny Elvis, that number is huge!\n"
00276      << iWARN << "We don't know how X-Plor represents over Z999 residues\n"
00277      << iWARN << "And you just tried " << s << " - so we'll fake it as " << num << "\n"
00278      << iWARN << "This is reversible, but only inside this program.\n" << endi;
00279                elvis_count = 1;
00280            }
00281        } else {
00282           elvis_count = 0;
00283        }
00284        residueseq( num);
00285     } else {
00286        residueseq(atoi(s));
00287     }
00288     // This used to be simply
00289 //    scan(data, len, SRESSEQ, LRESSEQ, &tempint, default_residueseq);
00290 //    residueseq( tempint);
00291  }
00292  
00293      // set the insertion code 
00294  scan(data, len, SINSERT, LINSERT, tempstr);
00295  insertioncode( tempstr);
00296  
00297      // set the X, Y, and Z coordinates
00298  scan(data, len, SX, LCOOR, &tempBigReal, default_coor);
00299  xcoor( tempBigReal);
00300  scan(data, len, SY, LCOOR, &tempBigReal, default_coor);
00301  ycoor( tempBigReal);
00302  scan(data, len, SZ, LCOOR, &tempBigReal, default_coor);
00303  zcoor( tempBigReal);
00304 
00305      // set the occupancy 
00306  scan(data, len, SOCC, LOCC, &tempBigReal, default_occupancy);
00307  occupancy( tempBigReal);
00308  
00309      // set the temperature factor
00310  scan(data, len, STEMPF, LTEMPF, &tempBigReal, default_temperaturefactor);
00311  temperaturefactor( tempBigReal);
00312  
00313      // set the footnote
00314  scan(data, len, SFOOT, LFOOT, &tempint, no_footnote);
00315  footnote( tempint);
00316 
00317    // this is for XPLOR style PDBs which have a segment name
00318  scan(data, len, SSEGNAME, LSEGNAME, tempstr);
00319  segmentname( tempstr);
00320 
00321  scan(data, len, SELEMENT, LELEMENT, tempstr);
00322  element( tempstr);
00323 }
00324 
00325 void PDBAtom::parse_field_data( const char *data)
00326 {
00327   char tempstr[100];
00328   // I already know that the first field starts with a '#' and that
00329   // the second is either ATOM or HETATM, so I'll start with the third
00330   field(data, 3, tempstr);
00331   serialnumber( tempstr[0] != '#' ? atoi(tempstr) : default_serial );
00332 
00333   field(data, 4, tempstr);
00334   name( tempstr[0] != '#' ? tempstr : "" );
00335 
00336   field(data, 5, tempstr);
00337   alternatelocation( tempstr[0] != '#' ? tempstr : "" );
00338 
00339   field(data, 6, tempstr);
00340   residuename( tempstr[0] != '#' ? tempstr : "" );
00341                 
00342   field(data, 7, tempstr);
00343   chain( tempstr[0] != '#' ? tempstr : "" );
00344 
00345   field(data, 8, tempstr);
00346   residueseq( tempstr[0] != '#' ? atoi(tempstr) : default_residueseq );
00347 
00348   field(data, 9, tempstr);
00349   insertioncode( tempstr[0] != '#' ? tempstr : "" );
00350 
00351   field(data, 10, tempstr);
00352   xcoor( tempstr[0] != '#' ?
00353         atof( tempstr) : default_coor);  // WARNING: assumes BigReal <= double
00354   field(data, 11, tempstr);
00355   ycoor( tempstr[0] != '#' ?
00356         atof( tempstr) : default_coor);  // WARNING: assumes BigReal <= double
00357   field(data, 12, tempstr);
00358   zcoor( tempstr[0] != '#' ?
00359         atof( tempstr) : default_coor);  // WARNING: assumes BigReal <= double
00360 
00361   field(data, 13, tempstr);
00362   occupancy( tempstr[0] != '#' ? 
00363         atof( tempstr) : default_occupancy );// WARNING: assumes BigReal <= double
00364 
00365   field(data, 14, tempstr);
00366   temperaturefactor( tempstr[0] != '#' ?
00367         atof( tempstr) : default_temperaturefactor ); // WARNING: ditto
00368 
00369   field(data, 15, tempstr);
00370   footnote( tempstr[0] != '#' ? atoi(tempstr) : no_footnote);
00371   
00372   field(data, 16, tempstr);
00373   segmentname( tempstr[0] != '#' ? tempstr : "");
00374   
00375   field(data, 17, tempstr);
00376   element( tempstr[0] != '#' ? tempstr : "");
00377 }
00378 
00379   // get/ set the serial number
00380 int PDBAtom:: serialnumber( void) 
00381 { return myserialnumber; }
00382 void PDBAtom:: serialnumber( int newserialnumber)
00383 { myserialnumber = newserialnumber; }
00384 
00385   // get/ set the serial number
00386 const char* PDBAtom:: name( void)
00387 { return myname; }
00388 void PDBAtom:: name( const char *newname)
00389 { strncpy(myname, newname, LNAME); myname[LNAME]=0; }
00390 
00391   // get/ set the alternate location
00392 const char* PDBAtom:: alternatelocation( void)
00393 { return myalternatelocation; }
00394 void PDBAtom:: alternatelocation( const char *newalternatelocation)
00395 { strncpy(myalternatelocation, newalternatelocation, LALT);
00396   myalternatelocation[LALT]=0;}
00397 
00398   // get/ set the residue name
00399 const char* PDBAtom:: residuename( void)
00400 { return myresiduename; }
00401 void PDBAtom:: residuename( const char *newresiduename)
00402 { strncpy(myresiduename, newresiduename, LRESNAME); myresiduename[LRESNAME]=0;}
00403 
00404   // get/ set the chain indentifier
00405 const char* PDBAtom:: chain( void)
00406 { return mychain; }
00407 void PDBAtom:: chain( const char *newchain)
00408 { strncpy(mychain, newchain, LCHAIN); mychain[LCHAIN]=0;}
00409 
00410   // get/ set the residue sequence number
00411 int PDBAtom:: residueseq( void)
00412 { return myresidueseq; }
00413 void PDBAtom:: residueseq( int newresidueseq)
00414 { myresidueseq = newresidueseq; }
00415 
00416   // get/ set the insertion code
00417 const char* PDBAtom:: insertioncode( void)
00418 { return myinsertioncode; }
00419 void PDBAtom:: insertioncode( const char *newinsertioncode)
00420 { strncpy(myinsertioncode, newinsertioncode, LINSERT); 
00421   myinsertioncode[LINSERT]=0;}
00422 
00423   // get/ set the different coordinates
00424   // either 1 by 1 ...
00425 BigReal PDBAtom:: xcoor( void)
00426 { return mycoor[0]; }
00427 void PDBAtom:: xcoor( BigReal newxcoor)
00428 { mycoor[0] = newxcoor; }
00429 BigReal PDBAtom:: ycoor( void)
00430 { return mycoor[1]; }
00431 void PDBAtom:: ycoor( BigReal newycoor)
00432 { mycoor[1] = newycoor; }
00433 BigReal PDBAtom:: zcoor( void)
00434 { return mycoor[2]; }
00435 void PDBAtom:: zcoor( BigReal newzcoor)
00436 { mycoor[2] = newzcoor; }
00437    // ...or all three at once
00438 const BigReal *PDBAtom:: coordinates( void)
00439 { return mycoor; }
00440 void PDBAtom:: coordinates(const BigReal *newcoordinates)
00441 { for (int i=0; i<3; i++) mycoor[i] = newcoordinates[i]; }
00442 
00443   // get/ set the occupancy
00444 BigReal PDBAtom:: occupancy( void)
00445 { return myoccupancy ;}
00446 void PDBAtom:: occupancy( BigReal newoccupancy)
00447 { myoccupancy = newoccupancy; }
00448 
00449   // get/ set the temperature factor
00450 BigReal PDBAtom:: temperaturefactor( void)
00451 { return mytemperaturefactor; }
00452 void PDBAtom:: temperaturefactor( BigReal newtemperaturefactor)
00453 { mytemperaturefactor = newtemperaturefactor; }
00454 
00455   // get/ set the footnote
00456 int PDBAtom:: footnote( void)
00457 { return myfootnote; }
00458 void PDBAtom:: footnote( int newfootnote)
00459 { myfootnote = newfootnote; }
00460 
00461   // get/ set the segment name
00462   // this is not part of the PDB format but is used by XPLOR instead of
00463   // the chain identifier (see XPLOR 3.1 manual, p 104)
00464 const char* PDBAtom:: segmentname( void)
00465 { return mysegmentname; }
00466 void PDBAtom:: segmentname( const char *newsegmentname)
00467 { strncpy(mysegmentname, newsegmentname, LSEGNAME); mysegmentname[LSEGNAME]=0;}
00468 
00469   // get/ set the element name
00470 const char* PDBAtom:: element( void)
00471 { return myelement; }
00472 void PDBAtom:: element( const char *newelement)
00473 { strncpy(myelement, newelement, LELEMENT); myelement[LELEMENT]=0;}
00474 
00475  
00476 // the function to print out an ATOM or HETATM
00477 // size or outstr must be >= 80!
00478 void PDBAtom::sprint_columns( char *outstr)
00479 { 
00480   int i;
00481   for (i=0; i<79; i++)   // dump spaces in outstr -- must be length > 80!!
00482     outstr[i] = 32;
00483   outstr[i] = 0;             // and terminate
00484 
00485   sprintcol(outstr, STYPE, LTYPE, PDBNames[type()] );
00486   sprintcol(outstr, SSERIAL, LSERIAL, serialnumber());
00487   {
00488      // For X-Plor compatability, if the name is 1, 2, or 3
00489      // characters, start it in the 2nd column of the field.
00490      // 4 letter names use the first column
00491      if (strlen(name()) == 4) {
00492         sprintcol(outstr, SNAME, LNAME, name());
00493      } else {
00494         sprintcol(outstr, SNAME+1, LNAME-1, name());
00495      }
00496   }
00497   sprintcol(outstr, SALT, LALT, alternatelocation());
00498   sprintcol(outstr, SRESNAME, LRESNAME, residuename());
00499   sprintcol(outstr, SCHAIN, LCHAIN, chain());
00500   {
00501      // Again, I may have to convert from a number > 9999 to
00502      // A000 or whatever (see the comments for the residueseq input)
00503      if (residueseq() <= 9999) {
00504         sprintcol(outstr, SRESSEQ, LRESSEQ, residueseq());
00505      } else {
00506         int val = residueseq() / 1000 - 10;  // integer arithmetic
00507         int modulo = residueseq() % 1000;
00508         char s[10];
00509         sprintf(s, "%c%03d", 'A' + val, modulo);
00510         sprintcol(outstr, SRESSEQ, LRESSEQ, s);
00511      }
00512      // This used to be just ...
00513 //  sprintcol(outstr, SRESSEQ, LRESSEQ, residueseq());
00514   }
00515   sprintcol(outstr, SINSERT, LINSERT, insertioncode());
00516   sprintcol(outstr, SX, LCOOR, LCOORPREC, xcoor());
00517   sprintcol(outstr, SY, LCOOR, LCOORPREC, ycoor());
00518   sprintcol(outstr, SZ, LCOOR, LCOORPREC, zcoor());
00519   sprintcol(outstr, SOCC, LOCC, LOCCPREC, occupancy());
00520   sprintcol(outstr, STEMPF, LTEMPF, LTEMPFPREC, temperaturefactor());
00521   if (footnote() == no_footnote)     // special case when no footnote
00522     sprintcol(outstr, SFOOT, LFOOT, "");
00523    else
00524     sprintcol(outstr, SFOOT, LFOOT, footnote() );         
00525   sprintcol(outstr, SSEGNAME, LSEGNAME, segmentname());
00526   // world's lamest right-justify
00527   int lelement = strlen(element());
00528   lelement = ( lelement > LELEMENT ? LELEMENT : lelement );
00529   sprintcol(outstr, SELEMENT+(LELEMENT-lelement), lelement, element());
00530 }
00531 
00533 void PDBAtom::sprint_fields( char *outstr)
00534 {
00535   char tmpstr[50];
00536   sprintf(outstr, "# %s", PDBNames[type()]);
00537   if (serialnumber() == default_serial)
00538       sprintf(tmpstr, " #");
00539      else
00540       sprintf(tmpstr, " %i", serialnumber());
00541   strcat(outstr, tmpstr);
00542   if (name()[0] == 0)
00543       sprintf(tmpstr, " #");
00544      else
00545       sprintf(tmpstr, " %s", name());
00546   strcat(outstr, tmpstr);
00547   if (alternatelocation()[0] == 0)
00548       sprintf(tmpstr, " #");
00549      else
00550       sprintf(tmpstr, " %s", alternatelocation());
00551   strcat(outstr, tmpstr);
00552   if (residuename()[0] == 0)
00553       sprintf(tmpstr, " #");
00554      else
00555       sprintf(tmpstr, " %s", residuename());
00556   strcat(outstr, tmpstr);
00557   if (chain()[0] == 0)
00558       sprintf(tmpstr, " #");
00559      else
00560       sprintf(tmpstr, " %s", chain());
00561   strcat(outstr, tmpstr);
00562   if (residueseq() == default_residueseq)
00563       sprintf(tmpstr, " #");
00564      else
00565       sprintf(tmpstr, " %d", residueseq());
00566   strcat(outstr, tmpstr);
00567   if (insertioncode()[0] == 0)
00568       sprintf(tmpstr, " #");
00569      else
00570       sprintf(tmpstr, " %s", insertioncode());
00571   strcat(outstr, tmpstr);
00572   if (xcoor() == default_coor)
00573       sprintf(tmpstr, " #");
00574      else
00575       sprintf(tmpstr, " %*.*f", LCOOR, LCOORPREC, xcoor());
00576   strcat(outstr, tmpstr);
00577   if (ycoor() == default_coor)
00578       sprintf(tmpstr, " #");
00579      else
00580       sprintf(tmpstr, " %*.*f", LCOOR, LCOORPREC, ycoor());
00581   strcat(outstr, tmpstr);
00582   if (zcoor() == default_coor)
00583       sprintf(tmpstr, " #");
00584      else
00585       sprintf(tmpstr, " %*.*f", LCOOR, LCOORPREC, zcoor());
00586   strcat(outstr, tmpstr);
00587 //  if (occupancy() == default_occupancy)  // no way to tell if the occ. is the default
00588 //      sprintf(tmpstr, " #");
00589 //     else
00590       sprintf(tmpstr, " %*.*f",  LOCC, LOCCPREC, occupancy());
00591   strcat(outstr, tmpstr);
00592 //  if (temperaturefactor() == default_temperaturefactor) // ditto previous
00593 //      sprintf(tmpstr, " #");
00594 //     else
00595       sprintf(tmpstr, " %*.*f", LTEMPF, LTEMPFPREC, temperaturefactor());
00596   strcat(outstr, tmpstr);
00597   if (footnote() == no_footnote)     // special case when no footnote
00598     sprintf(tmpstr, " #");
00599    else
00600     sprintf(tmpstr, " %d",  footnote() );         
00601   strcat(outstr, tmpstr);
00602   if (segmentname()[0] == 0)
00603       sprintf(tmpstr, " #");
00604      else
00605       sprintf(tmpstr, " %s", segmentname());
00606   strcat(outstr, tmpstr);
00607   if (element()[0] == 0)
00608       sprintf(tmpstr, " #");
00609      else
00610       sprintf(tmpstr, " %s", element());
00611   strcat(outstr, tmpstr);
00612 
00613 }
00614 
00615 void PDBAtom::sprint( char *outstr, PDBFormatStyle usestyle)
00616 {
00617  if (usestyle == PDBData::COLUMNS)
00618    sprint_columns( outstr);
00619  else
00620    sprint_fields( outstr);
00621 }
00622 
00623 //****************** The wrapper for all of the functions ************///
00624 PDBData *new_PDBData(const char *data)  // nasty
00625 {
00626   char temps1[160];
00627   char temps2[160];
00628   char *temps;
00629   sscanf(data, "%s %s ", temps1, temps2);
00630   if (temps1[0] == '#')
00631      temps = temps2;
00632     else
00633      temps = temps1;
00634 
00635      // go through the list of possible PDB data types
00636      //this _should_ be the same as: for(PDBTypes i=HEADER; i<UNKNOWN; i++)
00637   for (int i=0; i< (int)(sizeof(PDBData::PDBNames) /
00638                   sizeof(PDBData::PDBNames[0])); i++)
00639     if (!strcmp(temps, PDBData::PDBNames[i]))
00640        switch(i) {
00641          case PDBData::ATOM: return new PDBAtomRecord(data);
00642          case PDBData::HETATM: return new PDBHetatm(data);
00643          default: return new PDBUnknown(data);
00644        }
00645   // Now, if HETATM is right next to an aton number (like HETATM12345) then the above
00646   // test will fail, so I have to special case it:
00647   if (!strncmp(temps, PDBData::PDBNames[PDBData::HETATM], sizeof(PDBData::PDBNames[PDBData::HETATM]))) {
00648     return new PDBHetatm(data);
00649   }
00650    //  Hmm, looks like it isn't any data type, so I'll fake it
00651    return new PDBUnknown(data);
00652 }
00653 
00655 //#define TEST_PDBREADER
00656 #ifdef TEST_PDBREADER
00657 main()
00658 {
00659   char tempstr[100];
00660   PDBAtomRecord atom("ATOM   6312  CB TALA 3 235I     24.681  54.463 137.827  1.00 51.30      VP3");
00661   atom.sprint(tempstr, PDBData::COLUMNS);
00662   std::cout << tempstr << '\n';
00663   atom.sprint(tempstr, PDBData::FIELDS);
00664   std::cout << tempstr << '\n';
00665   std::cout << "Serial number : "  << atom.serialnumber()         <<  "\n";
00666   std::cout << "name          : '" << atom.name()                 <<  "'\n";
00667   std::cout << "alt. location : '" << atom.alternatelocation()    <<  "'\n";
00668   std::cout << "residue name  : '" << atom.residuename()          <<  "'\n";
00669   std::cout << "chain         : '" << atom.chain()                <<  "'\n";
00670   std::cout << "residue seq   : "  << atom.residueseq()           <<  "\n";
00671   std::cout << "insertion code: '" << atom.insertioncode()        << "'\n";
00672   std::cout << "X coordinate  : "  << atom.xcoor()                <<  "\n";
00673   std::cout << "Y coordinate  : "  << atom.ycoor()                << "\n";
00674   std::cout << "Z coordinate  : "  << atom.zcoor()                << "\n";
00675   std::cout << "occupancy     : "  << atom.occupancy()            << "\n";
00676   std::cout << "temperature factor: " << atom.temperaturefactor() << "\n";
00677   std::cout << "footnote      : " << atom.footnote()              << "\n";
00678   std::cout << "segment name  : '" << atom.segmentname()          << "'\n";
00679   std::cout << "element       : '" << atom.element()              << "'\n";
00680   std::cout << '\n';
00681   
00682   PDBAtomRecord atom2("# ATOM 6312 CB T ALA 3 235 I 24.681 54.463 137.827 1.00 51.30 # VP3");
00683   atom2.sprint(tempstr, PDBData::COLUMNS);
00684   std::cout << tempstr << '\n';
00685   atom2.sprint(tempstr, PDBData::FIELDS);
00686   std::cout << tempstr << '\n';
00687   std::cout << "Serial number : "  << atom2.serialnumber()         <<  "\n";
00688   std::cout << "name          : '" << atom2.name()                 <<  "'\n";
00689   std::cout << "alt. location : '" << atom2.alternatelocation()    <<  "'\n";
00690   std::cout << "residue name  : '" << atom2.residuename()          <<  "'\n";
00691   std::cout << "chain         : '" << atom2.chain()                <<  "'\n";
00692   std::cout << "residue seq   : "  << atom2.residueseq()           <<  "\n";
00693   std::cout << "insertion code: '" << atom2.insertioncode()        << "'\n";
00694   std::cout << "X coordinate  : "  << atom2.xcoor()                <<  "\n";
00695   std::cout << "Y coordinate  : "  << atom2.ycoor()                << "\n";
00696   std::cout << "Z coordinate  : "  << atom2.zcoor()                << "\n";
00697   std::cout << "occupancy     : "  << atom2.occupancy()            << "\n";
00698   std::cout << "temperature factor: " << atom2.temperaturefactor() << "\n";
00699   std::cout << "footnote      : " << atom2.footnote()              << "\n";
00700   std::cout << "segment name  : '" << atom2.segmentname()          << "'\n";
00701   std::cout << "element       : '" << atom2.element()              << "'\n";
00702   std::cout << '\n';
00703   
00704 
00705   PDBAtomRecord atom3("# ATOM # # Q WER # # # # 123.456 # # # 9 LAST anything?");
00706   atom3.sprint(tempstr, PDBData::COLUMNS);
00707   std::cout << tempstr << '\n';
00708   atom3.sprint(tempstr, PDBData::FIELDS);
00709   std::cout << tempstr << '\n';
00710   std::cout << "Serial number : "  << atom3.serialnumber()         <<  "\n";
00711   std::cout << "name          : '" << atom3.name()                 <<  "'\n";
00712   std::cout << "alt. location : '" << atom3.alternatelocation()    <<  "'\n";
00713   std::cout << "residue name  : '" << atom3.residuename()          <<  "'\n";
00714   std::cout << "chain         : '" << atom3.chain()                <<  "'\n";
00715   std::cout << "residue seq   : "  << atom3.residueseq()           <<  "\n";
00716   std::cout << "insertion code: '" << atom3.insertioncode()        << "'\n";
00717   std::cout << "X coordinate  : "  << atom3.xcoor()                <<  "\n";
00718   std::cout << "Y coordinate  : "  << atom3.ycoor()                << "\n";
00719   std::cout << "Z coordinate  : "  << atom3.zcoor()                << "\n";
00720   std::cout << "occupancy     : "  << atom3.occupancy()            << "\n";
00721   std::cout << "temperature factor: " << atom3.temperaturefactor() << "\n";
00722   std::cout << "footnote      : " << atom3.footnote()              << "\n";
00723   std::cout << "segment name  : '" << atom3.segmentname()          << "'\n";
00724   std::cout << "element       : '" << atom3.element()              << "'\n";
00725   std::cout << '\n';
00726   
00727 }
00728 #endif
00729 

Generated on Sun Nov 19 01:17:14 2017 for NAMD by  doxygen 1.4.7