31 const std::vector<Vector> & PDB::getPositions()
const{
35 const std::vector<double> & PDB::getOccupancy()
const{
39 const std::vector<double> & PDB::getBeta()
const{
43 const std::vector<std::string> & PDB::getRemark()
const{
47 const std::vector<AtomNumber> & PDB::getAtomNumbers()
const{
52 std::map<AtomNumber,unsigned>::const_iterator p;
53 p=number2index.find(a);
54 if(p==number2index.end())
return "";
55 else return atomsymb[p->second];
59 std::map<AtomNumber,unsigned>::const_iterator p;
60 p=number2index.find(a);
61 if(p==number2index.end())
return 0;
62 else return residue[p->second];
66 std::map<AtomNumber,unsigned>::const_iterator p;
67 p=number2index.find(a);
68 if(p==number2index.end())
return "";
69 else return residuenames[p->second];
72 unsigned PDB::size()
const{
73 return positions.size();
76 bool PDB::readFromFilepointer(FILE *fp,
bool naturalUnits,
double scale){
78 bool file_is_alive=
false;
79 if(naturalUnits) scale=1.0;
82 while(Tools::getline(fp,line)){
85 while(line.length()<80) line.push_back(
' ');
86 string record=line.substr(0,6);
87 string serial=line.substr(6,5);
88 string atomname=line.substr(12,4);
89 string residuename=line.substr(17,3);
90 string chainID=line.substr(21,1);
91 string resnum=line.substr(22,4);
92 string x=line.substr(30,8);
93 string y=line.substr(38,8);
94 string z=line.substr(46,8);
95 string occ=line.substr(54,6);
96 string bet=line.substr(60,6);
98 if(record==
"TER"){file_is_alive=
true;
break;}
99 if(record==
"END"){file_is_alive=
true;
break;}
100 if(record==
"ENDMDL"){file_is_alive=
true;
break;}
101 if(record==
"REMARK"){
102 vector<string> v1; v1=Tools::getWords(line.substr(6));
103 remark.insert(remark.begin(),v1.begin(),v1.end());
105 if(record==
"ATOM" || record==
"HETATM"){
109 Tools::convert(serial,a);
110 Tools::convert(resnum,resno);
111 Tools::convert(occ,o);
112 Tools::convert(bet,b);
113 Tools::convert(x,p[0]);
114 Tools::convert(y,p[1]);
115 Tools::convert(z,p[2]);
118 numbers.push_back(a);
119 number2index[
a]=positions.size();
120 std::size_t startpos=atomname.find_first_not_of(
" \t");
121 std::size_t endpos=atomname.find_last_not_of(
" \t");
122 atomsymb.push_back( atomname.substr(startpos, endpos-startpos+1) );
123 residue.push_back(resno);
124 chain.push_back(chainID);
125 occupancy.push_back(o);
127 positions.push_back(p);
128 residuenames.push_back(residuename);
131 return file_is_alive;
134 void PDB::setArgKeyword(
const std::string& new_args ){
136 for(
unsigned i=0;i<remark.size();++i){
137 if( remark[i].find(
"ARG=")!=std::string::npos){
138 remark[i]=new_args; replaced=
true;
141 plumed_assert( replaced );
144 bool PDB::read(
const std::string&file,
bool naturalUnits,
double scale){
145 FILE* fp=fopen(file.c_str(),
"r");
146 if(!fp)
return false;
147 readFromFilepointer(fp,naturalUnits,scale);
152 void PDB::getChainNames( std::vector<std::string>& chains )
const {
154 chains.push_back( chain[0] );
155 for(
unsigned i=1;i<size();++i){
156 if( chains[chains.size()-1]!=chain[i] ) chains.push_back( chain[i] );
160 void PDB::getResidueRange(
const std::string& chainname,
unsigned& res_start,
unsigned& res_end, std::string& errmsg )
const {
161 bool inres=
false, foundchain=
false;
162 for(
unsigned i=0;i<size();++i){
163 if( chain[i]==chainname ){
165 if(foundchain) errmsg=
"found second start of chain named " + chainname;
166 res_start=residue[i];
168 inres=
true; foundchain=
true;
169 }
else if( inres && chain[i]!=chainname ){
171 res_end=residue[i-1];
174 if(inres) res_end=residue[size()-1];
177 void PDB::getAtomRange(
const std::string& chainname,
AtomNumber& a_start,
AtomNumber& a_end, std::string& errmsg )
const {
178 bool inres=
false, foundchain=
false;
179 for(
unsigned i=0;i<size();++i){
180 if( chain[i]==chainname ){
182 if(foundchain) errmsg=
"found second start of chain named " + chainname;
185 inres=
true; foundchain=
true;
186 }
else if( inres && chain[i]!=chainname ){
191 if(inres) a_end=numbers[size()-1];
194 bool PDB::getBackbone(
const unsigned& resnumber,
const std::vector<std::string>& backnames, std::vector<AtomNumber>& backnumbers )
const {
195 if( backnames.size()!=backnumbers.size() ) plumed_merror(
"mismatch between number of backbone names and number of backbone numbers");
196 bool foundresidue=
false; std::vector<bool> found( backnames.size(), false );
197 for(
unsigned i=0;i<size();++i){
198 if( residue[i]==resnumber ){
200 for(
unsigned bno=0;bno<backnames.size();++bno){
201 if( backnames[bno]==atomsymb[i] ){ backnumbers[bno]=numbers[i]; found[bno]=
true; }
205 if(!foundresidue){
return false; }
206 for(
unsigned i=0;i<found.size();++i){
207 if( !found[i] )
return false;
212 std::string PDB::getChainID(
const unsigned& resnumber)
const {
213 for(
unsigned i=0;i<size();++i){
214 if(resnumber==residue[i])
return chain[i];
216 plumed_merror(
"Not enough residues in pdb input file");
219 void PDB::renameAtoms(
const std::string& old_name,
const std::string& new_name ){
220 for(
unsigned i=0;i<size();++i){
221 if( atomsymb[i]==old_name ) atomsymb[i]=new_name;
227 for(
unsigned i=0;i<pdb.
positions.size();i++){
Simple class to store the index of an atom.
Class implementing fixed size vectors of doubles.
Class containing the log stream.
std::vector< Vector > positions
Log & operator<<(Log &ostr, const PDB &pdb)
std::vector< AtomNumber > numbers