src/cpsp/PDBsupport.hh

Go to the documentation of this file.
00001 #ifndef PDBSUPPORT_HH_
00002 #define PDBSUPPORT_HH_
00003 
00004 #include <iostream>
00005 #include <vector>
00006 #include <string>
00007 #include <iomanip>
00008 #include <time.h>
00009 
00010 
00011 #include <biu/LatticeModel.hh>
00012 
00013 namespace cpsp {
00014 
00015 
00016 
00017 std::string getTimeString() {
00018       // return string
00019     std::string strTime = "";
00020       // get current time
00021     time_t tim = time(NULL);
00022     tm *curTime = localtime(&tim);
00023       // generate time string
00024     strTime += char(48+curTime->tm_mday/10);
00025     strTime += char(48+curTime->tm_mday%10);
00026     strTime += '-';
00027     switch (curTime->tm_mon) {
00028         case 0  : strTime += "JAN"; break;
00029         case 1  : strTime += "FEB"; break;
00030         case 2  : strTime += "MAR"; break;
00031         case 3  : strTime += "APR"; break;
00032         case 4  : strTime += "MAY"; break;
00033         case 5  : strTime += "JUN"; break;
00034         case 6  : strTime += "JUL"; break;
00035         case 7  : strTime += "AUG"; break;
00036         case 8  : strTime += "SEP"; break;
00037         case 9  : strTime += "OCT"; break;
00038         case 10 : strTime += "NOV"; break;
00039         case 11 : strTime += "DEC"; break;
00040     }
00041     strTime += '-';
00042     int year10 = curTime->tm_year%10;
00043     int year100 = (curTime->tm_year-year10)%100;
00044     strTime += char(48+year100);
00045     strTime += char(48+year10);
00046       // return time string
00047     return strTime;
00048 }
00049 
00050 
00051 void
00052 writePDBseq( const std::vector<std::string> & seq,
00053             char chainID,
00054             std::ostream & out) 
00055 {
00056       // print sequence information
00057     size_t seqResLine = 1;
00058     out <<"SEQRES "<<std::setw(3) <<seqResLine++ <<" " <<chainID <<" "<<std::setw(4) <<seq.size()<<" ";
00059     for (size_t i=1; i<=seq.size(); i++) {
00060         out <<" " <<seq[i-1];
00061         if (i%13==0 && (i+1) < seq.size()) {
00062             out <<"\nSEQRES "<<std::setw(3) <<seqResLine++ <<" " <<chainID <<" "<<std::setw(4) <<seq.size()<<" ";
00063         }
00064     }
00065     out <<"\n";
00066 }
00067 
00068 void
00069 writePDBhetatm( const std::vector<std::string> & seq,
00070             const biu::DPointVec & points, 
00071             char chainID,
00072             bool hasSideChain,
00073             std::ostream & out,
00074             size_t firstID) 
00075 {
00076     size_t atomID = 0;
00077     for (size_t i=0; i<seq.size(); i++) {
00078         
00079         out <<"HETATM" <<std::setw(5)<<(atomID+1+firstID) <<"  "<<"CA " <<" "<<seq[i]
00080             <<" " <<chainID <<std::setw(4) <<(i+1) 
00081             <<" " <<"   " <<std::fixed 
00082             <<std::setprecision(3) <<std::setw(8) <<points[atomID].getX() 
00083             <<std::setprecision(3) <<std::setw(8) <<points[atomID].getY() 
00084             <<std::setprecision(3) <<std::setw(8) <<points[atomID].getZ()
00085             <<std::setprecision(2) <<std::setw(6) <<1.0
00086             <<std::setprecision(2) <<std::setw(6) <<0.0
00087             <<"           C  "
00088             <<"\n";
00089         atomID++;
00090         if (hasSideChain) {
00091             out <<"HETATM" <<std::setw(5)<<(atomID+1+firstID) <<"  "<<"CB " <<" "<<seq[i]
00092                 <<" " <<chainID <<std::setw(4) <<(i+1) 
00093                 <<" " <<"   " <<std::fixed 
00094                 <<std::setw(8) <<std::setprecision(3) <<points[atomID].getX() 
00095                 <<std::setw(8) <<std::setprecision(3) <<points[atomID].getY() 
00096                 <<std::setw(8) <<std::setprecision(3) <<points[atomID].getZ()
00097                 <<std::setw(6) <<std::setprecision(2) <<1.0
00098                 <<std::setw(6) <<std::setprecision(2) <<0.0
00099                 <<"          Cl  "
00100                 <<"\n";
00101             atomID++;
00102         }
00103     }
00104 }
00105 
00106 void
00107 writePDBconect( const std::vector<std::string> & seq,
00108             const biu::DPointVec & points, 
00109             char chainID,
00110             bool hasSideChain,
00111             std::ostream & out,
00112             size_t firstID) 
00113 {
00114 //                  1         2         3         4         5         6         7         8         
00115 //         12345678901234567890123456789012345678901234567890123456789012345678901234567890
00116      // CA -> CB
00117      
00118     if (hasSideChain) {
00119           // give connections explicitly
00120         size_t i = firstID;
00121         out <<"CONECT" <<std::setw(5) <<(2*i)+1 
00122             <<std::setw(5) <<((2*i)+1)+1 <<std::setw(5) <<(2*(i+1))+1
00123             <<"                                                           \n";
00124         out <<"CONECT" <<std::setw(5) <<((2*i)+1)+1 
00125             <<std::setw(5) <<(2*i)+1
00126             <<"                                                                \n";
00127         for (i++; (i+1)< (firstID+seq.size()); i++) {
00128             out <<"CONECT" <<std::setw(5) <<(2*i)+1 
00129                 <<std::setw(5) <<(2*(i-1))+1 <<std::setw(5) <<((2*i)+1)+1 <<std::setw(5) <<(2*(i+1))+1
00130                 <<"                                                      \n";
00131             out <<"CONECT" <<std::setw(5) <<((2*i)+1)+1 
00132                 <<std::setw(5) <<(2*i)+1
00133                 <<"                                                                \n";
00134         }
00135         out <<"CONECT" <<std::setw(5) <<(2*i)+1 
00136             <<std::setw(5) <<(2*(i-1))+1 <<std::setw(5) <<((2*i)+1)+1
00137             <<"                                                           \n";
00138         out <<"CONECT" <<std::setw(5) <<((2*i)+1)+1 
00139             <<std::setw(5) <<(2*i)+1
00140             <<"                                                                \n";
00141     } else {
00142           // give connections explicitly
00143         size_t i = firstID;
00144         out <<"CONECT" <<std::setw(5) <<i+1 
00145             <<std::setw(5) <<(i+1)+1
00146             <<"                                                                \n";
00147         for (i++; (i+1)< (firstID+seq.size()); i++) {
00148             out <<"CONECT" <<std::setw(5) <<i+1 
00149                 <<std::setw(5) <<i <<std::setw(5) <<(i+1)+1
00150                 <<"                                                           \n";
00151         }
00152         out <<"CONECT" <<std::setw(5) <<i+1 
00153             <<std::setw(5) <<i
00154             <<"                                                                \n";
00155     }   
00156 }
00157 void 
00158 writePDB(   const std::string & sourceInfo,
00159             const std::vector<std::string> & seq,
00160             const biu::DPointVec & points, 
00161             const biu::LatticeModel & lattice,
00162             std::ostream & out,
00163             bool hasSideChain,
00164             const biu::LatticeModel & sclattice)  
00165 {
00166     const char chainID = 'L';
00167       // get lattice name
00168     std::string latName = lattice.getDescriptor()->getName();
00169     latName.resize(58,' ');
00170 
00171 //                  1         2         3         4         5         6         7         8         
00172 //         12345678901234567890123456789012345678901234567890123456789012345678901234567890
00173     out <<"HEADER    LATTICE PROTEIN STRUCTURE               "<<getTimeString()<<"                     \n"              
00174         <<"TITLE     HP LATTICE PROTEIN STRUCTURE                                          \n" 
00175         <<"COMPND    MOL_ID: 1;                                                            \n"
00176         <<"COMPND   2 MOLECULE: LATTICE PROTEIN;                                           \n"        
00177         <<"COMPND   3 CHAIN: "<<chainID<<";                                                            \n"
00178         <<"COMPND   4 ENGINEERED: YES                                                      \n"
00179         <<"SOURCE    MOL_ID: 1                                                             \n"
00180         <<"KEYWDS    LATTICE PROTEIN MODEL                                                 \n"
00181         <<"EXPDTA    THEORETICAL MODEL                                                     \n"
00182         <<"AUTHOR    CPSP-TOOLS SOFTWARE                                                   \n"               
00183         <<"REMARK  40                                                                      \n"
00184         <<"REMARK  40 GENERATED WITH CPSP-TOOLS (C) MARTIN MANN 2008                       \n";
00185     out <<"REMARK  40  LATTICE = "<<latName <<"\n";
00186     if (hasSideChain) {     
00187         std::string latName = sclattice.getDescriptor()->getName();
00188         latName.resize(47,' ');
00189     out <<"REMARK  40  SIDE CHAIN MODEL                                                    \n"
00190         <<"REMARK  40  SIDE CHAIN LATTICE = " <<latName <<"\n";
00191     }
00192     out <<"REMARK 220                                                                      \n"
00193         <<"REMARK 220 EXPERIMENTAL DETAILS                                                 \n"
00194         <<"REMARK 220 EXPERIMENT TYPE : THEORETICAL MODELLING                              \n";
00195     if (sourceInfo.size() > 0) {
00196         out <<"REMARK 220                                                                      \n"
00197             <<"REMARK 220 PARAMETER SET APPLIED :                                              \n"
00198             <<"REMARK 220                                                                      \n";
00199           // write parameter set out of infoString
00200         size_t cut = 0;
00201         while(cut < sourceInfo.size() && sourceInfo[cut]==' ')
00202         { cut++; }
00203         while(cut < sourceInfo.size() && cut != std::string::npos) {
00204             size_t start = cut;
00205             cut = sourceInfo.find(" ", start+1);
00206             std::string output = "REMARK 220  ";
00207             output += sourceInfo.substr(start,cut-start);
00208             output.resize(80,' ');
00209             out <<output <<"\n";
00210             if (cut!=std::string::npos)
00211                 cut++;
00212         }
00213     }
00214         
00215     out <<"REMARK 220                                                                      \n"
00216         <<"REMARK 225                                                                      \n"
00217         <<"REMARK 225 THEORETICAL MODEL                                                    \n"
00218         <<"REMARK 225 THE COORDINATES IN THIS ENTRY REPRESENT A LATTICE MODEL STRUCTURE.   \n"
00219         <<"REMARK 225                                                                      \n"
00220         <<"REMARK 225 THE LATTICE PROTEIN STRUCTURE IS GIVEN AS CHAIN "<<chainID<<".                   \n"
00221         <<"REMARK 225 THE BACKBONE IS REPRESENTED BY A CA-ATOMS AS A 'C'-elements.         \n"
00222         <<(hasSideChain?"REMARK 225 THE SIDE CHAIN IS REPRESENTED BY A CB-ATOMS AS A 'Cl'-elements.      \n":"")
00223         <<"REMARK 225                                                                      \n"
00224         ;
00225 
00226       // print sequence information
00227     writePDBseq( seq, chainID, out);
00228     
00229       // print coordinates
00230     writePDBhetatm( seq, points, chainID, hasSideChain, out, 0);
00231       // mark end of chain
00232     out <<"TER                                                                             \n";
00233     
00234       // print contacts
00235     writePDBconect( seq, points, chainID, hasSideChain, out, 0);
00236     
00237       // mark EOF
00238     out <<"END   "
00239         <<std::endl;
00240 
00241     out.flush();    
00242 }
00243 
00244 } // namespace cpsp
00245 
00246 #endif /*PDBSUPPORT_HH_*/