From 7101ae694c9cd77e1b228f11959e496005134cbd Mon Sep 17 00:00:00 2001
From: Marco Biasini <marco.biasini@unibas.ch>
Date: Wed, 13 Oct 2010 18:34:24 +0200
Subject: [PATCH] SaveCHARMMTraj saves the PDB files in CHARMM format
---
modules/io/pymod/__init__.py | 4 +-
modules/io/pymod/export_pdb_io.cc | 1 +
modules/io/src/mol/dcd_io.cc | 3 +-
modules/io/src/mol/pdb_writer.cc | 71 ++++++++++++++++++++-----------
modules/io/src/mol/pdb_writer.hh | 13 +++---
5 files changed, 58 insertions(+), 34 deletions(-)
diff --git a/modules/io/pymod/__init__.py b/modules/io/pymod/__init__.py
index 5534c540d..b43415cb6 100644
--- a/modules/io/pymod/__init__.py
+++ b/modules/io/pymod/__init__.py
@@ -145,7 +145,7 @@ def LoadPDB(filename, restrict_chains="", no_hetatms=False,
PDB.PopFlags()
raise
-def SavePDB(models, filename):
+def SavePDB(models, filename, dialect='PDB'):
"""
Save entity or list of entities to disk. If a list of entities is supplied the
PDB file will be saved as a multi PDB file. Each of the entities is wrapped
@@ -157,7 +157,7 @@ def SavePDB(models, filename):
"""
if not getattr(models, '__len__', None):
models=[models]
- writer=PDBWriter(filename)
+ writer=PDBWriter(filename, dialect=='CHARMM')
try:
if len(models)>1:
PDB.PushFlags(PDB.Flags() |PDB.WRITE_MULTIPLE_MODELS)
diff --git a/modules/io/pymod/export_pdb_io.cc b/modules/io/pymod/export_pdb_io.cc
index 1275ea7a6..27e62dc1c 100644
--- a/modules/io/pymod/export_pdb_io.cc
+++ b/modules/io/pymod/export_pdb_io.cc
@@ -64,6 +64,7 @@ void export_pdb_io()
;
class_<PDBWriter, boost::noncopyable>("PDBWriter", init<String>())
+ .def(init<String, bool>())
.def("Write", write_a)
.def("Write", write_b)
;
diff --git a/modules/io/src/mol/dcd_io.cc b/modules/io/src/mol/dcd_io.cc
index 86a59ed42..b7723de0a 100644
--- a/modules/io/src/mol/dcd_io.cc
+++ b/modules/io/src/mol/dcd_io.cc
@@ -242,6 +242,7 @@ public:
mol::CoordSource(atoms), filename_(filename),
stream_(filename.c_str(), std::ios::binary), loaded_(false), stride_(stride)
{
+ frame_count_=0;
this->SetMutable(false);
frame_=mol::CoordFramePtr(new mol::CoordFrame(atoms.size()));
}
@@ -373,7 +374,7 @@ void SaveCHARMMTraj(const mol::CoordGroupHandle& coord_group,
{
if(stepsize==0) stepsize=1;
if(!pdb_filename.empty()) {
- PDBWriter writer(pdb_filename);
+ PDBWriter writer(pdb_filename, true);
writer.Write(coord_group.GetAtomList());
}
std::ofstream out(dcd_filename.c_str(), std::ios::binary);
diff --git a/modules/io/src/mol/pdb_writer.cc b/modules/io/src/mol/pdb_writer.cc
index 580c40752..8b16c2cce 100644
--- a/modules/io/src/mol/pdb_writer.cc
+++ b/modules/io/src/mol/pdb_writer.cc
@@ -61,7 +61,7 @@ bool shift_left(const String& atom_name, bool is_hetatm,
void write_atom(std::ostream& ostr, FormattedLine& line,
const mol::AtomHandle& atom, int atomnum,
- bool is_pqr)
+ bool is_pqr, bool charmm_style)
{
mol::ResidueHandle res=atom.GetResidue();
char ins_code=res.GetNumber().GetInsCode();
@@ -83,20 +83,32 @@ void write_atom(std::ostream& ostr, FormattedLine& line,
} else {
line(13, 3)=fmt::RPadded(atom_name);
}
- if (res.GetKey().size()>3) {
- throw IOException("Residue name '"+res.GetQualifiedName()+
- "' is too long for PDB output. At most 3 characters are "
- "allowed");
+ if (charmm_style) {
+ if (res.GetKey().size()>4) {
+ throw IOException("Residue name '"+res.GetName()+
+ "' is too long for CHARMM-PDB output. At most 4 "
+ "characters are allowed");
+ }
+ line(17, 4)=fmt::LPadded(res.GetKey());
+ } else {
+ if (res.GetKey().size()>3) {
+ throw IOException("Residue name '"+res.GetName()+
+ "' is too long for PDB output. At most 3 characters "
+ "are allowed");
+ }
+ line(17, 3)=fmt::LPadded(res.GetKey());
}
- line(17, 3)=fmt::LPadded(res.GetKey());
+
String chain_name=res.GetChain().GetName();
- if (chain_name.size()>0) {
- if (chain_name.size()==1) {
- line[21]=chain_name[0];
- } else {
- throw IOException("PDB format only support single-letter chain names");
- }
+ if (!charmm_style) {
+ if (chain_name.size()>0) {
+ if (chain_name.size()==1) {
+ line[21]=chain_name[0];
+ } else {
+ throw IOException("PDB format only support single-letter chain names");
+ }
+ }
}
line(22, 4)=fmt::LPaddedInt(res.GetNumber().GetNum());
if (ins_code!=0) {
@@ -116,7 +128,9 @@ void write_atom(std::ostream& ostr, FormattedLine& line,
line(54, 6)=fmt::LPaddedFloat(atom.GetOccupancy(), 2);
line(60, 6)=fmt::LPaddedFloat(atom.GetBFactor(), 2);
}
-
+ if (charmm_style) {
+ line(72, 4)=fmt::RPadded(chain_name);
+ }
line(76, 2)=fmt::LPadded(atom.GetElement());
ostr << line;
} else {
@@ -142,7 +156,9 @@ void write_atom(std::ostream& ostr, FormattedLine& line,
line(54, 6)=fmt::LPaddedFloat(atom.GetOccupancy(), 2);
line(60, 6)=fmt::LPaddedFloat(atom.GetBFactor(), 2);
}
-
+ if (charmm_style) {
+ line(72, 4)=fmt::RPadded(chain_name);
+ }
line(76, 2)=fmt::LPadded(atom.GetElement());
ostr << line;
}
@@ -171,15 +187,16 @@ void write_conect(std::ostream& ostr, int atom_index,
class PDBWriterImpl : public mol::EntityVisitor {
public:
PDBWriterImpl(std::ostream& ostream, FormattedLine& line,
- std::map<long,int>& atom_indices)
+ std::map<long,int>& atom_indices, bool charmm_style)
: ostr_(ostream), counter_(0), is_pqr_(false),
- atom_indices_(atom_indices), line_(line), peptide_(false) {
+ atom_indices_(atom_indices), line_(line), peptide_(false),
+ charmm_style_(charmm_style) {
}
private:
public:
virtual bool VisitAtom(const mol::AtomHandle& atom) {
counter_++;
- write_atom(ostr_, line_, atom, counter_, is_pqr_);
+ write_atom(ostr_, line_, atom, counter_, is_pqr_, charmm_style_);
if (atom.GetAtomProps().is_hetatm) {
atom_indices_[atom.GetHashCode()]=counter_;
}
@@ -240,6 +257,7 @@ private:
FormattedLine& line_;
mol::ResidueHandle prev_;
bool peptide_;
+ bool charmm_style_;
};
class PDBConectWriterImpl : public mol::EntityVisitor {
@@ -287,19 +305,21 @@ struct ForcePOSIX {
}
-PDBWriter::PDBWriter(std::ostream& stream):
- outfile_(), outstream_(stream), mol_count_(0), line_(80)
+PDBWriter::PDBWriter(std::ostream& stream, bool charmm_style):
+ outfile_(), outstream_(stream), mol_count_(0), line_(80),
+ charmm_style_(charmm_style)
{
}
-PDBWriter::PDBWriter(const boost::filesystem::path& filename):
+PDBWriter::PDBWriter(const boost::filesystem::path& filename, bool charmm_style):
outfile_(filename.file_string().c_str()), outstream_(outfile_),
- mol_count_(0), line_(80)
+ mol_count_(0), line_(80), charmm_style_(charmm_style)
{}
-PDBWriter::PDBWriter(const String& filename):
- outfile_(filename.c_str()), outstream_(outfile_), mol_count_(0), line_(80)
+PDBWriter::PDBWriter(const String& filename, bool charmm_style):
+ outfile_(filename.c_str()), outstream_(outfile_), mol_count_(0), line_(80),
+ charmm_style_(charmm_style)
{}
void PDBWriter::WriteModelLeader()
@@ -325,7 +345,7 @@ void PDBWriter::WriteModel(H ent)
{
ForcePOSIX posix = ForcePOSIX();
this->WriteModelLeader();
- PDBWriterImpl writer(outstream_,line_, atom_indices_);
+ PDBWriterImpl writer(outstream_,line_, atom_indices_, charmm_style_);
if (PDB::Flags() & PDB::PQR_FORMAT) {
writer.SetIsPQR(true);
}
@@ -353,7 +373,8 @@ void PDBWriter::Write(const mol::AtomHandleList& atoms)
mol::ChainHandle last_chain;
for (mol::AtomHandleList::const_iterator i=atoms.begin(),
e=atoms.end(); i!=e; ++i, ++counter) {
- write_atom(outstream_, line_, *i, counter, (PDB::Flags() & PDB::PQR_FORMAT) != 0);
+ write_atom(outstream_, line_, *i, counter,
+ (PDB::Flags() & PDB::PQR_FORMAT) != 0, charmm_style_);
}
this->WriteModelTrailer();
}
diff --git a/modules/io/src/mol/pdb_writer.hh b/modules/io/src/mol/pdb_writer.hh
index 40ed1cbff..c89657215 100644
--- a/modules/io/src/mol/pdb_writer.hh
+++ b/modules/io/src/mol/pdb_writer.hh
@@ -40,9 +40,9 @@ namespace ost { namespace io {
class DLLEXPORT_OST_IO PDBWriter {
public:
- PDBWriter(const String& filename);
- PDBWriter(const boost::filesystem::path& filename);
- PDBWriter(std::ostream& outstream);
+ PDBWriter(const String& filename, bool charmm_style=false);
+ PDBWriter(const boost::filesystem::path& filename, bool charmm_style=false);
+ PDBWriter(std::ostream& outstream, bool charmm_style=false);
void Write(const mol::EntityView& ent);
void Write(const mol::EntityHandle& ent);
@@ -56,11 +56,12 @@ private:
void WriteModelLeader();
void WriteModelTrailer();
- std::ofstream outfile_;
- std::ostream& outstream_;
- int mol_count_;
+ std::ofstream outfile_;
+ std::ostream& outstream_;
+ int mol_count_;
std::map<long, int> atom_indices_;
FormattedLine line_;
+ bool charmm_style_;
};
}}
--
GitLab