implement ReducedPotential and ReducedStatistics for views

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>()))

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, 

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_;
 };
 
 }}}

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();
+}
+
 
 }}

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);

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

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);