ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/applications/atom2md/oopseformat.cpp
(Generate patch)

Comparing trunk/OOPSE-4/src/applications/atom2md/oopseformat.cpp (file contents):
Revision 3424 by gezelter, Tue Jul 1 13:28:23 2008 UTC vs.
Revision 3429 by gezelter, Thu Jul 3 12:45:30 2008 UTC

# Line 53 | Line 53 | namespace OpenBabel
53        return NOTREADABLE | WRITEONEONLY;
54      }
55      
56 –    //*** This section identical for most OBMol conversions ***
57 –    ////////////////////////////////////////////////////
58 –    /// The "API" interface functions
56      virtual bool WriteMolecule(OBBase* pOb, OBConversion* pConv);
57  
58    private:
59      bool AreSameFragments(OBMol& mol, vector<int>& frag1, vector<int>& frag2);
63 –    //void findAngles(OBMol& mol);
60      OBMol* createMolFromFragment(OBMol& mol, vector<int>& fragment);
61 <    void WriteMDFile(vector<OBMol*> mols, vector<int> numMols, ostream& os, OBMol& mol, vector<int>& indices);
61 >    void WriteMDFile(vector<OBMol*> mols, vector<int> numMols, ostream& os,
62 >                     OBMol& mol, vector<int>& indices);
63    };
64 <
64 >  
65    //Make an instance of the format class
66    OOPSEFormat theOOPSEFormat;
67 <  
68 <  /////////////////////////////////////////////////////////////////
72 <
73 <
74 <  bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConversion* pConv)
75 <  {
67 >
68 >  bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConversion* pConv) {
69      OBMol* pmol = dynamic_cast<OBMol*>(pOb);
70      if(pmol==NULL)
71        return false;
# Line 84 | Line 77 | namespace OpenBabel
77      vector<vector<int> > molecules;
78      vector<int> indices;
79      for(int i =0; i < used.size(); ++i) {
80 <      if (used[i])
81 <        {
82 <          continue;
90 <        }
80 >      if (used[i])
81 >        continue;
82 >
83        used[i] = true;
84        vector<int> sameMolTypes;
85        sameMolTypes.push_back(i);
86        indices.insert(indices.end(), fragmentLists[i].begin(),
87                       fragmentLists[i].end());
88 <      for (int j = i + 1;j < used.size(); ++j)
89 <        {
90 <          if (used[j])
91 <            {
92 <              continue;
93 <            }
94 <          
95 <          if (AreSameFragments(*pmol, fragmentLists[i], fragmentLists[j]))
96 <            {
105 <              sameMolTypes.push_back(j);
106 <              indices.insert(indices.end(), fragmentLists[j].begin(),
107 <                             fragmentLists[j].end());
108 <              used[j]=true;
109 <            }
88 >      for (int j = i + 1;j < used.size(); ++j) {
89 >        if (used[j])
90 >          continue;
91 >                
92 >        if (AreSameFragments(*pmol, fragmentLists[i], fragmentLists[j])) {
93 >          sameMolTypes.push_back(j);
94 >          indices.insert(indices.end(), fragmentLists[j].begin(),
95 >                         fragmentLists[j].end());
96 >          used[j]=true;
97          }
98 <      molecules.push_back(sameMolTypes);
99 <      
98 >      }
99 >      molecules.push_back(sameMolTypes);      
100      }
101      
115 –    //
102      vector<OBMol*> mdMols;    
103      vector<int> numMols;
104      for(vector<vector<int> >::iterator  i = molecules.begin();
105 <        i != molecules.end(); ++i)
106 <      {
107 <        mdMols.push_back(createMolFromFragment(*pmol,
108 <                                               fragmentLists[i->front()]));
109 <        numMols.push_back((*i).size());
110 <      }
105 >        i != molecules.end(); ++i) {
106 >      
107 >      mdMols.push_back(createMolFromFragment(*pmol,
108 >                                             fragmentLists[i->front()]));
109 >      numMols.push_back((*i).size());
110 >    }
111      
112      string OutputFileName = pConv->GetInFilename();
113      size_t pos = OutputFileName.rfind(".");
# Line 129 | Line 115 | namespace OpenBabel
115        OutputFileName = OutputFileName.substr(0, pos) + ".md";      
116      else
117        OutputFileName += ".md";
118 <
118 >    
119      ofstream ofs(OutputFileName.c_str());
120 <    if(!ofs)
135 <      {
120 >    if(!ofs) {
121          cerr << "Cannot write to " << OutputFileName <<endl;
122          return false;
123 <      }
123 >    }
124      
125      WriteMDFile(mdMols, numMols, ofs, *pmol, indices);
126      
127 <    for(vector<OBMol*>::iterator  i = mdMols.begin(); i != mdMols.end(); ++i)
128 <      {
129 <        delete *i;
145 <      }
127 >    for(vector<OBMol*>::iterator  i = mdMols.begin(); i != mdMols.end(); ++i) {
128 >      delete *i;
129 >    }
130      
147 –    //    
148 –    
131      return(true);
132    }
133    
134    bool OOPSEFormat::AreSameFragments(OBMol& mol, vector<int>& frag1,
135 <                                     vector<int>& frag2)
154 <  {
135 >                                     vector<int>& frag2) {
136      if (frag1.size() != frag2.size())
137 <      {
157 <        return false;
158 <      }
137 >      return false;
138      
139 <    //exact graph matching is a NP complete problem
140 <    /** @todo using sparse matrix to store the connectivities*/
141 <    for (unsigned int i =0 ; i < frag1.size(); ++i)
142 <      {
139 >    // Exact graph matching is an NP complete problem.
140 >    // This just matches all of the atom atomic numbers and may falsely
141 >    // detect identical fragments which aren't really identical.
142 >    // @todo using sparse matrix to store the connectivities
143 >
144 >    for (unsigned int i =0 ; i < frag1.size(); ++i) {
145          OBAtom* atom1 = mol.GetAtom(frag1[i]);
146          OBAtom* atom2 = mol.GetAtom(frag2[i]);
147          
148 <        if (atom1->GetAtomicNum() != atom2->GetAtomicNum())
149 <          {
150 <            return false;
151 <          }
171 <      }
148 >        if (atom1->GetAtomicNum() != atom2->GetAtomicNum())
149 >          return false;
150 >        
151 >    }
152      return true;
153    }
154    
155 <  struct SameAngle
176 <  {
155 >  struct SameAngle {
156      bool operator()(const triple<OBAtom*,OBAtom*,OBAtom*> t1,
157 <                    const triple<OBAtom*,OBAtom*,OBAtom*> t2) const
179 <    {
157 >                    const triple<OBAtom*,OBAtom*,OBAtom*> t2) const {
158        return (t1.second == t2.second) && ( (t1.first == t2.first && t1.third == t2.third) || (t1.first == t2.third && t1.third == t2.first));
159      }
160    };
161  
184 –  /*  
185 –    void OOPSEFormat::findAngles(OBMol& mol) {
162  
163 <    //if already has data return
164 <    if(mol.HasData(OBGenericDataType::AngleData))
189 <    return;
190 <    
191 <    vector<OBEdgeBase*>::iterator bi1,bi2;
192 <    OBBond* bond;
193 <    OBAtom *a,*b,*c;
163 >  OBMol* OOPSEFormat::createMolFromFragment(OBMol& mol,
164 >                                            vector<int>& fragment) {
165      
195 –    set<triple<OBAtom*,OBAtom*,OBAtom*>, SameAngle> uniqueAngles;
196 –    //loop through all bonds generating torsions
197 –    for(bond = mol.BeginBond(bi1);bond;bond = mol.NextBond(bi1))
198 –    {
199 –        b = bond->GetBeginAtom();
200 –        c = bond->GetEndAtom();
201 –        if(b->IsHydrogen())
202 –            continue;
203 –
204 –        for(a = b->BeginNbrAtom(bi2);a;a = b->NextNbrAtom(bi2))
205 –        {
206 –            if(a == c)
207 –                continue;          
208 –            
209 –            uniqueAngles.insert(triple<OBAtom*,OBAtom*,OBAtom*>(a, b, c));
210 –        }
211 –    }
212 –
213 –    //get new data and attach it to molecule
214 –    OBAngleData *angles = new OBAngleData;
215 –    mol.SetData(angles);
216 –    set<triple<OBAtom*,OBAtom*,OBAtom*>, SameAngle>::iterator i;
217 –
218 –    for (i = uniqueAngles.begin(); i != uniqueAngles.end(); ++i) {
219 –        OBAngle angle;
220 –        angle.SetAtoms(i->first, i->second, i->second);
221 –        angles->SetData(angle);
222 –    }
223 –
224 –  }
225 –  */
226 –  
227 –  OBMol* OOPSEFormat::createMolFromFragment(OBMol& mol, vector<int>& fragment) {
228 –    
166      OBMol* newMol = new OBMol();
230 –    OBChainsParser* chainParser = new OBChainsParser();
167      newMol->ReserveAtoms(fragment.size());
168      newMol->BeginModify();
169      for(vector<int>::iterator i = fragment.begin(); i != fragment.end(); ++i) {
170        OBAtom* newAtom = newMol->NewAtom();
171        *newAtom = *mol.GetAtom(*i);
172      }
173 +
174      newMol->EndModify();
175      newMol->ConnectTheDots();
176      newMol->PerceiveBondOrders();
240 –    newMol->FindAngles();
241 –    newMol->FindTorsions();
242 –    //chainParser->PerceiveChains(*newMol, false);
177  
178      return newMol;
179    }
180    
181 <  void OOPSEFormat::WriteMDFile(vector<OBMol*> mols, vector<int> numMols, ostream& os, OBMol& mol, vector<int>& indices) {
182 <    std::string indentLevel1("  ");
183 <    std::string indentLevel2("    ");
181 >  void OOPSEFormat::WriteMDFile(vector<OBMol*> mols, vector<int> numMols,
182 >                                ostream& os, OBMol& mol,
183 >                                vector<int>& indices) {
184 >    
185      std::string molPrefix("MolName");
186      std::string resName;
187      unsigned int i;
# Line 254 | Line 189 | namespace OpenBabel
189      char buffer[BUFFLEN];
190      string str, str1;
191      OBAtom *a, *b, *c, *d;    
192 <    OBChainsParser* chainParser = new OBChainsParser();
192 >    bool molIsWater = false;
193      OBResidue *r;
194      int resKey;
195      char type_name[10];
196      char *element_name;
197      int res_num;
198 +    OBChainsParser* chainParser = new OBChainsParser();
199 +
200      
264 –    std::cerr << "yo\n";
201      os << "<OOPSE version=4>" << endl;
202      os << "  <MetaData>" << endl << endl;
203      
204      for(i = 0; i < mols.size(); ++i) {
205        OBMol* pmol = mols[i];
270 –      std::cerr << "yo1\n";
271 –      pmol->ConnectTheDots();
272 –      pmol->PerceiveBondOrders();
273 –      chainParser->PerceiveChains(*pmol, false);
274 –      //pmol->FindSSSR();
275 –      //pmol->SetAromaticPerceived();
276 –      //pmol->Kekulize();
277 –      
206        map<OBAtom*, int> atomMap;
279 –      os << "molecule {\n";
280 –      sprintf(buffer, "%d", i);
281 –      os << indentLevel1 << "name = " << "\"" << molPrefix << buffer << "\"" << ";\n";
282 –      
207  
208 <      //atom
209 <      
210 <      int ai = 0;
211 <      FOR_RESIDUES_OF_MOL(res, *pmol) {
212 <
213 <        std::cerr << "found residue" << res->GetName() << "\n";
208 >      chainParser->PerceiveChains(*pmol, false);
209 >      molIsWater = false;
210 >      FOR_RESIDUES_OF_MOL(residue, *pmol) {
211 >        std::cerr << "residue = " << residue->GetName() << "\n";
212 >        if (residue->GetName().compare("HOH") == 0) {
213 >          molIsWater = true;
214 >        }
215        }
216        
217 < //       FOR_RESIDUES_OF_MOL(res, *pmol) {
218 <        
219 < //      resName = res->GetName();
220 < //      resKey = res->GetResKey();
217 >      if (molIsWater) {
218 >        // water include files define all of the known water types
219 >        os << "#include \"water.md\";\n";
220 >        pmol->SetTitle("HOH");
221 >      } else {
222  
223 <        
224 <
225 < //      std::cerr << "found residue " << resName << "\n";
226 < //      std::cerr << "  num = " << res->GetNum() << "\n";
227 < //      std::cerr << "  numAtoms = " << res->GetNumAtoms() << "\n";
228 < //      std::cerr << "  num = " << res->GetNum() << "\n";
229 < //      std::cerr << "  chain = " << res->GetChain() << "\n";
230 < //      std::cerr << "  chainnum = " << res->GetChainNum() << "\n";
231 < //      std::cerr << "  idx = " << res->GetIdx() << "\n";
232 < //      std::cerr << "  key = " << res->GetResKey() << "\n";
223 >        os << "molecule {\n";
224 >        sprintf(buffer, "%d", i);
225 >        os << "  name = \"" << molPrefix << buffer << "\";\n";
226 >        
227 >        int ai = 0;
228 >        FOR_ATOMS_OF_MOL(atom, *pmol ) {
229 >          str = atom->GetType();
230 >          r = atom->GetResidue();
231 >          
232 >          if (r == NULL)
233 >            resName = "NULL";
234 >          else
235 >            resName = r->GetName();
236 >          
237 >          if (resName.compare("NULL") ==0 ||
238 >              resName.compare("LIG") == 0 ||
239 >              resName.compare("UNK") == 0) {
240 >            // Either couldn't find a residue at all or couldn't find a
241 >            // reasonable residue name to use.  We'll punt and use
242 >            // OpenBabel's internal atom typing:
243 >            ttab.SetFromType("INT");
244 >            ttab.SetToType("INT");
245 >            ttab.Translate(str1, str);
246 >          } else {                        
247 >            
248 >            // If we know what residue we've got, the specific atom name can
249 >            // be used to help specify partial charges.
250 >            
251 >            str = r->GetAtomID(&*atom);
252 >            size_t startpos = str.find_first_not_of(" ");
253 >            size_t endpos = str.find_last_not_of(" ");
254 >            str = str.substr( startpos, endpos-startpos+1 );
255 >            str1 = resName + "-" + str;
256 >          }      
257 >          os << "  atom[" << ai << "] { ";
258 >          os << "type = " << "\"" << str1 << "\"" << "; ";
259 >          os << "}\n";
260 >          atomMap[&(*atom)] = ai++;
261 >        }        
262 >        os << "\n";
263 >        
264 >        //bond
265 >        
266 >        int b1, b2;
267 >        FOR_BONDS_OF_MOL(bond, *pmol ) {
268 >          b1 = atomMap[bond->GetBeginAtom()];
269 >          b2 = atomMap[bond->GetEndAtom()];
270 >          
271 >          os << "  bond { ";
272 >          
273 >          if (b1 < b2)
274 >            os << "members(" << b1 <<  ", " << b2 << "); ";
275 >          else
276 >            os << "members(" << b2 <<  ", " << b1 << "); ";
277 >          
278 >          os << "}\n";
279 >        }
280 >        
281 >        os << "\n";
282  
283 <
284 < //      FOR_ATOMS_OF_RESIDUE(atom, &*res) {
285 < //        str = atom->GetType();
311 < //        ttab.SetFromType("INT");
312 < //        ttab.SetToType("INT");
313 < //        ttab.Translate(str1,str);
314 < //        os << indentLevel1 << "atom[" << ai << "] {\n";
315 < //        os << indentLevel2 << "type = " << "\"" << resName << "-" << str1 << "\"" << ";\n";
316 < //        os << indentLevel1 << "}\n";
317 < //        atomMap[&(*atom)] = ai++;
318 < //      }        
319 < //      os << "\n";
320 <        
321 < //       }      
322 <
323 <      ai = 0;
324 <      FOR_ATOMS_OF_MOL(atom, *pmol ) {
325 <        str = atom->GetType();
326 <        r = atom->GetResidue();
327 <        
328 <        if (r == NULL)
329 <          resName = "NULL";
330 <        else
331 <          resName = r->GetName();
332 <                
333 <        if (resName.compare("NULL") ==0 || resName.compare("LIG") == 0 || resName.compare("UNK") == 0) {
334 <          // Either couldn't find a residue at all or couldn't find a reasonable
335 <          // residue name to use.  We'll punt and use OpenBabel's internal atom typing:
336 <          ttab.SetFromType("INT");
337 <          ttab.SetToType("INT");
338 <          ttab.Translate(str1,str);
339 <        } else {
340 <                
341 <
342 <          if (resName.compare("HOH") == 0) {
343 <            // HOH is a special residue name for water, and the standard atom types
344 <            // are OW and HW, so just append W to the string for the atom type:
345 <            ttab.SetFromType("INT");
346 <            ttab.SetToType("XYZ");
347 <            ttab.Translate(str1,str);
348 <            str1 += "W";
349 <          } else {
350 <
351 <            std::cerr << "found residue " << resName << "\n";
352 <
353 <            // If we know what residue we've got, the specific atom name can
354 <            // be used to help specify partial charges.
355 <
356 <            str = r->GetAtomID(&*atom);
357 <            size_t startpos = str.find_first_not_of(" ");
358 <            size_t endpos = str.find_last_not_of(" ");
359 <            str = str.substr( startpos, endpos-startpos+1 );
360 <            str1 = resName + "-" + str;
361 <          }
362 <
363 <        }
364 <
365 < //         os << indentLevel1 << "atom[" << ai << "] {\n";
366 < //         os << indentLevel2 << "type = " << "\"" << str1 << "\"" << ";\n";
367 < //         os << indentLevel1 << "}\n";
368 <        os << "atom[" << ai << "] { ";
369 <        os << "type = " << "\"" << str1 << "\"" << "; ";
370 <        os << "}\n";
371 <        atomMap[&(*atom)] = ai++;
372 <      }        
373 <      os << "\n";
374 <      
375 <      //bond
376 <      FOR_BONDS_OF_MOL(bond, *pmol ) {
377 < //         os << indentLevel1 << "bond {\n";
378 < //         os << indentLevel2 << "members(" << atomMap[bond->GetBeginAtom()] <<  ", " << atomMap[bond->GetEndAtom()] << ");\n";
379 < //         os << indentLevel1 << "}\n";
380 <        os << "bond { ";
381 <        os << "members(" << atomMap[bond->GetBeginAtom()] <<  ", " << atomMap[bond->GetEndAtom()] << "); ";
382 <        os << "}\n";
383 <      }
384 <
385 <      FOR_ANGLES_OF_MOL(angle, *pmol) {
386 <
387 <        // OpenBabel's getAtoms returns the 3 atom pointer for the
388 <        // angle with the vertex first.  These need to be reordered
389 <        // for vertex-second ordering for OOPSE.
390 <
391 <         b = pmol->GetAtom((*angle)[0] + 1);
392 <         a = pmol->GetAtom((*angle)[1] + 1);
393 <         c = pmol->GetAtom((*angle)[2] + 1);
394 <
395 < //       os << indentLevel1 << "bend {\n";
396 < //       os << indentLevel2 << "members(" << atomMap[a] <<  ", " << atomMap[b] << ", " << atomMap[c] << ");\n";
397 < //       os << indentLevel1 << "}\n";
398 <         os << "bend { ";
399 <         os << "members(" << atomMap[a] <<  ", " << atomMap[b] << ", " << atomMap[c] << "); ";
400 <         os << "}\n";
401 <      }
402 <
403 <      FOR_TORSIONS_OF_MOL(torsion, *pmol) {
404 <
405 <        // OpenBabel's getAtoms returns the 3 atom pointer for the
406 <        // angle with the vertex first.  These need to be reordered
407 <        // for vertex-second ordering for OOPSE.
408 <
409 <         a = pmol->GetAtom((*torsion)[0] + 1);
410 <         b = pmol->GetAtom((*torsion)[1] + 1);
411 <         c = pmol->GetAtom((*torsion)[2] + 1);
412 <         d = pmol->GetAtom((*torsion)[3] + 1);
413 <
414 < //         os << indentLevel1 << "torsion {\n";
415 < //         os << indentLevel2 << "members(" << atomMap[a] <<  ", " << atomMap[b] << ", " << atomMap[c] << ", " << atomMap[d] << ");\n";
416 < //         os << indentLevel1 << "}\n";
417 <
418 <        os << "torsion { ";
419 <        os << "members(" << atomMap[a] <<  ", " << atomMap[b] << ", " << atomMap[c] << ", " << atomMap[d] << "); ";
420 <        os << "}\n";
421 <      }
422 <
423 <      /*
424 <      //bend
425 <      OBGenericData* pGenericData = pmol->GetData(OBGenericDataType::AngleData);
426 <      OBAngleData* pAngleData = dynamic_cast<OBAngleData*>(pGenericData);
427 <      vector<OBAngle> angles = pAngleData->GetData();
428 <      
429 <      os << indentLevel1 << "nBends = " << angles.size() << ";\n";        
430 <      int bendIndex = 0;
431 <      for (vector<OBAngle>::iterator ti = angles.begin(); ti != angles.end(); ++ti)
432 <      {
433 <      triple<OBAtom*, OBAtom*, OBAtom*> bendAtoms = ti->getAtoms();
434 <      os << indentLevel1 << "bend[" << bendIndex++ << "] {\n";
435 <          os << indentLevel2 << "member(" << atomMap[bendAtoms.first] <<  ", " << atomMap[bendAtoms.second] << atomMap[bendAtoms.third] <<");\n";
436 <          os << indentLevel1 << "}\n";            
437 <          }
283 >        std::vector<int> possibleInversion;
284 >        FOR_ATOMS_OF_MOL(atom, *pmol) {
285 >          possibleInversion.clear();
286            
287 <          //torsion
288 <          pGenericData = pmol->GetData(OBGenericDataType::TorsionData);
441 <          OBTorsionData* pTorsionData = dynamic_cast<OBTorsionData*>(pGenericData);
442 <          vector<OBTorsion> torsions = pTorsionData->GetData();
443 <          vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > torsionArray;
444 <          for (vector<OBTorsion>::iterator ti = torsions.begin(); ti != torsions.end(); ++ti)
445 <          {
446 <          vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > tmpTorsions = ti->getTorsions();
447 <          torsionArray.insert(torsionArray.end(), tmpTorsions.begin(), tmpTorsions.end());            
287 >          FOR_NBORS_OF_ATOM(nbor, &*atom) {
288 >            possibleInversion.push_back(atomMap[&(*nbor)]);          
289            }
290            
291 <          os << indentLevel1 << "nTorsions = " << torsionArray.size() << ";\n";
292 <          int torsionIndex = 0;
293 <          for (vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> >::iterator ti = torsionArray.begin(); ti != torsionArray.end(); ++ti)
294 <          {
295 <          os << indentLevel1 << "torsion[" << torsionIndex++ << "] {\n";
455 <          os << indentLevel2 << "member(" << atomMap[ti->first] <<  ", " << atomMap[ti->second] <<", " << atomMap[ti->third] <<", " << atomMap[ti->forth] << ");\n";
456 <          os << indentLevel1 << "}\n";          
291 >          if (possibleInversion.size() == 3) {
292 >
293 >            os << "  inversion { ";
294 >            os << "center(" << atomMap[&(*atom)] << "); ";
295 >            os << "}\n";
296            }
297 <      */
298 <      os << "}" << endl;
299 <      os << endl;
300 <      
297 >          
298 >        }
299 >        os << "}" << endl;
300 >        os << endl;
301 >      }
302      }
303      
304      os << endl;
305      
306      
307 <    for(i=0; i < mols.size(); ++i) {      
307 >    for(i=0; i < mols.size(); ++i) {
308 >      OBMol* pmol = mols[i];      
309        os << "component{" << endl;
310 <      sprintf(buffer, "%d", i);
311 <      os << indentLevel1 << "type = " << molPrefix << buffer << ";" << endl;
312 <      os << indentLevel1 << "nMol = " << numMols[i]<< ";" << endl;
310 >      if (std::string(pmol->GetTitle()).compare("HOH") == 0) {
311 >        os << "  type = " << "HOH" << ";" << endl;
312 >      } else {
313 >        sprintf(buffer, "%d", i);
314 >        os << "  type = " << molPrefix << buffer << ";" << endl;
315 >      }
316 >      os << "  nMol = " << numMols[i]<< ";" << endl;
317        os << "}" << endl;
318      }
319      
# Line 480 | Line 325 | namespace OpenBabel
325      
326      os << buffer << endl;
327      
328 +
329 +    // should really compute a bounding box here:
330      sprintf(buffer, "        Hmat: {{ %.10g, %.10g, %.10g }, { %.10g, %.10g, %.10g }, { %.10g, %.10g, %.10g }}", 100.0, 0.0, 0.0, 0.0, 100.0, 0.0, 0.0, 0.0, 100.0);
331      
332      os << buffer << endl;
# Line 488 | Line 335 | namespace OpenBabel
335      
336      OBAtom *atom;
337      
338 +    // still to do: intercept waters and recompute pvqj lines
339 +
340      for(vector<int>::iterator i = indices.begin();i != indices.end(); ++i) {    
341        atom = mol.GetAtom(*i);
342        sprintf(buffer, "%10d %7s %18.10g %18.10g %18.10g %13e %13e %13e", *i - 1, "pv", atom->GetX(), atom->GetY(), atom->GetZ(), 0.0, 0.0, 0.0);

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines