diff --git a/modules/qa/src/packing_potential.cc b/modules/qa/src/packing_potential.cc index e7a473fd3ac496291c2c459da974ae930fec9548..10a7de4aeb20f1642f743915bc3cab7719d7062c 100644 --- a/modules/qa/src/packing_potential.cc +++ b/modules/qa/src/packing_potential.cc @@ -187,10 +187,7 @@ bool PackingPotential::VisitAtom(const mol::AtomHandle& atom) e=views_.end(); i!=e; ++i) { count+=i->FindWithin(atom.GetPos(), options_.cutoff).size(); } - - count=count>options_.max_counts ? options_.max_counts : count; - - energy_+=energies_.Get(aa, count/options_.bucket_size); + energy_+=this->GetPackingEnergy(aa, count); energy_counts_++; return false; } diff --git a/modules/qa/src/packing_potential.hh b/modules/qa/src/packing_potential.hh index df49d006c4d7a3774b1646cce3f67496243c0afc..6d11d2e6c6ff43535e725c1f4267910edbc4a4cf 100644 --- a/modules/qa/src/packing_potential.hh +++ b/modules/qa/src/packing_potential.hh @@ -88,6 +88,13 @@ public: int GetEnergyCounts(); + + Real GetPackingEnergy(AminoAcid aa, int count) const + { + count=count>options_.max_counts ? options_.max_counts : count; + return energies_.Get(aa, count/options_.bucket_size); + } + template <typename DS> void Serialize(DS& ds); diff --git a/modules/qa/src/torsion_potential.cc b/modules/qa/src/torsion_potential.cc index 651718a4f04bc0d53faf6fadf3dd1ab8eb328e0c..2ecf8c655f24de49220a00be1b09b24d1e14b856 100644 --- a/modules/qa/src/torsion_potential.cc +++ b/modules/qa/src/torsion_potential.cc @@ -16,6 +16,7 @@ // along with this library; if not, write to the Free Software Foundation, Inc., // 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA //------------------------------------------------------------------------------ + #include "torsion_potential.hh" #include "amino_acids.hh" @@ -39,10 +40,9 @@ namespace { class TorsionEnergyCalc : public mol::EntityVisitor { public: - TorsionEnergyCalc(TorsionPotential::TorsionEnergies& energies, + TorsionEnergyCalc(const TorsionPotentialPtr& pot, TorsionPotentialOpts opts): - energies_(energies), energy_(0.0), - opts_(opts), num_torsions_(0) + pot_(pot), energy_(0.0), num_torsions_(0) { } @@ -75,24 +75,16 @@ public: Real next_phi=next_phit.GetAngle()*180/M_PI - 0.001; if (next_phi<-180) next_phi=-180; Real next_psi=next_psit.GetAngle()*180/M_PI - 0.001; if (next_psi<-180) next_psi=-180; + // calculate position of the amino acid in the alphabet - int icenter=this->GetAAIndex(ca); - - if (icenter >-1) { - energy_+=energies_.Get(icenter, prev_phi, prev_psi, - central_phi, central_psi, - next_phi, next_psi); - num_torsions_++; - } - else { - LOG_INFO("Amino acid not found in alphabets..."); - } + energy_+=pot_->GetTorsionEnergy(ca, prev_phi, prev_psi, central_phi, + central_psi, next_phi, next_psi); + ++num_torsions_; return false; } Real GetEnergy() const { return energy_; - //return num_torsions_==0?0:energy_/num_torsions_; } int GetEnergyCounts() const { @@ -100,29 +92,11 @@ public: } private: - int GetAAIndex(AminoAcid aa) - { - return this->GetIndex(aa, opts_.alphabet); - } - - - int GetIndex(AminoAcid aa, AminoAcidAlphabet& alpha) { - AminoAcidAlphabet::iterator i=alpha.begin(); - for (int index=0; i!=alpha.end(); ++i, ++index) { - if ((*i).Contains(aa)) { - return index; - } - - } - return -1; - } -private: - TorsionPotential::TorsionEnergies& energies_; + TorsionPotentialPtr pot_; AminoAcid prev_; AminoAcid center_; mol::ResidueHandle cr_; Real energy_; - TorsionPotentialOpts opts_; int num_torsions_; }; @@ -153,6 +127,18 @@ TorsionPotentialOpts::TorsionPotentialOpts(): } +int TorsionPotential::GetAAIndex(AminoAcid aa) const +{ + AminoAcidAlphabet::const_iterator i=options_.alphabet.begin(); + for (int index=0; i!=options_.alphabet.end(); ++i, ++index) { + if ((*i).Contains(aa)) { + return index; + } + } + return -1; +} + + TorsionPotentialPtr TorsionPotential::Create(TorsionStatisticsPtr statistics, const TorsionPotentialOpts& opts, bool calculate_average_energy_flag) @@ -225,7 +211,7 @@ void TorsionPotential::SaveToFile(const String& path) Real TorsionPotential::GetTotalEnergy(mol::EntityHandle entity) { - TorsionEnergyCalc c(energies_, options_); + TorsionEnergyCalc c(this->shared_from_this(), options_); entity.Apply(c); num_torsions_ = c.GetEnergyCounts(); return c.GetEnergy(); @@ -233,7 +219,7 @@ Real TorsionPotential::GetTotalEnergy(mol::EntityHandle entity) Real TorsionPotential::GetTotalEnergy(mol::EntityView entity) { - TorsionEnergyCalc c(energies_, options_); + TorsionEnergyCalc c(this->shared_from_this(), options_); entity.Apply(c); num_torsions_ = c.GetEnergyCounts(); return c.GetEnergy(); diff --git a/modules/qa/src/torsion_potential.hh b/modules/qa/src/torsion_potential.hh index 68faeb24fc2ab71cc4c261edf6d719452dd25201..d0048e84ccce0ada40901e93dd70fea536d57d04 100644 --- a/modules/qa/src/torsion_potential.hh +++ b/modules/qa/src/torsion_potential.hh @@ -30,6 +30,7 @@ */ #include <ost/qa/torsion_statistics.hh> #include <boost/shared_ptr.hpp> +#include <boost/enable_shared_from_this.hpp> #include <vector> #include <boost/scoped_ptr.hpp> @@ -67,7 +68,8 @@ typedef boost::shared_ptr<TorsionPotential> TorsionPotentialPtr; /// \brief Torsion potential /// /// The torsion potential class is parametrisable by TorsionPotentialOpts. -class DLLEXPORT_OST_QA TorsionPotential { +class DLLEXPORT_OST_QA TorsionPotential : + public boost::enable_shared_from_this<TorsionPotential> { public: /// \brief create new torsion potential with the given torsion statistics /// and options @@ -89,6 +91,18 @@ public: /// \brief retrieve total number of energy local (i.e. valid residues) int GetEnergyCounts() const; + + + + Real GetTorsionEnergy(AminoAcid central_aa, Real prev_phi, Real prev_psi, + Real central_phi, Real central_psi, + Real next_phi, Real next_psi) const { + int icenter=this->GetAAIndex(central_aa); + return energies_.Get(icenter, prev_phi, prev_psi, + central_phi, central_psi, + next_phi, next_psi); + } + /// \brief save torsion potential /// @@ -106,6 +120,9 @@ public: typedef MultiClassifier<float, int, Real, Real, Real, Real, Real, Real> TorsionEnergies; private: + + int GetAAIndex(AminoAcid aa) const; + void Fill(const TorsionStatisticsPtr& stat, bool calculate_average_energy_flag);