diff --git a/modules/io/pymod/wrap_io.cc b/modules/io/pymod/wrap_io.cc
index 4e58c969f10f536f7047f56c933ad796c14a3d20..79920189ef60ba1cb71c122bb566907291646bef 100644
--- a/modules/io/pymod/wrap_io.cc
+++ b/modules/io/pymod/wrap_io.cc
@@ -30,7 +30,7 @@ using namespace boost::python;
#include <ost/io/mol/entity_io_crd_handler.hh>
#include <ost/io/mol/entity_io_sdf_handler.hh>
#include <ost/io/mol/pdb_reader.hh>
-#include <ost/io/mol/save_dcd.hh>
+#include <ost/io/mol/dcd_io.hh>
using namespace ost;
using namespace ost::io;
@@ -65,10 +65,14 @@ BOOST_PYTHON_FUNCTION_OVERLOADS(save_entity_handle_ov,
BOOST_PYTHON_FUNCTION_OVERLOADS(save_entity_view_ov,
save_ent_view, 2, 3)
-BOOST_PYTHON_FUNCTION_OVERLOADS(load_charmm_trj_ov,
- LoadCHARMMTraj, 2, 3)
BOOST_PYTHON_FUNCTION_OVERLOADS(save_charmm_trj_ov,
SaveCHARMMTraj, 3, 4)
+
+mol::CoordGroupHandle load_dcd1(const String& c, const String& t) {return LoadCHARMMTraj(c,t);}
+mol::CoordGroupHandle load_dcd2(const String& c, const String& t, unsigned int s) {return LoadCHARMMTraj(c,t,s);}
+mol::CoordGroupHandle load_dcd3(const mol::EntityHandle& e, const String& t) {return LoadCHARMMTraj(e,t);}
+mol::CoordGroupHandle load_dcd4(const mol::EntityHandle& e, const String& t, unsigned int s) {return LoadCHARMMTraj(e,t,s);}
+
}
void export_pdb_io();
@@ -111,7 +115,10 @@ BOOST_PYTHON_MODULE(_io)
def("LoadSDF", &LoadSDF);
def("LoadCRD", &LoadCRD);
- def("LoadCHARMMTraj",LoadCHARMMTraj,load_charmm_trj_ov());
+ def("LoadCHARMMTraj",load_dcd1);
+ def("LoadCHARMMTraj",load_dcd2);
+ def("LoadCHARMMTraj",load_dcd3);
+ def("LoadCHARMMTraj",load_dcd4);
def("SaveCHARMMTraj",SaveCHARMMTraj,save_charmm_trj_ov());
export_pdb_io();
#if OST_IMG_ENABLED
diff --git a/modules/io/src/mol/CMakeLists.txt b/modules/io/src/mol/CMakeLists.txt
index 76ebcfec1b051adf2e855c06e5a5a0edf52b5a63..b3fd346b6e87552242a509f853a0eb5ff093fd38 100644
--- a/modules/io/src/mol/CMakeLists.txt
+++ b/modules/io/src/mol/CMakeLists.txt
@@ -8,14 +8,14 @@ save_entity.cc
load_entity.cc
surface_io_msms_handler.cc
load_surface.cc
-save_dcd.cc
+dcd_io.cc
star_parser.cc
PARENT_SCOPE
)
set(OST_IO_MOL_HEADERS
star_parser.hh
-save_dcd.hh
+dcd_io.hh
entity_io_crd_handler.hh
pdb_io.hh
entity_io_handler.hh
diff --git a/modules/io/src/mol/dcd_io.cc b/modules/io/src/mol/dcd_io.cc
new file mode 100644
index 0000000000000000000000000000000000000000..efe9f1dfb368cc888c5785279af00ee8510a3fb7
--- /dev/null
+++ b/modules/io/src/mol/dcd_io.cc
@@ -0,0 +1,315 @@
+//------------------------------------------------------------------------------
+// This file is part of the OpenStructure project <www.openstructure.org>
+//
+// Copyright (C) 2008-2010 by the OpenStructure authors
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License as published by the Free
+// Software Foundation; either version 3.0 of the License, or (at your option)
+// any later version.
+// This library 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 Lesser General Public License for more
+// details.
+//
+// You should have received a copy of the GNU Lesser General Public License
+// along with this library; if not, write to the Free Software Foundation, Inc.,
+// 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
+//------------------------------------------------------------------------------
+
+/*
+ Authors: Ansgar Philippsen, Marco Biasini
+*/
+
+#include <fstream>
+#include <boost/iostreams/filter/gzip.hpp>
+#include <boost/filesystem/convenience.hpp>
+
+#include <ost/stdint.hh>
+#include <ost/log.hh>
+#include <ost/profile.hh>
+#include <ost/conop/conop.hh>
+#include <ost/mol/xcs_editor.hh>
+
+#include "dcd_io.hh"
+#include "pdb_reader.hh"
+#include "pdb_writer.hh"
+#include "entity_io_crd_handler.hh"
+#include <ost/io/io_exception.hh>
+#include <ost/io/swap_util.hh>
+
+namespace ost { namespace io {
+
+/*
+ icntrl[1]: number of coordinate sets in file
+ icntrl[2]: number of previous dynamic steps
+ icntrl[3]: frequency for saving coordinates
+ icntlr[4]: number of steps for creation run
+*/
+
+namespace {
+
+struct DCDHeader {
+ char hdrr[4];
+ int icntrl[20];
+ int ntitle;
+ char title[1024];
+ int num, istep, freq,nstep;
+ int t_atom_count,f_atom_count, atom_count;
+};
+
+mol::CoordGroupHandle load_dcd(const mol::AtomHandleList& alist,
+ const String& trj_fn,
+ unsigned int stride)
+{
+ Profile profile_load("LoadCHARMMTraj");
+
+ bool gap_flag = true;
+
+ boost::filesystem::path trj_f(trj_fn);
+ boost::filesystem::ifstream ff(trj_f);
+
+ DCDHeader header;
+ char dummy[4];
+ bool swap_flag=false;
+
+ LOGN_MESSAGE("importing trajectory data");
+
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(header.hdrr,sizeof(char)*4);
+ ff.read(reinterpret_cast<char*>(header.icntrl),sizeof(int)*20);
+
+ if(header.icntrl[1]<0 || header.icntrl[1]>1e6) {
+ // nonsense atom count, try swapping
+ swap_int(header.icntrl,20);
+ if(header.icntrl[1]<0 || header.icntrl[1]>1e6) {
+ throw(IOException("LoadCHARMMTraj: nonsense atom count in header"));
+ } else {
+ LOGN_MESSAGE("LoadCHARMMTraj: byte-swapping");
+ swap_flag=true;
+ }
+ }
+
+ bool skip_flag=false;
+
+ if(header.icntrl[19]!=0) { // CHARMM format
+ skip_flag=(header.icntrl[10]!=0);
+ if(skip_flag) {
+ LOGN_VERBOSE("LoadCHARMMTraj: using CHARMM format with per-frame header");
+ } else {
+ LOGN_VERBOSE("LoadCHARMMTraj: using CHARMM format");
+ }
+ } else {
+ // XPLOR format
+ LOGN_VERBOSE("LoadCHARMMTraj: using XPLOR format");
+ }
+
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(reinterpret_cast<char*>(&header.ntitle),sizeof(int));
+ if(swap_flag) swap_int(&header.ntitle,1);
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(header.title,sizeof(char)*header.ntitle);
+ header.title[header.ntitle]='\0';
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(reinterpret_cast<char*>(&header.t_atom_count),sizeof(int));
+ if(swap_flag) swap_int(&header.t_atom_count,1);
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ header.num=header.icntrl[0];
+ header.istep=header.icntrl[1];
+ header.freq=header.icntrl[2];
+ header.nstep=header.icntrl[3];
+ header.f_atom_count=header.icntrl[8];
+ header.atom_count=header.t_atom_count-header.f_atom_count;
+
+ LOGN_MESSAGE("LoadCHARMMTraj: " << header.num << " trajectories with "
+ << header.atom_count << " atoms (" << header.f_atom_count
+ << " fixed) each");
+
+ if(alist.size() != static_cast<size_t>(header.t_atom_count)) {
+ LOGN_ERROR("LoadCHARMMTraj: atom count missmatch: " << alist.size()
+ << " in coordinate file, " << header.t_atom_count
+ << " in each traj frame");
+ throw(IOException("invalid trajectory"));
+ }
+
+ mol::CoordGroupHandle cg=CreateCoordGroup(alist);
+ std::vector<geom::Vec3> clist(header.t_atom_count);
+ std::vector<float> xlist(header.t_atom_count);
+
+ for(int i=0;i<header.num;++i) {
+ if(skip_flag) ff.seekg(14*4,std::ios_base::cur);
+ // read each frame
+ if(!ff) {
+ /* premature EOF */
+ LOGN_ERROR("LoadCHARMMTraj: premature end of file, " << i
+ << " frames read");
+ break;
+ }
+ // x coord
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(reinterpret_cast<char*>(&xlist[0]),sizeof(float)*xlist.size());
+ if(swap_flag) swap_float(&xlist[0],xlist.size());
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ for(uint j=0;j<clist.size();++j) {
+ clist[j].SetX(xlist[j]);
+ }
+ // y coord
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(reinterpret_cast<char*>(&xlist[0]),sizeof(float)*xlist.size());
+ if(swap_flag) swap_float(&xlist[0],xlist.size());
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ for(uint j=0;j<clist.size();++j) {
+ clist[j].SetY(xlist[j]);
+ }
+ // z coord
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ ff.read(reinterpret_cast<char*>(&xlist[0]),sizeof(float)*xlist.size());
+ if(swap_flag) swap_float(&xlist[0],xlist.size());
+ if(gap_flag) ff.read(dummy,sizeof(dummy));
+ for(uint j=0;j<clist.size();++j) {
+ clist[j].SetZ(xlist[j]);
+ }
+ if(i%stride) {
+ cg.AddFrame(clist);
+ }
+ }
+
+ ff.get();
+ if(!ff.eof()) {
+ LOGN_VERBOSE("LoadCHARMMTraj: unexpected trailing file data, bytes read: "
+ << ff.tellg());
+ }
+ return cg;
+}
+
+} // anon ns
+
+mol::CoordGroupHandle LoadCHARMMTraj(const String& crd_fn,
+ const String& trj_fn,
+ unsigned int stride)
+{
+ boost::filesystem::path crd_f(crd_fn);
+
+ // first the coordinate reference file
+ conop::BuilderP builder = conop::Conopology::Instance().GetBuilder();
+ mol::EntityHandle ent=mol::CreateEntity();
+ mol::XCSEditor editor=ent.RequestXCSEditor(mol::BUFFERED_EDIT);
+ std::vector<mol::AtomHandle> alist;
+ if(boost::filesystem::extension(crd_f)==".pdb") {
+ PDBReader reader(crd_f);
+ reader.SetFlags(PDB::SEQUENTIAL_ATOM_IMPORT);
+ LOGN_MESSAGE("importing coordinate data");
+ reader.Import(ent);
+ alist = reader.GetSequentialAtoms();
+ } else if (boost::filesystem::extension(crd_f)==".crd") {
+ CRDReader reader(crd_f);
+ LOGN_MESSAGE("importing coordinate data");
+ reader.Import(ent);
+ alist = reader.GetSequentialAtoms();
+ } else {
+ throw(IOException("unsupported coordinate file format"));
+ }
+ conop::Conopology::Instance().ConnectAll(builder,ent);
+
+ return load_dcd(alist,trj_fn,stride);
+}
+
+
+mol::CoordGroupHandle DLLEXPORT_OST_IO LoadCHARMMTraj(const mol::EntityHandle& e,
+ const String& trj_fn,
+ unsigned int stride)
+{
+ return load_dcd(e.GetAtomList(),trj_fn,stride);
+}
+
+
+namespace {
+
+void write_dcd_hdr(std::ofstream& out,
+ const mol::CoordGroupHandle& coord_group,
+ unsigned int stepsize)
+{
+ int32_t magic_number=84;
+ char crd[]={'C', 'O', 'R', 'D'};
+ out.write(reinterpret_cast<char*>(&magic_number), 4);
+ out.write(crd, 4);
+ int32_t num_frames=coord_group.GetFrameCount()/stepsize;
+ out.write(reinterpret_cast<char*>(&num_frames), 4);
+ int32_t zero=0;
+ // write zero for istart
+ out.write(reinterpret_cast<char*>(&zero), 4);
+ int32_t one=1;
+ // number of timesteps between dcd saves.
+ out.write(reinterpret_cast<char*>(&one), 4);
+ // write spacer of 5 blank integers
+ for (int i=0; i<5; ++i) {
+ out.write(reinterpret_cast<char*>(&zero), 4);
+ }
+ // write number of fixed atoms
+
+ out.write(reinterpret_cast<char*>(&zero), 4);
+ float delta=1.0;
+ out.write(reinterpret_cast<char*>(&delta), 4);
+ // write spacer of 10 blank integers
+ for (int i=0; i <10; ++i) {
+ out.write(reinterpret_cast<char*>(&zero), 4);
+ }
+ out.write(reinterpret_cast<char*>(&magic_number), 4);
+ // we don't write any titles for now. This means that the block has only to
+ // accomodate one int, harr, harr, harr.
+ int32_t title_block_size=4;
+ out.write(reinterpret_cast<char*>(&title_block_size), 4);
+ out.write(reinterpret_cast<char*>(&zero), 4);
+ int32_t four=4;
+ // again block size for titles?
+ out.write(reinterpret_cast<char*>(&four), 4);
+ // has to be 4
+ out.write(reinterpret_cast<char*>(&four), 4);
+ int32_t atom_count=coord_group.GetAtomCount();
+ out.write(reinterpret_cast<char*>(&atom_count), 4);
+ out.write(reinterpret_cast<char*>(&four), 4);
+}
+
+}
+
+void SaveCHARMMTraj(const mol::CoordGroupHandle& coord_group,
+ const String& pdb_filename, const String& dcd_filename,
+ unsigned int stepsize)
+{
+ if(stepsize==0) stepsize=1;
+ if(!pdb_filename.empty()) {
+ PDBWriter writer(pdb_filename);
+ writer.Write(coord_group.GetAtomList());
+ }
+ std::ofstream out(dcd_filename.c_str(), std::ios::binary);
+ write_dcd_hdr(out, coord_group,stepsize);
+ int frame_count=coord_group.GetFrameCount();
+ int atom_count=coord_group.GetAtomCount();
+
+ std::vector<float> x(atom_count);
+ std::vector<float> y(atom_count);
+ std::vector<float> z(atom_count);
+ int32_t out_n=atom_count*4;
+ for (int i=0; i<frame_count; i+=stepsize) {
+ mol::CoordFramePtr frame=coord_group.GetFrame(i);
+ int k=0;
+ for (mol::CoordFrame::iterator j=frame->begin(),
+ e=frame->end(); j!=e; ++j, ++k) {
+ //geom::Vec3 v=*j;
+ x[k]=float((*j)[0]);
+ y[k]=float((*j)[1]);
+ z[k]=float((*j)[2]);
+ }
+ out.write(reinterpret_cast<char*>(&out_n), 4);
+ out.write(reinterpret_cast<char*>(&x[0]), out_n);
+ out.write(reinterpret_cast<char*>(&out_n), 4);
+ out.write(reinterpret_cast<char*>(&out_n), 4);
+ out.write(reinterpret_cast<char*>(&y[0]), out_n);
+ out.write(reinterpret_cast<char*>(&out_n), 4);
+ out.write(reinterpret_cast<char*>(&out_n), 4);
+ out.write(reinterpret_cast<char*>(&z[0]), out_n);
+ out.write(reinterpret_cast<char*>(&out_n), 4);
+ }
+}
+
+}} // ns
diff --git a/modules/io/src/mol/save_dcd.hh b/modules/io/src/mol/dcd_io.hh
similarity index 53%
rename from modules/io/src/mol/save_dcd.hh
rename to modules/io/src/mol/dcd_io.hh
index 7b3c4ce390182e33cb1f409b71fe875d078155ce..b20a3fb0f4b7bcbb2ccf6b47b1b206fdeb3b3e57 100644
--- a/modules/io/src/mol/save_dcd.hh
+++ b/modules/io/src/mol/dcd_io.hh
@@ -16,23 +16,47 @@
// along with this library; if not, write to the Free Software Foundation, Inc.,
// 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
//------------------------------------------------------------------------------
-#ifndef OST_SAVE_DCD_HH
-#define OST_SAVE_DCD_HH
+#ifndef OST_IO_ENTITY_DCD_IO_HH
+#define OST_IO_ENTITY_DCD_IO_HH
/*
- Author: Marco Biasini
+ CHARMM trajectory IO
+ Authors: Ansgar Philippsen, Marco Biasini
*/
-
+
#include <ost/io/module_config.hh>
#include <ost/mol/coord_group.hh>
namespace ost { namespace io {
-
-/// \brief export coord group as PDB file and DCD trajectory
+
+/*! \brief import a CHARMM trajectory in dcd format
+ requires the coordinate and the trajectory file; the format
+ of the coordinate file will be automatically deduced from the extension
+ the optional stride parameter will cause only every nth frame to be loaded
+*/
+mol::CoordGroupHandle DLLEXPORT_OST_IO LoadCHARMMTraj(const String& coord,
+ const String& trj,
+ unsigned int stride=1);
+
+/*! \brief import a CHARMM trajectory in dcd format with an existing entity
+ requires the existing entity and the trajectory file - obviously the
+ atom layout of the entity must match the trajectory file
+*/
+mol::CoordGroupHandle DLLEXPORT_OST_IO LoadCHARMMTraj(const mol::EntityHandle& e,
+ const String& trj,
+ unsigned int stride=1);
+
+
+/*! \brief export coord group as PDB file and DCD trajectory
+ if the pdb filename is an empty string, it won't be exported
+ the optional stride parameter will cause every nth frame to be exported
+ */
void DLLEXPORT_OST_IO SaveCHARMMTraj(const mol::CoordGroupHandle& coord_group,
const String& pdb_filename,
const String& dcd_filename,
- unsigned int stepsize=1);
-}}
+ unsigned int stride=1);
+
+
+}} // ns
#endif
diff --git a/modules/io/src/mol/entity_io_crd_handler.cc b/modules/io/src/mol/entity_io_crd_handler.cc
index a7870b33ed239d98bee19519cdfbabb7a73fa7ca..ce4a1679b765cc8aa13fbddb873b82c33f427871 100644
--- a/modules/io/src/mol/entity_io_crd_handler.cc
+++ b/modules/io/src/mol/entity_io_crd_handler.cc
@@ -35,10 +35,10 @@
#include <ost/profile.hh>
#include <ost/io/io_exception.hh>
-#include "entity_io_crd_handler.hh"
-#include "pdb_reader.hh"
#include <ost/io/swap_util.hh>
+#include "entity_io_crd_handler.hh"
+
namespace ost { namespace io {
using boost::format;
@@ -219,173 +219,6 @@ mol::EntityHandle LoadCRD(const String& file_name)
return ent;
}
-namespace {
-
-struct TrjHeader {
- char hdrr[4];
- int icntrl[20];
- int ntitle;
- char title[1024];
- int num, istep, freq,nstep;
- int t_atom_count,f_atom_count, atom_count;
-};
-
-}
-
-/*
- icntrl[1]: number of coordinate sets in file
- icntrl[2]: number of previous dynamic steps
- icntrl[3]: frequency for saving coordinates
- icntlr[4]: number of steps for creation run
-*/
-
-mol::CoordGroupHandle LoadCHARMMTraj(const String& crd_fn,
- const String& trj_fn,
- int flags)
-{
- Profile profile_load("LoadCHARMMTraj");
-
- bool gap_flag = true;
-
- boost::filesystem::path crd_f(crd_fn);
- boost::filesystem::path trj_f(trj_fn);
-
-
-
- // first the coordinate reference file
- conop::BuilderP builder = conop::Conopology::Instance().GetBuilder();
- mol::EntityHandle ent=mol::CreateEntity();
- mol::XCSEditor editor=ent.RequestXCSEditor(mol::BUFFERED_EDIT);
- std::vector<mol::AtomHandle> alist;
- if(boost::filesystem::extension(crd_f)==".pdb") {
- PDBReader reader(crd_f);
- reader.SetFlags(PDB::SEQUENTIAL_ATOM_IMPORT);
- LOGN_MESSAGE("importing coordinate data");
- reader.Import(ent);
- alist = reader.GetSequentialAtoms();
- } else if (boost::filesystem::extension(crd_f)==".crd") {
- CRDReader reader(crd_f);
- LOGN_MESSAGE("importing coordinate data");
- reader.Import(ent);
- alist = reader.GetSequentialAtoms();
- } else {
- throw(IOException("unsupported coordinate file format"));
- }
- conop::Conopology::Instance().ConnectAll(builder,ent);
-
- // now the trajectory
- boost::filesystem::ifstream ff(trj_f);
-
- TrjHeader header;
- char dummy[4];
- bool swap_flag=false;
-
- LOGN_MESSAGE("importing trajectory data");
-
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(header.hdrr,sizeof(char)*4);
- ff.read(reinterpret_cast<char*>(header.icntrl),sizeof(int)*20);
-
- if(header.icntrl[1]<0 || header.icntrl[1]>1e6) {
- // nonsense atom count, try swapping
- swap_int(header.icntrl,20);
- if(header.icntrl[1]<0 || header.icntrl[1]>1e6) {
- throw(IOException("LoadCHARMMTraj: nonsense atom count in header"));
- } else {
- LOGN_MESSAGE("LoadCHARMMTraj: byte-swapping");
- swap_flag=true;
- }
- }
-
- bool skip_flag=false;
-
- if(header.icntrl[19]!=0) { // CHARMM format
- skip_flag=(header.icntrl[10]!=0);
- if(skip_flag) {
- LOGN_VERBOSE("LoadCHARMMTraj: using CHARMM format with per-frame header");
- } else {
- LOGN_VERBOSE("LoadCHARMMTraj: using CHARMM format");
- }
- } else {
- // XPLOR format
- LOGN_VERBOSE("LoadCHARMMTraj: using XPLOR format");
- }
-
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(reinterpret_cast<char*>(&header.ntitle),sizeof(int));
- if(swap_flag) swap_int(&header.ntitle,1);
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(header.title,sizeof(char)*header.ntitle);
- header.title[header.ntitle]='\0';
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(reinterpret_cast<char*>(&header.t_atom_count),sizeof(int));
- if(swap_flag) swap_int(&header.t_atom_count,1);
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- header.num=header.icntrl[0];
- header.istep=header.icntrl[1];
- header.freq=header.icntrl[2];
- header.nstep=header.icntrl[3];
- header.f_atom_count=header.icntrl[8];
- header.atom_count=header.t_atom_count-header.f_atom_count;
-
- LOGN_MESSAGE("LoadCHARMMTraj: " << header.num << " trajectories with "
- << header.atom_count << " atoms (" << header.f_atom_count
- << " fixed) each");
-
- if(alist.size() != static_cast<size_t>(header.t_atom_count)) {
- LOGN_ERROR("LoadCHARMMTraj: atom count missmatch: " << alist.size()
- << " in coordinate file, " << header.t_atom_count
- << " in each traj frame");
- throw(IOException("invalid trajectory"));
- }
-
- mol::CoordGroupHandle cg=CreateCoordGroup(alist);
- std::vector<geom::Vec3> clist(header.t_atom_count);
- std::vector<float> xlist(header.t_atom_count);
-
- for(int i=0;i<header.num;++i) {
- if(skip_flag) ff.seekg(14*4,std::ios_base::cur);
- // read each frame
- if(!ff) {
- /* premature EOF */
- LOGN_ERROR("LoadCHARMMTraj: premature end of file, " << i
- << " frames read");
- break;
- }
- // x coord
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(reinterpret_cast<char*>(&xlist[0]),sizeof(float)*xlist.size());
- if(swap_flag) swap_float(&xlist[0],xlist.size());
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- for(uint j=0;j<clist.size();++j) {
- clist[j].SetX(xlist[j]);
- }
- // y coord
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(reinterpret_cast<char*>(&xlist[0]),sizeof(float)*xlist.size());
- if(swap_flag) swap_float(&xlist[0],xlist.size());
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- for(uint j=0;j<clist.size();++j) {
- clist[j].SetY(xlist[j]);
- }
- // z coord
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- ff.read(reinterpret_cast<char*>(&xlist[0]),sizeof(float)*xlist.size());
- if(swap_flag) swap_float(&xlist[0],xlist.size());
- if(gap_flag) ff.read(dummy,sizeof(dummy));
- for(uint j=0;j<clist.size();++j) {
- clist[j].SetZ(xlist[j]);
- }
- cg.AddFrame(clist);
- }
-
- ff.get();
- if(!ff.eof()) {
- LOGN_VERBOSE("LoadCHARMMTraj: unexpected trailing file data, bytes read: "
- << ff.tellg());
- }
- return cg;
-}
void EntityIOCRDHandler::Import(mol::EntityHandle& ent,
std::istream& stream)
diff --git a/modules/io/src/mol/entity_io_crd_handler.hh b/modules/io/src/mol/entity_io_crd_handler.hh
index 72984150e6684309972309f6fc19b36b520ff8f0..ceb635f81a2723900b70a0e6215e468c8fbe356d 100644
--- a/modules/io/src/mol/entity_io_crd_handler.hh
+++ b/modules/io/src/mol/entity_io_crd_handler.hh
@@ -25,8 +25,6 @@
*/
#include <ost/io/mol/entity_io_handler.hh>
-#include <ost/mol/coord_group.hh>
-
#include <boost/iostreams/filtering_stream.hpp>
#include <boost/filesystem/fstream.hpp>
@@ -80,14 +78,6 @@ typedef EntityIOHandlerFactory<EntityIOCRDHandler> EntityIOCRDHandlerFactory;
mol::EntityHandle DLLEXPORT_OST_IO LoadCRD(const String& file_name);
-/*
- the flag can tweak the behavior of the importer; right now, the zero bit
- can be used to switch old format compatibility on
-*/
-mol::CoordGroupHandle DLLEXPORT_OST_IO LoadCHARMMTraj(const String& crd,
- const String& trj,
- int flags=0);
-
}} // ns
#endif
diff --git a/modules/io/src/mol/save_dcd.cc b/modules/io/src/mol/save_dcd.cc
deleted file mode 100644
index 4ed551260ea5efa95c712484996d987dc0b10216..0000000000000000000000000000000000000000
--- a/modules/io/src/mol/save_dcd.cc
+++ /dev/null
@@ -1,118 +0,0 @@
-//------------------------------------------------------------------------------
-// This file is part of the OpenStructure project <www.openstructure.org>
-//
-// Copyright (C) 2008-2010 by the OpenStructure authors
-//
-// This library is free software; you can redistribute it and/or modify it under
-// the terms of the GNU Lesser General Public License as published by the Free
-// Software Foundation; either version 3.0 of the License, or (at your option)
-// any later version.
-// This library 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 Lesser General Public License for more
-// details.
-//
-// You should have received a copy of the GNU Lesser General Public License
-// along with this library; if not, write to the Free Software Foundation, Inc.,
-// 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
-//------------------------------------------------------------------------------
-
-#include <fstream>
-
-#include <ost/io/mol/pdb_writer.hh>
-#include <ost/io/mol/save_dcd.hh>
-#include <ost/stdint.hh>
-/*
- Author: Marco Biasini
- */
-namespace ost { namespace io {
-
-
-namespace {
-
-void write_dcd_hdr(std::ofstream& out,
- const mol::CoordGroupHandle& coord_group,
- unsigned int stepsize)
-{
- int32_t magic_number=84;
- char crd[]={'C', 'O', 'R', 'D'};
- out.write(reinterpret_cast<char*>(&magic_number), 4);
- out.write(crd, 4);
- int32_t num_frames=coord_group.GetFrameCount()/stepsize;
- out.write(reinterpret_cast<char*>(&num_frames), 4);
- int32_t zero=0;
- // write zero for istart
- out.write(reinterpret_cast<char*>(&zero), 4);
- int32_t one=1;
- // number of timesteps between dcd saves.
- out.write(reinterpret_cast<char*>(&one), 4);
- // write spacer of 5 blank integers
- for (int i=0; i<5; ++i) {
- out.write(reinterpret_cast<char*>(&zero), 4);
- }
- // write number of fixed atoms
-
- out.write(reinterpret_cast<char*>(&zero), 4);
- float delta=1.0;
- out.write(reinterpret_cast<char*>(&delta), 4);
- // write spacer of 10 blank integers
- for (int i=0; i <10; ++i) {
- out.write(reinterpret_cast<char*>(&zero), 4);
- }
- out.write(reinterpret_cast<char*>(&magic_number), 4);
- // we don't write any titles for now. This means that the block has only to
- // accomodate one int, harr, harr, harr.
- int32_t title_block_size=4;
- out.write(reinterpret_cast<char*>(&title_block_size), 4);
- out.write(reinterpret_cast<char*>(&zero), 4);
- int32_t four=4;
- // again block size for titles?
- out.write(reinterpret_cast<char*>(&four), 4);
- // has to be 4
- out.write(reinterpret_cast<char*>(&four), 4);
- int32_t atom_count=coord_group.GetAtomCount();
- out.write(reinterpret_cast<char*>(&atom_count), 4);
- out.write(reinterpret_cast<char*>(&four), 4);
-}
-
-}
-
-void SaveCHARMMTraj(const mol::CoordGroupHandle& coord_group,
- const String& pdb_filename, const String& dcd_filename,
- unsigned int stepsize)
-{
- if(stepsize==0) stepsize=1;
- PDBWriter writer(pdb_filename);
- writer.Write(coord_group.GetAtomList());
- std::ofstream out(dcd_filename.c_str(), std::ios::binary);
- write_dcd_hdr(out, coord_group,stepsize);
- int frame_count=coord_group.GetFrameCount();
- int atom_count=coord_group.GetAtomCount();
-
- std::vector<float> x(atom_count);
- std::vector<float> y(atom_count);
- std::vector<float> z(atom_count);
- int32_t out_n=atom_count*4;
- for (int i=0; i<frame_count; i+=stepsize) {
- mol::CoordFramePtr frame=coord_group.GetFrame(i);
- int k=0;
- for (mol::CoordFrame::iterator j=frame->begin(),
- e=frame->end(); j!=e; ++j, ++k) {
- //geom::Vec3 v=*j;
- x[k]=float((*j)[0]);
- y[k]=float((*j)[1]);
- z[k]=float((*j)[2]);
- }
- out.write(reinterpret_cast<char*>(&out_n), 4);
- out.write(reinterpret_cast<char*>(&x[0]), out_n);
- out.write(reinterpret_cast<char*>(&out_n), 4);
- out.write(reinterpret_cast<char*>(&out_n), 4);
- out.write(reinterpret_cast<char*>(&y[0]), out_n);
- out.write(reinterpret_cast<char*>(&out_n), 4);
- out.write(reinterpret_cast<char*>(&out_n), 4);
- out.write(reinterpret_cast<char*>(&z[0]), out_n);
- out.write(reinterpret_cast<char*>(&out_n), 4);
- }
-}
-
-}}