ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/openbabel/oopseformat.cpp
(Generate patch)

Comparing trunk/src/openbabel/oopseformat.cpp (file contents):
Revision 741 by tim, Wed Nov 16 19:42:11 2005 UTC vs.
Revision 1080 by gezelter, Thu Oct 19 19:51:34 2006 UTC

# Line 13 | Line 13 | GNU General Public License for more details.
13   GNU General Public License for more details.
14   ***********************************************************************/
15  
16 < #include <set>
17 <
18 < #include "mol.hpp"
19 < #include "obconversion.hpp"
20 < #include "obmolecformat.hpp"
21 <
22 < #ifdef HAVE_SSTREAM
23 < #include <sstream>
24 < #else
25 < #include <strstream>
26 < #endif
27 <
28 < using namespace std;
16 > #include "oopseformat.hpp"
17 > #include <fstream>
18   namespace OpenBabel
19   {
20  
32 – class OOPSEFormat : public OBMoleculeFormat
33 – {
34 – public:
35 –  //Register this format type ID
36 –  OOPSEFormat()
37 –  {
38 –    OBConversion::RegisterFormat("in", this, "chemical/x-in");
39 –  }
40 –
41 –  virtual const char* Description() //required
42 –  {
43 –    return
44 –      "oopse cartesian coordinates format\n";
45 –  };
46 –
47 –  virtual const char* SpecificationURL()
48 –  { return "";}; //optional
49 –
50 –  virtual const char* GetMIMEType()
51 –  { return "chemical/x-in"; };
52 –
53 –   virtual unsigned int Flags() { return WRITEONEONLY;}
54 –
55 –  //*** This section identical for most OBMol conversions ***
56 –  ////////////////////////////////////////////////////
57 –  /// The "API" interface functions
58 –  //virtual bool ReadMolecule(OBBase* pOb, OBConversion* pConv);
59 –  virtual bool WriteMolecule(OBBase* pOb, OBConversion* pConv);
60 –
61 –  private:
62 –      bool AreSameFragments(OBMol& mol, vector<int>& frag1, vector<int>& frag2);
63 –      void findAngles(OBMol& mol);
64 –      OBMol* createMolFromFragment(OBMol& mol, vector<int>& fragment);
65 –      void WriteMDFile(vector<OBMol*> mols, vector<int> numMols);
66 –      void WriteINFile(OBMol& mol, ostream& ofs, vector<int>& indices);
67 – };
68 –
69 –
70 – //Make an instance of the format class
71 – OOPSEFormat theOOPSEFormat;
72 –
73 – ////////////////////////////////////////////////////////////////
74 –
21   bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConversion* pConv)
22   {
23      OBMol* pmol = dynamic_cast<OBMol*>(pOb);
# Line 120 | Line 66 | bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConvers
66          numMols.push_back((*i).size());
67      }
68  
69 <    WriteMDFile(mdMols, numMols);    
69 >    string OutputFileName = pConv->GetInFilename();
70 >    unsigned int pos = OutputFileName.rfind(".");
71 >    if(pos==string::npos)
72 >        OutputFileName += ".md";
73 >    else
74 >        OutputFileName = OutputFileName.substr(0, pos) + ".md";      
75 >    ofstream ofs(OutputFileName.c_str());
76 >    if(!ofs)
77 >    {
78 >         cerr << "Cannot write to " << OutputFileName <<endl;
79 >         return false;
80 >    }
81 >                
82 >    WriteMDFile(mdMols, numMols, ofs, *pmol, indices);
83  
84      for(vector<OBMol*>::iterator  i = mdMols.begin(); i != mdMols.end(); ++i)
85      {
# Line 128 | Line 87 | bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConvers
87      }
88  
89      //    
131 –    WriteINFile(*pmol, *pConv->GetOutStream(), indices);
90  
133 –
134 –
91      return(true);
92   }
93  
# Line 142 | Line 98 | bool OOPSEFormat::AreSameFragments(OBMol& mol, vector<
98          return false;
99      }
100  
101 <    //exact graph matching is a NP complete problem,
101 >    //exact graph matching is a NP complete problem
102 >    /** @todo using sparse matrix to store the connectivities*/
103      for (unsigned int i =0 ; i < frag1.size(); ++i)
104      {
105 <        if (strcmp(mol.GetAtom(frag1[i])->GetType(), mol.GetAtom(frag2[i])->GetType()) )
105 >        OBAtom* atom1 = mol.GetAtom(frag1[i]);
106 >        OBAtom* atom2 = mol.GetAtom(frag2[i]);
107 >        
108 >        if (atom1->GetAtomicNum() != atom2->GetAtomicNum())
109          {
110              return false;
111          }
# Line 161 | Line 121 | struct SameAngle
121    }
122   };
123  
124 < void OOPSEFormat::findAngles(OBMol& mol)
165 < {
124 >  void OOPSEFormat::findAngles(OBMol& mol) {
125      /*
126      //if already has data return
127      if(mol.HasData(OBGenericDataType::AngleData))
# Line 201 | Line 160 | void OOPSEFormat::findAngles(OBMol& mol)
160          angles->SetData(angle);
161      }
162      */
163 < }
164 < OBMol* OOPSEFormat::createMolFromFragment(OBMol& mol, vector<int>& fragment)
165 < {
163 >  }
164 >
165 >  OBMol* OOPSEFormat::createMolFromFragment(OBMol& mol, vector<int>& fragment) {
166 >    
167      OBMol* newMol = new OBMol();
168      newMol->ReserveAtoms(fragment.size());
169      newMol->BeginModify();
170 <    for(vector<int>::iterator i = fragment.begin(); i != fragment.end(); ++i)
171 <    {
172 <        OBAtom* newAtom = newMol->NewAtom();
213 <        *newAtom = *mol.GetAtom(*i);
170 >    for(vector<int>::iterator i = fragment.begin(); i != fragment.end(); ++i) {
171 >      OBAtom* newAtom = newMol->NewAtom();
172 >      *newAtom = *mol.GetAtom(*i);
173      }
174      newMol->EndModify();
175      newMol->ConnectTheDots();
176      findAngles(*newMol);
177      newMol->FindTorsions();
178      return newMol;
179 < }
180 < void OOPSEFormat::WriteMDFile(vector<OBMol*> mols, vector<int> numMols)
181 < {
182 <    std::string identLevel1("\t");
183 <    std::string identLevel2("\t\t");
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("    ");
184      std::string molPrefix("MolName");
185 <    ostream& ofs = std::cout;
185 >    unsigned int i;
186      const int BUFFLEN = 1024;
187      char buffer[BUFFLEN];
188 +    string str, str1;
189      
190 <    for(unsigned int i = 0; i < mols.size(); ++i)
191 <    {
192 <        OBMol* pmol = mols[i];
193 <        map<OBAtom*, int> atomMap;
194 <        ofs << "molecule {\n";
195 <        sprintf(buffer, "%d", i);
236 <        ofs << identLevel1 << "name = " << "\"" << molPrefix << buffer << "\"" << ";\n";
190 >    
191 >    os << "<OOPSE version=4>" << endl;
192 >    os << "  <MetaData>" << endl << endl;
193 >    
194 >    for(i = 0; i < mols.size(); ++i) {
195 >      OBMol* pmol = mols[i];
196  
197 <        
198 <        //atom
199 <        ofs << identLevel1 << "nAtoms = " << pmol->NumAtoms() << ";\n";
200 <        int ai = 0;
201 <        FOR_ATOMS_OF_MOL(atom, *pmol ) {
243 <            ofs << identLevel1 << "atom[" << ai << "] {\n";
244 <            ofs << identLevel2 << "type = " << "\"" << atom->GetType() << "\"" << ";\n";
245 <            ofs << identLevel2 << "position(" << atom->GetX() << ", " << atom->GetY() << ", " << atom->GetZ() << ");\n";
246 <            ofs << identLevel1 << "}\n";
247 <            atomMap[&(*atom)] = ai++;
248 <        }        
197 >      pmol->ConnectTheDots();
198 >      pmol->PerceiveBondOrders();
199 >      //pmol->FindSSSR();
200 >      //pmol->SetAromaticPerceived();
201 >      //pmol->Kekulize();
202  
203 <        //bond
204 <        ofs << identLevel1 << "nBonds = " << pmol->NumBonds() << ";\n";
205 <        int bi = 0;
206 <        FOR_BONDS_OF_MOL(bond, *pmol ) {
207 <            ofs << identLevel1 << "bond[" << bi++ << "] {\n";
208 <            ofs << identLevel2 << "member(" << atomMap[bond->GetBeginAtom()] <<  ", " << atomMap[bond->GetEndAtom()] << ");\n";
209 <            ofs << identLevel1 << "}\n";
210 <        }  
211 <        /*
212 <        //bend
213 <        OBGenericData* pGenericData = pmol->GetData(OBGenericDataType::AngleData);
214 <        OBAngleData* pAngleData = dynamic_cast<OBAngleData*>(pGenericData);
215 <        vector<OBAngle> angles = pAngleData->GetData();
216 <
217 <        ofs << identLevel1 << "nBends = " << angles.size() << ";\n";        
218 <        int bendIndex = 0;
219 <        for (vector<OBAngle>::iterator ti = angles.begin(); ti != angles.end(); ++ti)
203 >      map<OBAtom*, int> atomMap;
204 >      os << "molecule {\n";
205 >      sprintf(buffer, "%d", i);
206 >      os << indentLevel1 << "name = " << "\"" << molPrefix << buffer << "\"" << ";\n";
207 >          
208 >      //atom
209 >      int ai = 0;
210 >      FOR_ATOMS_OF_MOL(atom, *pmol ) {
211 >        str = atom->GetType();
212 >        ttab.SetFromType("INT");
213 >        ttab.SetToType("INT");
214 >        ttab.Translate(str1,str);
215 >        os << indentLevel1 << "atom[" << ai << "] {\n";
216 >        // os << indentLevel2 << "type = " << "\"" << etab.GetSymbol(atom->GetAtomicNum()) << "\"" << ";\n";
217 >        os << indentLevel2 << "type = " << "\"" << str1 << "\"" << ";\n";
218 >        os << indentLevel1 << "}\n";
219 >        atomMap[&(*atom)] = ai++;
220 >      }        
221 >      os << "\n";
222 >      
223 >      //bond
224 >      FOR_BONDS_OF_MOL(bond, *pmol ) {
225 >        os << indentLevel1 << "bond {\n";
226 >        os << indentLevel2 << "members(" << atomMap[bond->GetBeginAtom()] <<  ", " << atomMap[bond->GetEndAtom()] << ");\n";
227 >        os << indentLevel1 << "}\n";
228 >      }  
229 >      /*
230 >      //bend
231 >      OBGenericData* pGenericData = pmol->GetData(OBGenericDataType::AngleData);
232 >      OBAngleData* pAngleData = dynamic_cast<OBAngleData*>(pGenericData);
233 >      vector<OBAngle> angles = pAngleData->GetData();
234 >      
235 >      os << indentLevel1 << "nBends = " << angles.size() << ";\n";        
236 >      int bendIndex = 0;
237 >      for (vector<OBAngle>::iterator ti = angles.begin(); ti != angles.end(); ++ti)
238          {
239 <            triple<OBAtom*, OBAtom*, OBAtom*> bendAtoms = ti->getAtoms();
240 <            ofs << identLevel1 << "bend[" << bendIndex++ << "] {\n";
241 <            ofs << identLevel2 << "member(" << atomMap[bendAtoms.first] <<  ", " << atomMap[bendAtoms.second] << atomMap[bendAtoms.third] <<");\n";
242 <            ofs << identLevel1 << "}\n";            
239 >          triple<OBAtom*, OBAtom*, OBAtom*> bendAtoms = ti->getAtoms();
240 >          os << indentLevel1 << "bend[" << bendIndex++ << "] {\n";
241 >          os << indentLevel2 << "member(" << atomMap[bendAtoms.first] <<  ", " << atomMap[bendAtoms.second] << atomMap[bendAtoms.third] <<");\n";
242 >          os << indentLevel1 << "}\n";            
243          }
244 <        
245 <        //torsion
246 <        pGenericData = pmol->GetData(OBGenericDataType::TorsionData);
247 <        OBTorsionData* pTorsionData = dynamic_cast<OBTorsionData*>(pGenericData);
248 <        vector<OBTorsion> torsions = pTorsionData->GetData();
249 <        vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > torsionArray;
250 <        for (vector<OBTorsion>::iterator ti = torsions.begin(); ti != torsions.end(); ++ti)
244 >      
245 >      //torsion
246 >      pGenericData = pmol->GetData(OBGenericDataType::TorsionData);
247 >      OBTorsionData* pTorsionData = dynamic_cast<OBTorsionData*>(pGenericData);
248 >      vector<OBTorsion> torsions = pTorsionData->GetData();
249 >      vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > torsionArray;
250 >      for (vector<OBTorsion>::iterator ti = torsions.begin(); ti != torsions.end(); ++ti)
251          {
252 <            vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > tmpTorsions = ti->getTorsions();
253 <            torsionArray.insert(torsionArray.end(), tmpTorsions.begin(), tmpTorsions.end());            
252 >          vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > tmpTorsions = ti->getTorsions();
253 >          torsionArray.insert(torsionArray.end(), tmpTorsions.begin(), tmpTorsions.end());            
254          }
255 <
256 <        ofs << identLevel1 << "nTorsions = " << torsionArray.size() << ";\n";
257 <        int torsionIndex = 0;
258 <        for (vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> >::iterator ti = torsionArray.begin(); ti != torsionArray.end(); ++ti)
255 >      
256 >      os << indentLevel1 << "nTorsions = " << torsionArray.size() << ";\n";
257 >      int torsionIndex = 0;
258 >      for (vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> >::iterator ti = torsionArray.begin(); ti != torsionArray.end(); ++ti)
259          {
260 <            ofs << identLevel1 << "torsion[" << torsionIndex++ << "] {\n";
261 <            ofs << identLevel2 << "member(" << atomMap[ti->first] <<  ", " << atomMap[ti->second] <<", " << atomMap[ti->third] <<", " << atomMap[ti->forth] << ");\n";
262 <            ofs << identLevel1 << "}\n";          
260 >          os << indentLevel1 << "torsion[" << torsionIndex++ << "] {\n";
261 >          os << indentLevel2 << "member(" << atomMap[ti->first] <<  ", " << atomMap[ti->second] <<", " << atomMap[ti->third] <<", " << atomMap[ti->forth] << ");\n";
262 >          os << indentLevel1 << "}\n";          
263          }
264 <        */
265 <        ofs << "}\n";
264 >      */
265 >      os << "}" << endl;
266 >      os << endl;
267 >      
268      }
269 +    
270 +    os << endl;
271  
297 –    ofs << "nComponents = " << mols.size() << ";\n";
272      
273 <    for(unsigned int i =0; i < mols.size(); ++i)
274 <    {
275 <        ofs << "component{\n";
276 <        sprintf(buffer, "%d", i);
277 <        ofs << "type = " << molPrefix << buffer << ";\n";
278 <        ofs << "nMol = " << numMols[i]<< ";\n";
305 <        ofs << "}\n";
273 >    for(i=0; i < mols.size(); ++i) {      
274 >      os << "component{" << endl;
275 >      sprintf(buffer, "%d", i);
276 >      os << indentLevel1 << "type = " << molPrefix << buffer << ";" << endl;
277 >      os << indentLevel1 << "nMol = " << numMols[i]<< ";" << endl;
278 >      os << "}" << endl;
279      }
307 – }
308 – void OOPSEFormat::WriteINFile(OBMol& mol, ostream& ofs, vector<int>& indices)
309 – {
310 –    unsigned int i;
311 –    char buffer[BUFF_SIZE];
280      
281 <    sprintf(buffer,"%d", mol.NumAtoms());
282 <    ofs << buffer << endl;
283 <    //sprintf(buffer,"0;%f, 0, 0; 0, %15.7f, 0; 0, 0, %15.7f;", boxx, boxy, boxz);
284 <    sprintf(buffer, "0;%f, 0, 0; 0, %15.7f, 0; 0, 0, %15.7f;", 100.0, 100.0, 100.0);
285 <    ofs << buffer << endl;
281 >    os << "  </MetaData>" << endl;
282 >    os << "  <Snapshot>" << endl;
283 >    os << "    <FrameData>" << endl;
284 >    
285 >    sprintf(buffer, "        Time: %.10g", 0.0);
286 >    
287 >    os << buffer << endl;
288  
289 +    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);
290 +
291 +    os << buffer << endl;
292 +    os << "    </FrameData>" << endl;
293 +    os << "    <StuntDoubles>" << endl;
294 +
295      OBAtom *atom;
320 –    string str,str1;
296      
297 <    for(vector<int>::iterator i = indices.begin();i != indices.end(); ++i)
298 <    {
299 <        atom = mol.GetAtom(*i);
300 <        sprintf(buffer,"%15s%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f%15.5f",
326 <                atom->GetType(),
327 <                atom->GetX(), atom->GetY(), atom->GetZ(),
328 <                0.0, 0.0, 0.0,
329 <                0.0, 0.0, 0.0, 0.0,
330 <                0.0, 0.0, 0.0
331 <                );
332 <        ofs << buffer << endl;
297 >    for(vector<int>::iterator i = indices.begin();i != indices.end(); ++i) {    
298 >      atom = mol.GetAtom(*i);
299 >      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);
300 >      os << buffer << endl;
301      }
302 <    
302 >    os << "    </StuntDoubles>" << endl;
303 >    os << "  </Snapshot>" << endl;
304 >    os << "</OOPSE>" << endl;
305   }
306 <
306 >  
307   } //namespace OpenBabel
308  

Diff Legend

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