diff --git a/modules/io/src/img/map_io_tiff_handler.cc b/modules/io/src/img/map_io_tiff_handler.cc index 233a43c2d4b0069df0b4543c930ac3f174984463..10c1a5c3f754a2a9b70cfd9da68f704b8e43cb45 100644 --- a/modules/io/src/img/map_io_tiff_handler.cc +++ b/modules/io/src/img/map_io_tiff_handler.cc @@ -40,6 +40,7 @@ #include <ost/io/io_exception.hh> #include "map_io_tiff_handler.hh" +#include "tiff_util.hh" namespace ost { namespace io { diff --git a/modules/io/src/img/map_io_tiff_handler.hh b/modules/io/src/img/map_io_tiff_handler.hh index 7fbfb901e8b44cc6d6616648a70abf842ec26505..4fced407477c08c4ab62f28f7169cbe34408894f 100644 --- a/modules/io/src/img/map_io_tiff_handler.hh +++ b/modules/io/src/img/map_io_tiff_handler.hh @@ -17,21 +17,21 @@ // along with this library; if not, write to the Free Software Foundation, Inc., // 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA //------------------------------------------------------------------------------ -#ifndef DX_IO_TIFF_IO_PLUGIN_OST_HH -#define DX_IO_TIFF_IO_PLUGIN_OST_HH +#ifndef DX_IO_tiff_IO_PLUGIN_OST_HH +#define DX_IO_tiff_IO_PLUGIN_OST_HH /* Authors: Ansgar Philippsen, Andreas Schenk */ #include <boost/logic/tribool.hpp> -#include "tiff_util.hh" #include "map_io_handler.hh" - +struct tiff; + namespace ost { namespace io { class DLLEXPORT_OST_IO TIF: public ImageFormatBase @@ -65,7 +65,7 @@ class DLLEXPORT_OST_IO TIF: public ImageFormatBase protected: TIF(String format_string, boost::logic::tribool normalize_on_save = boost::logic::indeterminate, Format bit_depth = OST_DEFAULT_FORMAT, bool sign = false, bool phasecolor = false, int subimage = -1); - void do_export(const img::MapHandle& image,TIFF* tfile,const TIF& formattif) const; + void do_export(const img::MapHandle& image,tiff* tfile,const TIF& formattif) const; private: @@ -84,10 +84,15 @@ class DLLEXPORT_OST_IO MapIOTiffHandler: public MapIOHandler /// \brief Map IO handler to read/write Tiff files /// /// This map IO handler reads and writes tiff-formatted files. - virtual void Import(img::MapHandle& sh, const boost::filesystem::path& loc,const ImageFormatBase& formatstruct); - virtual void Import(img::MapHandle& sh, std::istream& loc,const ImageFormatBase& formatstruct); - virtual void Export(const img::MapHandle& sh, const boost::filesystem::path& loc,const ImageFormatBase& formatstruct) const; - virtual void Export(const img::MapHandle& sh, std::ostream& loc,const ImageFormatBase& formatstruct) const; + virtual void Import(img::MapHandle& sh, const boost::filesystem::path& loc, + const ImageFormatBase& formatstruct); + virtual void Import(img::MapHandle& sh, std::istream& loc, + const ImageFormatBase& formatstruct); + virtual void Export(const img::MapHandle& sh, + const boost::filesystem::path& loc, + const ImageFormatBase& formatstruct) const; + virtual void Export(const img::MapHandle& sh, std::ostream& loc, + const ImageFormatBase& formatstruct) const; static bool MatchContent(unsigned char* header); static bool MatchType(const ImageFormatBase& type); static bool MatchSuffix(const String& loc); @@ -96,10 +101,13 @@ class DLLEXPORT_OST_IO MapIOTiffHandler: public MapIOHandler protected: - TIFF* open_subimage_file(const boost::filesystem::path& location,const TIF& formattif); - TIFF* open_subimage_stream(std::istream& location,const TIF& formattif); - void load_image_data(TIFF* tfile, img::ImageHandle& image, const TIF& formattif); - virtual void do_export(const img::MapHandle& sh,TIFF* tfile,TIF& formatstruct) const; + tiff* open_subimage_file(const boost::filesystem::path& location, + const TIF& formattif); + tiff* open_subimage_stream(std::istream& location,const TIF& formattif); + void load_image_data(tiff* tfile, img::ImageHandle& image, + const TIF& formattif); + virtual void do_export(const img::MapHandle& sh,tiff* tfile, + TIF& formatstruct) const; }; diff --git a/modules/mol/alg/src/svd_superpose.cc b/modules/mol/alg/src/svd_superpose.cc index 6608372c13a2372196041a0d6dfee8f2e312e647..2f770167f236fa5c1a7dc02aef7016c77823fc90 100644 --- a/modules/mol/alg/src/svd_superpose.cc +++ b/modules/mol/alg/src/svd_superpose.cc @@ -21,13 +21,15 @@ * Author Juergen Haas */ #include <stdexcept> - +#include <iostream> #include <boost/bind.hpp> +#include <Eigen/Core> +#include <Eigen/Array> #include <Eigen/SVD> #include <Eigen/LU> -#include <iostream> + #include <ost/base.hh> #include <ost/geom/vec3.hh> #include <ost/mol/alg/svd_superpose.hh> @@ -36,7 +38,13 @@ namespace ost { namespace mol { namespace alg { + using boost::bind; +typedef Eigen::Matrix<Real,3,1> EVec3; +typedef Eigen::Matrix<Real,3,3> EMat3; +typedef Eigen::Matrix<Real,1,3> ERVec3; +typedef Eigen::Matrix<Real,Eigen::Dynamic,Eigen::Dynamic> EMatX; +typedef Eigen::Matrix<Real,1,Eigen::Dynamic> ERVecX; Real CalculateRMSD(const mol::EntityView& ev1, const mol::EntityView& ev2, @@ -93,22 +101,55 @@ Real calc_rmsd_for_atom_lists(const mol::AtomViewList& atoms1, return rmsd; } +class SuperposerSVDImpl { +public: + SuperposerSVDImpl(int natoms, bool alloc_atoms): + natoms_(natoms), alloc_atoms_(alloc_atoms) + { -SuperposerSVD::SuperposerSVD(int natoms, bool alloc_atoms): - natoms_(natoms), alloc_atoms_(alloc_atoms) -{ + ERVec3 avg1_=EVec3::Zero(); + ERVec3 avg2_=EVec3::Zero(); - ERVec3 avg1_=EVec3::Zero(); - ERVec3 avg2_=EVec3::Zero(); + if (alloc_atoms_) { + atoms1_=EMatX::Zero(natoms,3); + atoms2_=EMatX::Zero(natoms,3); + } - if (alloc_atoms_) { - atoms1_=EMatX::Zero(natoms,3); - atoms2_=EMatX::Zero(natoms,3); } +private: + int natoms_; + bool alloc_atoms_; + EMatX atoms1_; + EMatX atoms2_; +public: + void EntToMatrices(const mol::EntityView& ev1, const mol::EntityView& ev2); + geom::Vec3 EigenVec3ToVec3(const EVec3 &vec); + geom::Mat3 EigenMat3ToMat3(const EMat3 &mat); + EVec3 Vec3ToEigenRVec(const geom::Vec3 &vec); + EVec3 Vec3ToEigenVec(const geom::Vec3 &vec); + EMatX SubtractVecFromMatrixRows(EMatX Mat, + ERVecX Vec); + SuperpositionResult Run(const mol::AtomViewList& atoms1, + const mol::AtomViewList& atoms2); + SuperpositionResult DoSuperposition(); + EIGEN_MAKE_ALIGNED_OPERATOR_NEW //needed only if *fixed* sized matrices are + //members of classes +}; + +SuperposerSVD::~SuperposerSVD() +{ + assert(impl_); + delete impl_; +} + +SuperposerSVD::SuperposerSVD(int natoms, bool alloc_atoms): + impl_(new SuperposerSVDImpl(natoms, alloc_atoms)) +{ + } -geom::Vec3 SuperposerSVD::EigenVec3ToVec3(const EVec3 &vec){ +geom::Vec3 SuperposerSVDImpl::EigenVec3ToVec3(const EVec3 &vec){ #if defined _MSC_VER geom::Vec3 myvec3; for (int i=0; i<3; ++i) { @@ -120,7 +161,8 @@ geom::Vec3 SuperposerSVD::EigenVec3ToVec3(const EVec3 &vec){ #endif } -geom::Mat3 SuperposerSVD::EigenMat3ToMat3(const EMat3 &mat){ +geom::Mat3 SuperposerSVDImpl::EigenMat3ToMat3(const EMat3 &mat) +{ #if defined _MSC_VER geom::Mat3 mymat3; for (int i=0; i<3; ++i) { @@ -135,7 +177,8 @@ geom::Mat3 SuperposerSVD::EigenMat3ToMat3(const EMat3 &mat){ } -EVec3 SuperposerSVD::Vec3ToEigenRVec(const geom::Vec3 &vec){ +EVec3 SuperposerSVDImpl::Vec3ToEigenRVec(const geom::Vec3 &vec) +{ #if defined _MSC_VER EVec3 myvec3=EVec3::Zero(); for (int counter=0; counter<3; ++counter) { @@ -147,7 +190,8 @@ EVec3 SuperposerSVD::Vec3ToEigenRVec(const geom::Vec3 &vec){ #endif } -EVec3 SuperposerSVD::Vec3ToEigenVec(const geom::Vec3 &vec){ +EVec3 SuperposerSVDImpl::Vec3ToEigenVec(const geom::Vec3 &vec) +{ #if defined _MSC_VER EVec3 myvec3=EVec3::Zero(); for (int counter=0; counter<3; ++counter) { @@ -159,8 +203,9 @@ EVec3 SuperposerSVD::Vec3ToEigenVec(const geom::Vec3 &vec){ #endif } -EMatX SuperposerSVD::SubtractVecFromMatrixRows(EMatX Mat, - ERVecX Vec){ +EMatX SuperposerSVDImpl::SubtractVecFromMatrixRows(EMatX Mat, + ERVecX Vec) +{ for (int row=0; row<Mat.rows();++row) { Mat.row(row)-=Vec; @@ -168,8 +213,8 @@ EMatX SuperposerSVD::SubtractVecFromMatrixRows(EMatX Mat, return Mat; } -void SuperposerSVD::EntToMatrices(const mol::EntityView& ev1, - const mol::EntityView& ev2) +void SuperposerSVDImpl::EntToMatrices(const mol::EntityView& ev1, + const mol::EntityView& ev2) { //Now iterate over atoms in entities and extract coords into Nx3 matrices //for superposition @@ -183,23 +228,23 @@ void SuperposerSVD::EntToMatrices(const mol::EntityView& ev1, for (int counter=0; a_ev1_end!=a_ev1; ++a_ev1, ++a_ev2, ++counter) { av1=*a_ev1; pos=av1.GetPos(); - e=SuperposerSVD::Vec3ToEigenVec(pos); + e=SuperposerSVDImpl::Vec3ToEigenVec(pos); atoms1_.row(counter)=e ; av2=*a_ev2; pos=av2.GetPos(); - e=SuperposerSVD::Vec3ToEigenVec(pos); + e=SuperposerSVDImpl::Vec3ToEigenVec(pos); atoms2_.row(counter)=e; } } -SuperpositionResult SuperposerSVD::DoSuperposition() +SuperpositionResult SuperposerSVDImpl::DoSuperposition() { ERVec3 avg1_=atoms1_.colwise().sum()/atoms1_.rows(); ERVec3 avg2_=atoms2_.colwise().sum()/atoms2_.rows(); //center the structures - atoms1_=this->SuperposerSVD::SubtractVecFromMatrixRows(atoms1_, avg1_); - atoms2_=this->SuperposerSVD::SubtractVecFromMatrixRows(atoms2_, avg2_); + atoms1_=this->SubtractVecFromMatrixRows(atoms1_, avg1_); + atoms2_=this->SubtractVecFromMatrixRows(atoms2_, avg2_); EMatX atoms2=atoms2_.transpose(); //single value decomposition @@ -221,9 +266,9 @@ SuperpositionResult SuperposerSVD::DoSuperposition() SuperpositionResult res; // prepare rmsd calculation - geom::Vec3 shift=SuperposerSVD::EigenVec3ToVec3(avg2_); - geom::Vec3 com_vec=-SuperposerSVD::EigenVec3ToVec3(avg1_); - geom::Mat3 rot=SuperposerSVD::EigenMat3ToMat3(ERot.transpose()); + geom::Vec3 shift=SuperposerSVDImpl::EigenVec3ToVec3(avg2_); + geom::Vec3 com_vec=-SuperposerSVDImpl::EigenVec3ToVec3(avg1_); + geom::Mat3 rot=SuperposerSVDImpl::EigenMat3ToMat3(ERot.transpose()); geom::Mat4 mat4_com, mat4_rot, mat4_shift; mat4_rot.PasteRotation(rot); mat4_shift.PasteTranslation(shift); @@ -238,14 +283,20 @@ SuperpositionResult SuperposerSVD::DoSuperposition() SuperpositionResult SuperposerSVD::Run(const mol::EntityView& ev1, const mol::EntityView& ev2) { - this->EntToMatrices(ev1, ev2); - SuperpositionResult r=this->DoSuperposition(); + impl_->EntToMatrices(ev1, ev2); + SuperpositionResult r=impl_->DoSuperposition(); r.rmsd=CalculateRMSD(ev1, ev2, r.transformation); return r; } SuperpositionResult SuperposerSVD::Run(const mol::AtomViewList& atoms1, const mol::AtomViewList& atoms2) +{ + return impl_->Run(atoms1, atoms2); +} + +SuperpositionResult SuperposerSVDImpl::Run(const mol::AtomViewList& atoms1, + const mol::AtomViewList& atoms2) { mol::AtomViewList::const_iterator a_ev1=atoms1.begin(); mol::AtomViewList::const_iterator a_ev1_end=atoms1.end(); @@ -254,23 +305,17 @@ SuperpositionResult SuperposerSVD::Run(const mol::AtomViewList& atoms1, for (int counter=0; a_ev1_end!=a_ev1; ++a_ev1, ++a_ev2, ++counter) { mol::AtomView av=*a_ev1; geom::Vec3 pos=av.GetPos(); - ERVec3 e=SuperposerSVD::Vec3ToEigenVec(pos); + ERVec3 e=SuperposerSVDImpl::Vec3ToEigenVec(pos); atoms1_.row(counter)=e ; av=*a_ev2; pos=av.GetPos(); - e=SuperposerSVD::Vec3ToEigenVec(pos); + e=SuperposerSVDImpl::Vec3ToEigenVec(pos); atoms2_.row(counter)=e; } return this->DoSuperposition(); } -SuperpositionResult SuperposerSVD::Run(EMatX atoms1, - EMatX atoms2) -{ - atoms1_=atoms1; - atoms2_=atoms2; - return this->DoSuperposition(); -} + SuperpositionResult SuperposeAtoms(const mol::AtomViewList& atoms1, const mol::AtomViewList& atoms2, diff --git a/modules/mol/alg/src/svd_superpose.hh b/modules/mol/alg/src/svd_superpose.hh index de7403bc1c83658da4c3a2eb196256cff63630b3..80f0762785df0d71bbd90a588c14ed79f30f37e1 100644 --- a/modules/mol/alg/src/svd_superpose.hh +++ b/modules/mol/alg/src/svd_superpose.hh @@ -23,9 +23,6 @@ #ifndef OST_SUPERPOSITION_HH #define OST_SUPERPOSITION_HH -#include <Eigen/Core> -#include <Eigen/Array> - #include <ost/base.hh> #include <ost/geom/geom.hh> #include <ost/mol/mol.hh> @@ -48,46 +45,22 @@ struct DLLEXPORT_OST_MOL_ALG SuperpositionResult { mol::EntityView entity_view2; }; - -typedef Eigen::Matrix<Real,3,1> EVec3; -typedef Eigen::Matrix<Real,3,3> EMat3; -typedef Eigen::Matrix<Real,1,3> ERVec3; -typedef Eigen::Matrix<Real,Eigen::Dynamic,Eigen::Dynamic> EMatX; -typedef Eigen::Matrix<Real,1,Eigen::Dynamic> ERVecX; +class SuperposerSVDImpl; /// \brief effiently superpose a bunch of models with the same number of atoms -/// Choose either two EntityViews, two AtomViewLists or two EigenMatrices of the -/// same size to superpose. +/// Choose either two EntityViews or two AtomViewLists. class DLLEXPORT_OST_MOL_ALG SuperposerSVD { public: SuperposerSVD(int natoms, bool alloc_atoms); - + ~SuperposerSVD(); SuperpositionResult Run(const mol::EntityView& ev1, const mol::EntityView& ev2); SuperpositionResult Run(const mol::AtomViewList& atoms1, const mol::AtomViewList& atoms2); - SuperpositionResult Run(const EMatX atoms1, - const EMatX atoms2); private: - - - int natoms_; - bool alloc_atoms_; - EMatX atoms1_; - EMatX atoms2_; - void EntToMatrices(const mol::EntityView& ev1, const mol::EntityView& ev2); - geom::Vec3 EigenVec3ToVec3(const EVec3 &vec); - geom::Mat3 EigenMat3ToMat3(const EMat3 &mat); - EVec3 Vec3ToEigenRVec(const geom::Vec3 &vec); - EVec3 Vec3ToEigenVec(const geom::Vec3 &vec); - EMatX SubtractVecFromMatrixRows(EMatX Mat, - ERVecX Vec); - EIGEN_MAKE_ALIGNED_OPERATOR_NEW //needed only if *fixed* sized matrices are - //members of classes - SuperpositionResult DoSuperposition(); - + SuperposerSVDImpl* impl_; }; /// \brief takes the corresponding atoms and superposes them diff --git a/modules/qa/pymod/CMakeLists.txt b/modules/qa/pymod/CMakeLists.txt index 5923ef71e1869586fbd6c9ca121a385783b75c52..c4fff8b16c9a6f8fc579717f9d2ccd988c59002c 100644 --- a/modules/qa/pymod/CMakeLists.txt +++ b/modules/qa/pymod/CMakeLists.txt @@ -1,5 +1,4 @@ set(OST_QA_PYMOD_SOURCES - export_ca_distance.cc export_interaction.cc export_torsion.cc export_packing.cc diff --git a/modules/qa/pymod/export_ca_distance.cc b/modules/qa/pymod/export_ca_distance.cc deleted file mode 100644 index ef6bd2846dd5ad617f7c3af9ee8496ed4571bd7c..0000000000000000000000000000000000000000 --- a/modules/qa/pymod/export_ca_distance.cc +++ /dev/null @@ -1,71 +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 <boost/python.hpp> - - -#include <ost/qa/ca_distance.hh> - -using namespace ost::qa; -using namespace ost; - -///\todo put this into a header? seems like other might need this as well. -namespace boost{ - - template<class T> - inline T* get_pointer( boost::shared_ptr<const T> const& p ){ - return const_cast< T* >( p.get() ); - } - -} - -namespace boost{ namespace python{ - - template<class T> - struct pointee< boost::shared_ptr<T const> >{ - typedef T type; - }; - -} } //boost::python - -namespace utils{ - - template< class T > - void register_shared_ptrs_to_python(){ - namespace bpl = boost::python; - bpl::register_ptr_to_python< boost::shared_ptr< T > >(); - bpl::register_ptr_to_python< boost::shared_ptr< const T > >(); - bpl::implicitly_convertible< boost::shared_ptr< T >, - boost::shared_ptr< const T > >(); - } - -} -using namespace boost::python; - -void export_CaDistance() -{ - def("CalculateCalphaDistance", &CalculateCalphaDistance); - - class_<CalphaDistance, boost::noncopyable>("CalphaDistance", no_init) - .def("IsDistanceDefined", &CalphaDistance::IsDistanceDefined, - args("number_a", "number_b")) - .def("GetDistance", &CalphaDistance::GetDistance, - args("number_a", "number_b")) - ; - utils::register_shared_ptrs_to_python<CalphaDistance>(); -} diff --git a/modules/qa/pymod/wrap_qa.cc b/modules/qa/pymod/wrap_qa.cc index 543ccc9f518d46b58f5352273643ac737eb69628..bb6519155ef12476d31c1e1fac7752069ca182d2 100644 --- a/modules/qa/pymod/wrap_qa.cc +++ b/modules/qa/pymod/wrap_qa.cc @@ -21,13 +21,11 @@ void export_Torsion(); void export_Interaction(); void export_Packing(); -void export_CaDistance(); void export_Clash(); BOOST_PYTHON_MODULE(_qa) { export_Torsion(); export_Interaction(); - export_Packing(); - export_CaDistance(); + export_Packing(); export_Clash(); } diff --git a/modules/qa/src/CMakeLists.txt b/modules/qa/src/CMakeLists.txt index 088906a2f1b57d6e972a6e4e90f5928e9fb3e3db..cc75ea3b85e8c3a6de3f98aeb75732ccf14f0849 100644 --- a/modules/qa/src/CMakeLists.txt +++ b/modules/qa/src/CMakeLists.txt @@ -10,7 +10,6 @@ solis_torsion_potential.hh interaction_statistics.hh all_atom_potential.hh packing_statistics.hh -ca_distance.hh atom_types.hh packing_potential.hh clash_score.hh @@ -24,7 +23,6 @@ solis_torsion_statistics.cc solis_torsion_potential.cc interaction_statistics.cc packing_statistics.cc -ca_distance.cc amino_acids.cc atom_types.cc all_atom_potential.cc diff --git a/modules/qa/src/ca_distance.cc b/modules/qa/src/ca_distance.cc deleted file mode 100644 index 1b9da37ed6de0046a88866eddbe7ce5f1a074d12..0000000000000000000000000000000000000000 --- a/modules/qa/src/ca_distance.cc +++ /dev/null @@ -1,103 +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 "ca_distance.hh" -#include <ost/geom/geom.hh> -#include <Eigen/Array> -/* - Author: Marco Biasini - */ - -namespace ost { namespace qa { - -namespace { -class DistanceCalculator : public mol::EntityVisitor { -public: - DistanceCalculator(const mol::AtomHandle& atom_handle, - Eigen::MatrixXd mat) - : atom_a_(atom_handle), distance_mat_(mat) { - - } - virtual bool VisitAtom(const mol::AtomHandle& atom) { - Real d=geom::Length(atom_a_.GetPos()-atom.GetPos()); - int i=atom_a_.GetResidue().GetNumber().GetNum(); - int j=atom.GetResidue().GetNumber().GetNum(); - distance_mat_(i-1, j-1)= d; - return true; - } -private: - mol::AtomHandle atom_a_; - Eigen::MatrixXd distance_mat_; -}; - -} - -mol::ResNum CalphaDistance::GetBiggestResNumber(const mol::EntityView& e) const { - assert(e.GetChainCount()==1); - return e.GetChainList().front().GetResidueList().back().GetNumber(); -} - -CalphaDistancePtr CalculateCalphaDistance(const mol::EntityView& a, - const mol::EntityView& b) { - if (a.GetChainCount()!=1 || b.GetChainCount()!=1) { - throw Error("One or more views have chain count not equal to 1"); - } - return CalphaDistancePtr(new CalphaDistance(a, b)); -} - -CalphaDistance::CalphaDistance(const mol::EntityView& a, const mol::EntityView& b) - : a_(a), b_(b), distance_mat_() { - assert(a_.GetChainCount()==1); - assert(b_.GetChainCount()==1); - mol::ResNum last_a=this->GetBiggestResNumber(a_); - mol::ResNum last_b=this->GetBiggestResNumber(b_); - distance_mat_=Eigen::MatrixXd::Constant(last_a.GetNum(), last_b.GetNum(), -1); - this->CalculateDistances(); -} - -void CalphaDistance::CalculateDistances() { - a_.Apply(*this); -} - -bool CalphaDistance::IsDistanceDefined(const mol::ResNum& number_a, - const mol::ResNum& number_b) const { - int i=number_a.GetNum(); - int j=number_b.GetNum(); - return (1<=i && 1<=j && static_cast<int>(distance_mat_.cols())>=i && - static_cast<int>(distance_mat_.rows())>=j && - distance_mat_(i-1, j-1)>=0.0); - -} - - -bool CalphaDistance::VisitAtom(const mol::AtomHandle& atom) { - DistanceCalculator dist_calc(atom, distance_mat_); - b_.Apply(dist_calc); - return true; -} - -Real CalphaDistance::GetDistance(const mol::ResNum& number_a, - const mol::ResNum& number_b) const { - if (this->IsDistanceDefined(number_a, number_b)) - return distance_mat_(number_a.GetNum()-1, - number_b.GetNum()-1); - throw Error("distance between C-alpha atoms is undefined"); -} - -}} - diff --git a/modules/qa/src/ca_distance.hh b/modules/qa/src/ca_distance.hh deleted file mode 100644 index 4690ad934848cdeef098c872be3cc678f285b312..0000000000000000000000000000000000000000 --- a/modules/qa/src/ca_distance.hh +++ /dev/null @@ -1,92 +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 -//------------------------------------------------------------------------------ -#ifndef OST_QA_CA_DISTANCE_HH -#define OST_QA_CA_DISTANCE_HH - -/* - Author: Marco Biasini - */ - -#include <ost/mol/mol.hh> -#include <ost/message.hh> -#include <Eigen/Core> -#include <boost/scoped_ptr.hpp> -#include <ost/qa/module_config.hh> - -namespace ost { namespace qa { - -class CalphaDistance; - -typedef boost::shared_ptr<CalphaDistance> CalphaDistancePtr; - -/// \brief calculate C-alpha distance matrix -/// -/// \param a first entity view. The entity view is required to contain only -/// C-alpha atoms. -/// -/// \param b second entity view. The entity view is required to contain only -/// C-alpha atoms. -/// -/// \throws Error, if one or both of the two views do not -/// contain exactly one chain -CalphaDistancePtr DLLEXPORT_OST_QA CalculateCalphaDistance(const mol::EntityView& a, - const mol::EntityView& b); - - -class DLLEXPORT_OST_QA CalphaDistance : public mol::EntityVisitor { -public: - CalphaDistance(const mol::EntityView& a, const mol::EntityView& b); - - virtual bool VisitAtom(const mol::AtomHandle& atom); - - /// \brief Check if distance between two C-alpha is defined - /// - /// Check whether the distance between C-alpha of residue with \c number_a - /// in the first entity and residue with \c number_b in the second entity is - /// defined, i.e. if both C-alphas are part of the view. - /// - /// \param number_a residue number in first entity view. - /// - /// \param number_b residue number in second entity view. - bool IsDistanceDefined(const mol::ResNum& number_a, - const mol::ResNum& number_b) const; - - /// \brief Get distance between two C-alpha atoms. - /// - /// \throws Error when distance between the two C-alpha atoms is - /// not defined. - /// - /// \sa IsDistanceDefined - Real GetDistance(const mol::ResNum& number_a, - const mol::ResNum& number_b) const; - - -private: - CalphaDistance(const CalphaDistance&); - mol::ResNum GetBiggestResNumber(const mol::EntityView& e) const; - void CalculateDistances(); - - mol::EntityView a_; - mol::EntityView b_; - Eigen::MatrixXd distance_mat_; -}; - -}} - -#endif