src/openbabel/oopseformat.cpp

/**********************************************************************
Copyright (C) 2000 by OpenEye Scientific Software, Inc.
Some portions Copyright (C) 2001-2005 by Geoffrey R. Hutchison
Some portions Copyright (C) 2004 by Chris Morley
 
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation version 2 of the License.
 
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.
***********************************************************************/

#include "oopseformat.hpp"
#include <fstream>
namespace OpenBabel
{

bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConversion* pConv)
{
    OBMol* pmol = dynamic_cast<OBMol*>(pOb);
    if(pmol==NULL)
        return false;
    
    vector<vector<int> > fragmentLists;
    pmol->ContigFragList(fragmentLists);
    OBBitVec unused;
    vector<bool> used(fragmentLists.size(), 0);
    vector<vector<int> > molecules;
    vector<int> indices;
    for(int i =0; i < used.size(); ++i) {
        if (used[i])
        {
            continue;
        }
        used[i] = true;
        vector<int> sameMolTypes;
        sameMolTypes.push_back(i);
        indices.insert(indices.end(), fragmentLists[i].begin(), fragmentLists[i].end());
        for (int j = i + 1;j < used.size(); ++j)
        {
            if (used[j])
            {
                continue;
            }
            
            if (AreSameFragments(*pmol, fragmentLists[i], fragmentLists[j]))
            {
                sameMolTypes.push_back(j);
                indices.insert(indices.end(), fragmentLists[j].begin(), fragmentLists[j].end());
                used[j]=true;
            }
        }
        molecules.push_back(sameMolTypes);
        
    }

    //
    vector<OBMol*> mdMols;    
    vector<int> numMols;
    for(vector<vector<int> >::iterator  i = molecules.begin(); i != molecules.end(); ++i) 
    {
        mdMols.push_back(createMolFromFragment(*pmol, fragmentLists[i->front()]));
        numMols.push_back((*i).size());
    }

    string OutputFileName = pConv->GetInFilename();
    unsigned int pos = OutputFileName.rfind(".");
    if(pos==string::npos)
        OutputFileName += ".md";
    else
        OutputFileName = OutputFileName.substr(0, pos) + ".md";       
    ofstream ofs(OutputFileName.c_str());
    if(!ofs)
    {
         cerr << "Cannot write to " << OutputFileName <<endl;
         return false;
    }
                
    WriteMDFile(mdMols, numMols, ofs, *pmol, indices);

    for(vector<OBMol*>::iterator  i = mdMols.begin(); i != mdMols.end(); ++i) 
    {
        delete *i;
    }

    //    

    return(true);
}

bool OOPSEFormat::AreSameFragments(OBMol& mol, vector<int>& frag1, vector<int>& frag2)
{
    if (frag1.size() != frag2.size())
    {
        return false;
    }

    //exact graph matching is a NP complete problem
    /** @todo using sparse matrix to store the connectivities*/
    for (unsigned int i =0 ; i < frag1.size(); ++i)
    {
        OBAtom* atom1 = mol.GetAtom(frag1[i]);
        OBAtom* atom2 = mol.GetAtom(frag2[i]);
        
        if (atom1->GetAtomicNum() != atom2->GetAtomicNum())
        {
            return false;
        }
    }
    return true;
}

struct SameAngle
{
  bool operator()(const triple<OBAtom*,OBAtom*,OBAtom*> t1, const triple<OBAtom*,OBAtom*,OBAtom*> t2) const
  {
    return (t1.second == t2.second) && ( (t1.first == t2.first && t1.third == t2.third) || (t1.first == t2.third && t1.third == t2.first));
  }
};

  void OOPSEFormat::findAngles(OBMol& mol) {
    /*
    //if already has data return
    if(mol.HasData(OBGenericDataType::AngleData))
        return;

    vector<OBEdgeBase*>::iterator bi1,bi2;
    OBBond* bond;
    OBAtom *a,*b,*c;

    set<triple<OBAtom*,OBAtom*,OBAtom*>, SameAngle> uniqueAngles;
    //loop through all bonds generating torsions
    for(bond = mol.BeginBond(bi1);bond;bond = mol.NextBond(bi1))
    {
        b = bond->GetBeginAtom();
        c = bond->GetEndAtom();
        if(b->IsHydrogen())
            continue;

        for(a = b->BeginNbrAtom(bi2);a;a = b->NextNbrAtom(bi2))
        {
            if(a == c)
                continue;          
            
            uniqueAngles.insert(triple<OBAtom*,OBAtom*,OBAtom*>(a, b, c));
        }
    }

    //get new data and attach it to molecule
    OBAngleData *angles = new OBAngleData;
    mol.SetData(angles);
    set<triple<OBAtom*,OBAtom*,OBAtom*>, SameAngle>::iterator i;

    for (i = uniqueAngles.begin(); i != uniqueAngles.end(); ++i) {
        OBAngle angle;
        angle.SetAtoms(i->first, i->second, i->second);
        angles->SetData(angle);
    }
    */
  }

  OBMol* OOPSEFormat::createMolFromFragment(OBMol& mol, vector<int>& fragment) {
    
    OBMol* newMol = new OBMol();
    newMol->ReserveAtoms(fragment.size());
    newMol->BeginModify();
    for(vector<int>::iterator i = fragment.begin(); i != fragment.end(); ++i) {
      OBAtom* newAtom = newMol->NewAtom();
      *newAtom = *mol.GetAtom(*i);
    }
    newMol->EndModify();
    newMol->ConnectTheDots();
    findAngles(*newMol);
    newMol->FindTorsions();
    return newMol;
  }
  
  void OOPSEFormat::WriteMDFile(vector<OBMol*> mols, vector<int> numMols, ostream& os, OBMol& mol, vector<int>& indices) {
    std::string indentLevel1("  ");
    std::string indentLevel2("    ");
    std::string molPrefix("MolName");
    unsigned int i;
    const int BUFFLEN = 1024;
    char buffer[BUFFLEN];
    string str, str1;
    
    
    os << "<OOPSE version=4>" << endl;
    os << "  <MetaData>" << endl << endl;
    
    for(i = 0; i < mols.size(); ++i) {
      OBMol* pmol = mols[i];

      pmol->ConnectTheDots();
      pmol->PerceiveBondOrders();
      //pmol->FindSSSR();
      //pmol->SetAromaticPerceived();
      //pmol->Kekulize();

      map<OBAtom*, int> atomMap;
      os << "molecule {\n";
      sprintf(buffer, "%d", i);
      os << indentLevel1 << "name = " << "\"" << molPrefix << buffer << "\"" << ";\n";
           
      //atom
      int ai = 0;
      FOR_ATOMS_OF_MOL(atom, *pmol ) {
        str = atom->GetType();
        ttab.SetFromType("INT");
        ttab.SetToType("INT");
        ttab.Translate(str1,str);
        os << indentLevel1 << "atom[" << ai << "] {\n";
        // os << indentLevel2 << "type = " << "\"" << etab.GetSymbol(atom->GetAtomicNum()) << "\"" << ";\n";
        os << indentLevel2 << "type = " << "\"" << str1 << "\"" << ";\n";
        os << indentLevel1 << "}\n";
        atomMap[&(*atom)] = ai++;
      }        
      os << "\n";
      
      //bond
      FOR_BONDS_OF_MOL(bond, *pmol ) {
        os << indentLevel1 << "bond {\n";
        os << indentLevel2 << "members(" << atomMap[bond->GetBeginAtom()] <<  ", " << atomMap[bond->GetEndAtom()] << ");\n";
        os << indentLevel1 << "}\n";
      }  
      /*
      //bend
      OBGenericData* pGenericData = pmol->GetData(OBGenericDataType::AngleData);
      OBAngleData* pAngleData = dynamic_cast<OBAngleData*>(pGenericData);
      vector<OBAngle> angles = pAngleData->GetData();
      
      os << indentLevel1 << "nBends = " << angles.size() << ";\n";        
      int bendIndex = 0;
      for (vector<OBAngle>::iterator ti = angles.begin(); ti != angles.end(); ++ti)
        {
          triple<OBAtom*, OBAtom*, OBAtom*> bendAtoms = ti->getAtoms();
          os << indentLevel1 << "bend[" << bendIndex++ << "] {\n";
          os << indentLevel2 << "member(" << atomMap[bendAtoms.first] <<  ", " << atomMap[bendAtoms.second] << atomMap[bendAtoms.third] <<");\n";
          os << indentLevel1 << "}\n";            
        }
      
      //torsion
      pGenericData = pmol->GetData(OBGenericDataType::TorsionData);
      OBTorsionData* pTorsionData = dynamic_cast<OBTorsionData*>(pGenericData);
      vector<OBTorsion> torsions = pTorsionData->GetData();
      vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > torsionArray;
      for (vector<OBTorsion>::iterator ti = torsions.begin(); ti != torsions.end(); ++ti)
        {
          vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> > tmpTorsions = ti->getTorsions();
          torsionArray.insert(torsionArray.end(), tmpTorsions.begin(), tmpTorsions.end());            
        }
      
      os << indentLevel1 << "nTorsions = " << torsionArray.size() << ";\n";
      int torsionIndex = 0;
      for (vector<quad<OBAtom*,OBAtom*,OBAtom*,OBAtom*> >::iterator ti = torsionArray.begin(); ti != torsionArray.end(); ++ti)
        {
          os << indentLevel1 << "torsion[" << torsionIndex++ << "] {\n";
          os << indentLevel2 << "member(" << atomMap[ti->first] <<  ", " << atomMap[ti->second] <<", " << atomMap[ti->third] <<", " << atomMap[ti->forth] << ");\n";
          os << indentLevel1 << "}\n";          
        }
      */
      os << "}" << endl;
      os << endl;
      
    }
    
    os << endl;

    
    for(i=0; i < mols.size(); ++i) {      
      os << "component{" << endl;
      sprintf(buffer, "%d", i);
      os << indentLevel1 << "type = " << molPrefix << buffer << ";" << endl;
      os << indentLevel1 << "nMol = " << numMols[i]<< ";" << endl;
      os << "}" << endl;
    }
    
    os << "  </MetaData>" << endl;
    os << "  <Snapshot>" << endl;
    os << "    <FrameData>" << endl;
    
    sprintf(buffer, "        Time: %.10g", 0.0);
    
    os << buffer << endl;

    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);

    os << buffer << endl;
    os << "    </FrameData>" << endl;
    os << "    <StuntDoubles>" << endl;

    OBAtom *atom;
    
    for(vector<int>::iterator i = indices.begin();i != indices.end(); ++i) {     
      atom = mol.GetAtom(*i);
      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);
      os << buffer << endl;
    }
    os << "    </StuntDoubles>" << endl;
    os << "  </Snapshot>" << endl;
    os << "</OOPSE>" << endl;
}
  
} //namespace OpenBabel

Revision:	3056
Committed:	Thu Oct 19 19:51:34 2006 UTC (17 years, 11 months ago) by gezelter
File size:	10226 byte(s)
Log Message:	Slight changes to oopseformat to allow importing of openbabel's perceived atom types into OOPSE md files
#	User	Rev	Content
1	tim	2440	/**********************************************************************
2			Copyright (C) 2000 by OpenEye Scientific Software, Inc.
3			Some portions Copyright (C) 2001-2005 by Geoffrey R. Hutchison
4			Some portions Copyright (C) 2004 by Chris Morley
5
6			This program is free software; you can redistribute it and/or modify
7			it under the terms of the GNU General Public License as published by
8			the Free Software Foundation version 2 of the License.
9
10			This program is distributed in the hope that it will be useful,
11			but WITHOUT ANY WARRANTY; without even the implied warranty of
12			MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13			GNU General Public License for more details.
14			***********************************************************************/
15
16	tim	2445	#include "oopseformat.hpp"
17	tim	2455	#include <fstream>
18	tim	2440	namespace OpenBabel
19			{
20
21			bool OOPSEFormat::WriteMolecule(OBBase* pOb, OBConversion* pConv)
22			{
23			OBMol* pmol = dynamic_cast<OBMol*>(pOb);
24			if(pmol==NULL)
25			return false;
26
27			vector<vector<int> > fragmentLists;
28			pmol->ContigFragList(fragmentLists);
29			OBBitVec unused;
30			vector<bool> used(fragmentLists.size(), 0);
31			vector<vector<int> > molecules;
32			vector<int> indices;
33			for(int i =0; i < used.size(); ++i) {
34			if (used[i])
35			{
36			continue;
37			}
38			used[i] = true;
39			vector<int> sameMolTypes;
40			sameMolTypes.push_back(i);
41			indices.insert(indices.end(), fragmentLists[i].begin(), fragmentLists[i].end());
42			for (int j = i + 1;j < used.size(); ++j)
43			{
44			if (used[j])
45			{
46			continue;
47			}
48
49			if (AreSameFragments(*pmol, fragmentLists[i], fragmentLists[j]))
50			{
51			sameMolTypes.push_back(j);
52			indices.insert(indices.end(), fragmentLists[j].begin(), fragmentLists[j].end());
53			used[j]=true;
54			}
55			}
56			molecules.push_back(sameMolTypes);
57
58			}
59
60			//
61			vector<OBMol*> mdMols;
62			vector<int> numMols;
63			for(vector<vector<int> >::iterator i = molecules.begin(); i != molecules.end(); ++i)
64			{
65			mdMols.push_back(createMolFromFragment(*pmol, fragmentLists[i->front()]));
66			numMols.push_back((*i).size());
67			}
68
69	tim	2455	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	gezelter	2970	WriteMDFile(mdMols, numMols, ofs, *pmol, indices);
83	tim	2440
84			for(vector<OBMol*>::iterator i = mdMols.begin(); i != mdMols.end(); ++i)
85			{
86			delete *i;
87			}
88
89			//
90
91			return(true);
92			}
93
94			bool OOPSEFormat::AreSameFragments(OBMol& mol, vector<int>& frag1, vector<int>& frag2)
95			{
96			if (frag1.size() != frag2.size())
97			{
98			return false;
99			}
100
101	tim	2469	//exact graph matching is a NP complete problem
102			/** @todo using sparse matrix to store the connectivities*/
103	tim	2440	for (unsigned int i =0 ; i < frag1.size(); ++i)
104			{
105	tim	2469	OBAtom* atom1 = mol.GetAtom(frag1[i]);
106			OBAtom* atom2 = mol.GetAtom(frag2[i]);
107
108			if (atom1->GetAtomicNum() != atom2->GetAtomicNum())
109	tim	2440	{
110			return false;
111			}
112			}
113			return true;
114			}
115
116			struct SameAngle
117			{
118			bool operator()(const triple<OBAtom,OBAtom,OBAtom> t1, const triple<OBAtom,OBAtom,OBAtom> t2) const
119			{
120			return (t1.second == t2.second) && ( (t1.first == t2.first && t1.third == t2.third) \|\| (t1.first == t2.third && t1.third == t2.first));
121			}
122			};
123
124	gezelter	2970	void OOPSEFormat::findAngles(OBMol& mol) {
125	tim	2440	/*
126			//if already has data return
127			if(mol.HasData(OBGenericDataType::AngleData))
128			return;
129
130			vector<OBEdgeBase*>::iterator bi1,bi2;
131			OBBond* bond;
132			OBAtom a,b,*c;
133
134			set<triple<OBAtom,OBAtom,OBAtom*>, SameAngle> uniqueAngles;
135			//loop through all bonds generating torsions
136			for(bond = mol.BeginBond(bi1);bond;bond = mol.NextBond(bi1))
137			{
138			b = bond->GetBeginAtom();
139			c = bond->GetEndAtom();
140			if(b->IsHydrogen())
141			continue;
142
143			for(a = b->BeginNbrAtom(bi2);a;a = b->NextNbrAtom(bi2))
144			{
145			if(a == c)
146			continue;
147
148			uniqueAngles.insert(triple<OBAtom,OBAtom,OBAtom*>(a, b, c));
149			}
150			}
151
152			//get new data and attach it to molecule
153			OBAngleData *angles = new OBAngleData;
154			mol.SetData(angles);
155			set<triple<OBAtom,OBAtom,OBAtom*>, SameAngle>::iterator i;
156
157			for (i = uniqueAngles.begin(); i != uniqueAngles.end(); ++i) {
158			OBAngle angle;
159			angle.SetAtoms(i->first, i->second, i->second);
160			angles->SetData(angle);
161			}
162			*/
163	gezelter	2970	}
164
165			OBMol* OOPSEFormat::createMolFromFragment(OBMol& mol, vector<int>& fragment) {
166
167	tim	2440	OBMol* newMol = new OBMol();
168			newMol->ReserveAtoms(fragment.size());
169			newMol->BeginModify();
170	gezelter	2970	for(vector<int>::iterator i = fragment.begin(); i != fragment.end(); ++i) {
171			OBAtom* newAtom = newMol->NewAtom();
172			newAtom = mol.GetAtom(*i);
173	tim	2440	}
174			newMol->EndModify();
175			newMol->ConnectTheDots();
176			findAngles(*newMol);
177			newMol->FindTorsions();
178			return newMol;
179	gezelter	2970	}
180
181			void OOPSEFormat::WriteMDFile(vector<OBMol*> mols, vector<int> numMols, ostream& os, OBMol& mol, vector<int>& indices) {
182	tim	2982	std::string indentLevel1(" ");
183			std::string indentLevel2(" ");
184	tim	2440	std::string molPrefix("MolName");
185	gezelter	2970	unsigned int i;
186	tim	2440	const int BUFFLEN = 1024;
187			char buffer[BUFFLEN];
188	gezelter	3056	string str, str1;
189	tim	2440
190	gezelter	2970
191			os << "<OOPSE version=4>" << endl;
192	tim	2982	os << " <MetaData>" << endl << endl;
193	gezelter	2970
194			for(i = 0; i < mols.size(); ++i) {
195			OBMol* pmol = mols[i];
196	gezelter	3056
197			pmol->ConnectTheDots();
198			pmol->PerceiveBondOrders();
199			//pmol->FindSSSR();
200			//pmol->SetAromaticPerceived();
201			//pmol->Kekulize();
202
203	gezelter	2970	map<OBAtom*, int> atomMap;
204			os << "molecule {\n";
205			sprintf(buffer, "%d", i);
206			os << indentLevel1 << "name = " << "\"" << molPrefix << buffer << "\"" << ";\n";
207	gezelter	3056
208	gezelter	2970	//atom
209			int ai = 0;
210			FOR_ATOMS_OF_MOL(atom, *pmol ) {
211	gezelter	3056	str = atom->GetType();
212			ttab.SetFromType("INT");
213			ttab.SetToType("INT");
214			ttab.Translate(str1,str);
215	gezelter	2970	os << indentLevel1 << "atom[" << ai << "] {\n";
216	gezelter	3056	// os << indentLevel2 << "type = " << "\"" << etab.GetSymbol(atom->GetAtomicNum()) << "\"" << ";\n";
217			os << indentLevel2 << "type = " << "\"" << str1 << "\"" << ";\n";
218	gezelter	2970	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	tim	2440	{
239	gezelter	2970	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	tim	2440	}
244	gezelter	2970
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	tim	2440	{
252	gezelter	2970	vector<quad<OBAtom,OBAtom,OBAtom,OBAtom> > tmpTorsions = ti->getTorsions();
253			torsionArray.insert(torsionArray.end(), tmpTorsions.begin(), tmpTorsions.end());
254	tim	2440	}
255	gezelter	2970
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	tim	2440	{
260	gezelter	2970	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	tim	2440	}
264	gezelter	2970	*/
265			os << "}" << endl;
266			os << endl;
267
268	tim	2440	}
269
270	gezelter	2970	os << endl;
271	tim	2982
272	gezelter	2970
273			for(i=0; i < mols.size(); ++i) {
274			os << "component{" << endl;
275			sprintf(buffer, "%d", i);
276	tim	2982	os << indentLevel1 << "type = " << molPrefix << buffer << ";" << endl;
277			os << indentLevel1 << "nMol = " << numMols[i]<< ";" << endl;
278			os << "}" << endl;
279	tim	2440	}
280
281	gezelter	2970	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	tim	2440
289	gezelter	2970	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	tim	2440	OBAtom *atom;
296
297	gezelter	2970	for(vector<int>::iterator i = indices.begin();i != indices.end(); ++i) {
298			atom = mol.GetAtom(*i);
299	gezelter	2984	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	gezelter	2970	os << buffer << endl;
301	tim	2440	}
302	gezelter	2970	os << " </StuntDoubles>" << endl;
303			os << " </Snapshot>" << endl;
304			os << "</OOPSE>" << endl;
305	tim	2440	}
306	gezelter	2970
307	tim	2440	} //namespace OpenBabel
308