From 071f7b6f4428b94d68ae060d89b2b097c8625a50 Mon Sep 17 00:00:00 2001 From: Marco Biasini <marco.biasini@unibas.ch> Date: Mon, 3 Jan 2011 16:53:59 +0100 Subject: [PATCH] implement ReducedPotential and ReducedStatistics for views --- modules/qa/pymod/export_reduced.cc | 18 +++++-- modules/qa/src/impl/reduced_impl.cc | 74 +++++++++++++++++++--------- modules/qa/src/impl/reduced_impl.hh | 10 +++- modules/qa/src/reduced_potential.cc | 13 ++++- modules/qa/src/reduced_potential.hh | 2 +- modules/qa/src/reduced_statistics.cc | 19 +++++++ modules/qa/src/reduced_statistics.hh | 5 +- 7 files changed, 109 insertions(+), 32 deletions(-) diff --git a/modules/qa/pymod/export_reduced.cc b/modules/qa/pymod/export_reduced.cc index 0c230fb7d..3384c4529 100644 --- a/modules/qa/pymod/export_reduced.cc +++ b/modules/qa/pymod/export_reduced.cc @@ -23,10 +23,17 @@ #include <ost/qa/reduced_potential.hh> using namespace ost::qa; - +using namespace ost::mol; using namespace boost::python; -uint64_t (ReducedStatistics::*get_count)(AminoAcid,AminoAcid,int,int)=&ReducedStatistics::GetCount; +uint64_t (ReducedStatistics::*get_count)(AminoAcid, + AminoAcid, + int,int)=&ReducedStatistics::GetCount; +void (ReducedStatistics::*ex_h)(EntityHandle)=&ReducedStatistics::Extract; +void (ReducedStatistics::*ex_v)(EntityView)=&ReducedStatistics::Extract; + +Real (ReducedPotential::*get_h)(EntityHandle,bool)=&ReducedPotential::GetTotalEnergy; +Real (ReducedPotential::*get_v)(EntityView,bool)=&ReducedPotential::GetTotalEnergy; void export_Reduced() { class_<ReducedStatOptions>("ReducedStatOptions", init<>()) @@ -43,7 +50,8 @@ void export_Reduced() class_<ReducedStatistics, ReducedStatisticsPtr>("ReducedStatistics", no_init) .def(init<Real,Real,uint,uint,uint>((arg("l_cutoff"), arg("u_cutoff"), arg("num_ang_bins"), arg("num_dist_bins"), arg("sequence_sep")))) - .def("Extract", &ReducedStatistics::Extract) + .def("Extract", ex_h, arg("ent")) + .def("Extract", ex_v, arg("ent")) .add_property("options", make_function(&ReducedStatistics::GetOptions, return_value_policy<copy_const_reference>())) .def("Save", &ReducedStatistics::Save) @@ -55,7 +63,9 @@ void export_Reduced() .def("Create", &ReducedPotential::Create).staticmethod("Create") .def("Save", &ReducedPotential::Save) .def("GetEnergy", &ReducedPotential::GetEnergy) - .def("GetTotalEnergy", &ReducedPotential::GetTotalEnergy, + .def("GetTotalEnergy", get_h, + (arg("ent"), arg("normalize")=true)) + .def("GetTotalEnergy", get_v, (arg("ent"), arg("normalize")=true)) .add_property("options", make_function(&ReducedPotential::GetOptions, return_value_policy<copy_const_reference>())) diff --git a/modules/qa/src/impl/reduced_impl.cc b/modules/qa/src/impl/reduced_impl.cc index 8994fada4..a3778cb39 100644 --- a/modules/qa/src/impl/reduced_impl.cc +++ b/modules/qa/src/impl/reduced_impl.cc @@ -20,39 +20,65 @@ bool ReducedPotentialImpl::VisitResidue(const mol::ResidueHandle& res) if (!this->GetCAlphaCBetaPos(res, ca_pos_one, cb_pos_one)) { return false; } + + if (ent_) { + // we got a full entity handle. + mol::AtomHandleList within=ent_.FindWithin(ca_pos_one, + opts_.upper_cutoff-0.00001); + for (mol::AtomHandleList::const_iterator + i=within.begin(), e=within.end(); i!=e; ++i) { + if (i->GetName()=="CA") { + this->HandleResidue(aa_one, ca_pos_one, cb_pos_one, index, *i); + } - mol::AtomHandleList within=ent_.FindWithin(ca_pos_one, + } + return false; + } + mol::AtomViewList within=view_.FindWithin(ca_pos_one, opts_.upper_cutoff-0.00001); - for (mol::AtomHandleList::const_iterator + for (mol::AtomViewList::const_iterator i=within.begin(), e=within.end(); i!=e; ++i) { - mol::ResidueHandle res_two=i->GetResidue(); - if (res_two.GetIndex()-index<opts_.sequence_sep) { - continue; - } - if (i->GetName()!="CA" || !res_two.IsPeptideLinking()) { - continue; + if (i->GetName()=="CA") { + this->HandleResidue(aa_one, ca_pos_one, cb_pos_one, index, i->GetHandle()); } + } + return false; +} - AminoAcid aa_two=OneLetterCodeToAminoAcid(res_two.GetOneLetterCode()); - if (aa_two==Xxx) { - continue; - } - geom::Vec3 ca_pos_two; - geom::Vec3 cb_pos_two; - if (!this->GetCAlphaCBetaPos(res_two, ca_pos_two, cb_pos_two)) { - continue; - } +void ReducedPotentialImpl::HandleResidue(AminoAcid aa_one, + const geom::Vec3& ca_pos_one, + const geom::Vec3& cb_pos_one, + uint index_one, + const mol::AtomHandle& ca_two) +{ + + mol::ResidueHandle res_two=ca_two.GetResidue(); + if (res_two.GetIndex()-index_one<opts_.sequence_sep) { + return; + } + if (!res_two.IsPeptideLinking()) { + return; + } - Real dist=geom::Length(ca_pos_one-ca_pos_two); - if (dist<opts_.lower_cutoff) { - continue; - } - Real angle=geom::Angle(cb_pos_one-ca_pos_one, cb_pos_two-ca_pos_two); - this->OnInteraction(aa_one, aa_two, dist, angle); + AminoAcid aa_two=OneLetterCodeToAminoAcid(res_two.GetOneLetterCode()); + if (aa_two==Xxx) { + return; } - return false; + geom::Vec3 ca_pos_two; + geom::Vec3 cb_pos_two; + if (!this->GetCAlphaCBetaPos(res_two, ca_pos_two, cb_pos_two)) { + return; + } + + Real dist=geom::Length(ca_pos_one-ca_pos_two); + if (dist<opts_.lower_cutoff) { + return; + } + Real angle=geom::Angle(cb_pos_one-ca_pos_one, cb_pos_two-ca_pos_two); + + this->OnInteraction(aa_one, aa_two, dist, angle); } bool ReducedPotentialImpl::GetCAlphaCBetaPos(const mol::ResidueHandle& res, diff --git a/modules/qa/src/impl/reduced_impl.hh b/modules/qa/src/impl/reduced_impl.hh index 98651360b..43c312ea2 100644 --- a/modules/qa/src/impl/reduced_impl.hh +++ b/modules/qa/src/impl/reduced_impl.hh @@ -35,6 +35,10 @@ public: ReducedPotentialImpl(const ReducedStatOptions& opts, mol::EntityHandle ent): opts_(opts), ent_(ent) { } + + ReducedPotentialImpl(const ReducedStatOptions& opts, mol::EntityView ent): + opts_(opts), view_(ent) + { } virtual bool VisitResidue(const mol::ResidueHandle& res); @@ -44,9 +48,13 @@ public: virtual void OnInteraction(AminoAcid aa_one, AminoAcid aa_two, Real dist, Real angle)=0; -private: +private: + void HandleResidue(AminoAcid aa_one, const geom::Vec3& ca_pos_one, + const geom::Vec3& cb_pos_one, + uint index_one, const mol::AtomHandle& ca_two); ReducedStatOptions opts_; mol::EntityHandle ent_; + mol::EntityView view_; }; }}} diff --git a/modules/qa/src/reduced_potential.cc b/modules/qa/src/reduced_potential.cc index 78df5c902..1799c18e4 100644 --- a/modules/qa/src/reduced_potential.cc +++ b/modules/qa/src/reduced_potential.cc @@ -35,7 +35,11 @@ public: impl::ReducedPotentialImpl(opts, ent), energies_(energies), energy_(0.0), norm_(norm), count_(0) { } - + ReducedEnergiesCalc(const ReducedStatOptions& opts, ReducedEnergies& energies, + mol::EntityView ent, bool norm): + impl::ReducedPotentialImpl(opts, ent), + energies_(energies), energy_(0.0), norm_(norm), count_(0) + { } virtual void OnInteraction(AminoAcid aa_one, AminoAcid aa_two, Real dist, Real angle) { @@ -138,5 +142,12 @@ Real ReducedPotential::GetTotalEnergy(ost::mol::EntityHandle ent, bool norm) return calc.GetEnergy(); } +Real ReducedPotential::GetTotalEnergy(ost::mol::EntityView ent, bool norm) +{ + ReducedEnergiesCalc calc(opts_, energies_, ent, norm); + ent.Apply(calc); + return calc.GetEnergy(); +} + }} diff --git a/modules/qa/src/reduced_potential.hh b/modules/qa/src/reduced_potential.hh index 534892678..0f90debce 100644 --- a/modules/qa/src/reduced_potential.hh +++ b/modules/qa/src/reduced_potential.hh @@ -49,7 +49,7 @@ public: Real GetTotalEnergy(ost::mol::EntityHandle ent, bool normalize=true); - + Real GetTotalEnergy(ost::mol::EntityView ent, bool normalize=true); static bool GetCAlphaCBetaPos(const mol::ResidueHandle& res, geom::Vec3& ca_pos, geom::Vec3& cb_pos); diff --git a/modules/qa/src/reduced_statistics.cc b/modules/qa/src/reduced_statistics.cc index 9c20a9cd2..723fe29f9 100644 --- a/modules/qa/src/reduced_statistics.cc +++ b/modules/qa/src/reduced_statistics.cc @@ -41,6 +41,12 @@ public: histo_(histo) { } + ReducedStatExtractor(const ReducedStatOptions& opts, ReducedHistogram& histo, + mol::EntityView ent): + impl::ReducedPotentialImpl(opts, ent), + histo_(histo) + { } + virtual void OnInteraction(AminoAcid aa_one, AminoAcid aa_two, Real dist, Real angle) { @@ -130,4 +136,17 @@ void ReducedStatistics::Extract(mol::EntityHandle ent) ent.Apply(extractor); } +void ReducedStatistics::Extract(mol::EntityView ent) +{ + if (ent.GetChainCount()!=1) { + std::stringstream ss; + ss << "Expected exactly one chain, but entity has " + << ent.GetChainCount() << "chains"; + throw std::runtime_error(ss.str()); + } + ReducedStatExtractor extractor(opts_, histo_, ent); + ent.Apply(extractor); +} + + }} \ No newline at end of file diff --git a/modules/qa/src/reduced_statistics.hh b/modules/qa/src/reduced_statistics.hh index 419ac550f..7d1a24800 100644 --- a/modules/qa/src/reduced_statistics.hh +++ b/modules/qa/src/reduced_statistics.hh @@ -82,9 +82,12 @@ public: const ReducedStatOptions& GetOptions() const { return opts_; } - /// \brief extract the statistics from the given entity + /// \brief extract the statistics from the given entity (handle) void Extract(mol::EntityHandle ent); + /// \brief extract statistics from given entity (view) + void Extract(mol::EntityView ent); + void Save(const String& filename); -- GitLab