diff --git a/core/src/tetrahedral_polytope.hh b/core/src/tetrahedral_polytope.hh index bb12cab98c8cae3017d46c748b94cdbe44a3b6ec..c1cf892fe3b7cfdaf687ca82e23de7d356f720aa 100644 --- a/core/src/tetrahedral_polytope.hh +++ b/core/src/tetrahedral_polytope.hh @@ -8,10 +8,10 @@ namespace promod3{ namespace core{ -// Tetrahedral polytope for rapid collision detection. -// Initialize with a position and a collision distance (radius) +// Tetrahedral polytope for collision detection. +// Initialize it with a position and a collision distance (radius) // The Overlap function guarantees to return true if two polytopes are closer -// than the sum of their radii but the it might already return true if they're +// than the sum of their radii but it might already return true if they're // a bit further away... struct TetrahedralPolytope{ diff --git a/doc/tests/scripts/sidechain_steps.py b/doc/tests/scripts/sidechain_steps.py index 06b44e96c60de4599d5dad64c01367fc6e84c936..be9f5c7558e9e14d9f90c0507a0766933ed5daaa 100644 --- a/doc/tests/scripts/sidechain_steps.py +++ b/doc/tests/scripts/sidechain_steps.py @@ -55,7 +55,7 @@ for i,r in enumerate(prot.residues): r, rotamer_ids[i], i, library, torsion_angles[i][0], torsion_angles[i][1]) - frame.SetFrameEnergy(rot_group) + rot_group.SetFrameEnergy(frame) # remove super unlikely rotamer in rotamer group # e.g. those who clash with the frame rot_group.ApplySelfEnergyThresh() diff --git a/sidechain/pymod/_reconstruct_sidechains.py b/sidechain/pymod/_reconstruct_sidechains.py index 3c15dd49e454f887561aa3bea4af40b9d284819b..d402500c1f650893e0c4ce4ed8345a335b7cead4 100644 --- a/sidechain/pymod/_reconstruct_sidechains.py +++ b/sidechain/pymod/_reconstruct_sidechains.py @@ -289,7 +289,7 @@ def _GetRotamerGroups(res_list, rot_ids, indices, rot_lib, rot_constructor, except: continue # keep best ones - frame.SetFrameEnergy(rot_group) + rot_group.SetFrameEnergy(frame) rot_group.ApplySelfEnergyThresh() rotamer_groups.append(rot_group) residues_with_rotamer_group.append(i) @@ -324,7 +324,7 @@ def _GetDisulfidBridges(frame_residues, cystein_indices, res_list, rotamer_libra rot_group = _GetRotamerGroup(r.handle, sidechain.CYD, i, rotamer_library, rotamer_constructor, phi_angles[i], psi_angles[i], use_frm, bbdep) - frame.AddFrameEnergy(rot_group) + rot_group.AddFrameEnergy(frame) cystein_rotamers.append(rot_group) # get CYS with disulfid bonds and the chosen rotamers diff --git a/sidechain/pymod/export_frame.cc b/sidechain/pymod/export_frame.cc index 7cff9995e8fafd06688b2d39adcbb66391ec0274..f9b1e1389164bb150f1d46555c5e41220e8f781d 100644 --- a/sidechain/pymod/export_frame.cc +++ b/sidechain/pymod/export_frame.cc @@ -32,32 +32,12 @@ ParticlePtr WrapGetItem(FrameResiduePtr p, uint index){ return ret_particle; } -void SetFrameEnergyOne(FramePtr frame, RRMRotamerGroupPtr p){ - frame->SetFrameEnergy(p); -} - -void SetFrameEnergyTwo(FramePtr frame, FRMRotamerGroupPtr p){ - frame->SetFrameEnergy(p); -} - -void AddFrameEnergyOne(FramePtr frame, RRMRotamerGroupPtr p){ - frame->AddFrameEnergy(p); -} - -void AddFrameEnergyTwo(FramePtr frame, FRMRotamerGroupPtr p){ - frame->AddFrameEnergy(p); -} - } void export_Frame() { class_<Frame, boost::noncopyable>("Frame", no_init) .def("__init__",make_constructor(&WrapFrameInit)) - .def("SetFrameEnergy",&SetFrameEnergyOne,(arg("rrm_rotamer_group"))) - .def("SetFrameEnergy",&SetFrameEnergyTwo,(arg("frm_rotamer_group"))) - .def("AddFrameEnergy",&AddFrameEnergyOne,(arg("rrm_rotamer_group"))) - .def("AddFrameEnergy",&AddFrameEnergyTwo,(arg("frm_rotamer_group"))) ; register_ptr_to_python<FramePtr>(); diff --git a/sidechain/pymod/export_rotamer.cc b/sidechain/pymod/export_rotamer.cc index 4d88fa70e6eca5304cbc105fd2f1867db969efe1..da649f0e9defd2ab150c9235c8f879110f0f8de2 100644 --- a/sidechain/pymod/export_rotamer.cc +++ b/sidechain/pymod/export_rotamer.cc @@ -194,6 +194,8 @@ void export_Rotamer() .def("ApplyOnResidue", WrapRRMGApplyOnResAA, (arg("index"), arg("all_atom"), arg("res_idx"))) .def("Merge",&RRMRotamerGroup::Merge,(arg("other"))) + .def("SetFrameEnergy",&RRMRotamerGroup::SetFrameEnergy,(arg("frame"))) + .def("AddFrameEnergy",&RRMRotamerGroup::AddFrameEnergy,(arg("frame"))) .def("ApplySelfEnergyThresh", &RRMRotamerGroup::ApplySelfEnergyThresh, (arg("thresh")=30)) ; @@ -211,6 +213,8 @@ void export_Rotamer() .def("ApplyOnResidue", WrapFRMGApplyOnResAA, (arg("index"), arg("all_atom"), arg("res_idx"))) .def("Merge",&FRMRotamerGroup::Merge,(arg("other"))) + .def("SetFrameEnergy",&FRMRotamerGroup::SetFrameEnergy,(arg("frame"))) + .def("AddFrameEnergy",&FRMRotamerGroup::AddFrameEnergy,(arg("frame"))) .def("ApplySelfEnergyThresh", &FRMRotamerGroup::ApplySelfEnergyThresh, (arg("thresh")=30)) ; diff --git a/sidechain/src/frame.cc b/sidechain/src/frame.cc index 57d798e26ec1111caa2d4a0f2e0b3df3d09af7c1..52a526216e00f2e666177976d2972409b7b60595 100644 --- a/sidechain/src/frame.cc +++ b/sidechain/src/frame.cc @@ -1,38 +1,6 @@ #include <promod3/sidechain/frame.hh> #include <promod3/core/runtime_profiling.hh> - -namespace{ - -Real LogSumExp(double* values, int n){ - - double min = std::numeric_limits<double>::max(); - double max = -std::numeric_limits<double>::max(); - - for(int i = 0; i < n; ++i){ - min = std::min(min, values[i]); - max = std::max(max, values[i]); - } - - double range = max - min; - - if(range > 100){ - // The log sum exp function (LSE) is bound by - // max(x1, x2, ..., xn) <= LSE(x1, x2, ..., xn) <= - // max(x1, x2, ..., xn) + log(n) - // the maximum value is therefore a rather good approximation - // this gets returned if the range get's too large because - // even the LSE trick implemented below gets numerical problems... - return max; - } - else{ - double a = (max + min) / 2; - double sum = 0.0; - for(int i = 0; i < n; ++i) sum += std::exp(values[i]-a); - return std::log(sum) + a; - } -} - -} +#include <ost/mol/xcs_editor.hh> namespace promod3{ namespace sidechain{ @@ -87,136 +55,4 @@ Frame::Frame(const std::vector<FrameResiduePtr>& frame_residues): collision_tree_.ResetTree(particle_positions, collision_distances); } -void Frame::SetFrameEnergy(RRMRotamerGroupPtr p) const{ - - for(RRMRotamerGroup::iterator i = p->begin(); i != p->end(); ++i){ - (*i)->SetFrameEnergy(0.0); - } - this->AddFrameEnergy(p); -} - -void Frame::AddFrameEnergy(RRMRotamerGroupPtr p) const{ - - core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( - "Frame::AddFrameEnergy", 2); - - overlapping_particles_.clear(); - collision_tree_.OverlappingPolytopes(p->GetCollisionTree(), - overlapping_particles_); - - const std::vector<int>& p_rotamer_indices = p->GetRotamerIndices(); - const std::vector<Particle*>& p_particles = p->GetParticles(); - uint p_residue_index = p->GetResidueIndex(); - - Real frame_energies[p->size()]; - memset(frame_energies,0,p->size()*sizeof(Real)); - Real energy; - - for(std::vector<std::pair<int,int> >::iterator i = overlapping_particles_.begin(), - e = overlapping_particles_.end(); i != e; ++i){ - - if(residue_indices_[i->first] == p_residue_index) continue; - - energy = frame_particles_[i->first]->PairwiseEnergy(p_particles[i->second]); - frame_energies[p_rotamer_indices[i->second]] += energy; - } - - Real* e_ptr = frame_energies; - for(RRMRotamerGroup::iterator i = p->begin(), e = p->end(); i != e; ++i){ - (*i)->AddFrameEnergy(*e_ptr); - ++e_ptr; - } -} - -void Frame::SetFrameEnergy(FRMRotamerGroupPtr p) const{ - - for(FRMRotamerGroup::iterator i = p->begin(); i != p->end(); ++i){ - (*i)->SetFrameEnergy(0); - for(uint j = 0, e = (*i)->subrotamer_size(); j < e; ++j){ - (*i)->SetFrameEnergy(0,j); - } - } - this->AddFrameEnergy(p); -} - -void Frame::AddFrameEnergy(FRMRotamerGroupPtr p) const{ - - core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( - "Frame::AddFrameEnergy", 2); - - overlapping_particles_.clear(); - collision_tree_.OverlappingPolytopes(p->GetCollisionTree(), - overlapping_particles_); - - const std::vector<int>& p_rotamer_indices = p->GetRotamerIndices(); - const std::vector<FRMRotamer::subrotamer_association_iterator>& p_ass_ind_b = p->GetSubrotamerAssociationIndB(); - const std::vector<FRMRotamer::subrotamer_association_iterator>& p_ass_ind_e = p->GetSubrotamerAssociationIndE(); - const std::vector<Particle*>& p_particles = p->GetParticles(); - uint p_residue_index = p->GetResidueIndex(); - - - //reserve the memory - Real** frame_energies = new Real*[p->size()]; - int counter = 0; - int subrotamer_size; - int max_subrotamer_size = 0; - for(FRMRotamerGroup::iterator i = p->begin(); i != p->end(); ++i){ - subrotamer_size = (*i)->subrotamer_size(); - frame_energies[counter] = new Real[subrotamer_size]; - //set the subrotamer frame energies to the already existing - //energies - for(int j = 0; j < subrotamer_size; ++j){ - frame_energies[counter][j] = (*i)->GetFrameEnergy(j); - } - max_subrotamer_size = std::max(max_subrotamer_size,subrotamer_size); - ++counter; - } - - Real energy; - Real* e_ptr; - for(std::vector<std::pair<int,int> >::iterator i = overlapping_particles_.begin(), - e = overlapping_particles_.end(); i != e ; ++i){ - - - if(residue_indices_[i->first] == p_residue_index){ - continue; - } - - energy = frame_particles_[i->first]->PairwiseEnergy(p_particles[i->second]); - if(energy == 0) continue; - - e_ptr = frame_energies[p_rotamer_indices[i->second]]; - - for(FRMRotamer::subrotamer_association_iterator j = p_ass_ind_b[i->second]; - j != p_ass_ind_e[i->second]; ++j){ - e_ptr[*j] += energy; - } - } - - double exp_buffer[max_subrotamer_size]; - counter = 0; - Real T; - Real one_over_T; - int actual_subrotamer_size; - - for(FRMRotamerGroup::iterator i = p->begin(), e = p->end(); i != e; ++i){ - actual_subrotamer_size = (*i)->subrotamer_size(); - e_ptr = frame_energies[counter]; - T = (*i)->GetTemperature(); - one_over_T = 1.0/T; - for(int j = 0; j < actual_subrotamer_size; ++j){ - (*i)->SetFrameEnergy(e_ptr[j],j); - exp_buffer[j] = -e_ptr[j]*one_over_T; - - } - (*i)->SetFrameEnergy(-T*LogSumExp(exp_buffer,actual_subrotamer_size)); - ++counter; - } - - for(uint i = 0; i < p->size(); ++i){ - delete [] frame_energies[i]; - } - delete [] frame_energies; -} - }}//ns diff --git a/sidechain/src/frame.hh b/sidechain/src/frame.hh index aca991606b7614e7f144388481eb0f11ea50e825..abd444060c4a4eb4073a9f579e597be297cfacd2 100644 --- a/sidechain/src/frame.hh +++ b/sidechain/src/frame.hh @@ -9,7 +9,8 @@ #include <promod3/core/tetrahedral_polytope.hh> #include <promod3/sidechain/particle.hh> #include <promod3/core/message.hh> -#include <promod3/sidechain/rotamer_group.hh> + +#include <ost/mol/residue_handle.hh> #include <boost/shared_ptr.hpp> @@ -58,25 +59,24 @@ public: Frame(const std::vector<FrameResiduePtr>& frame_residues); - void SetFrameEnergy(RRMRotamerGroupPtr p) const; - - void SetFrameEnergy(FRMRotamerGroupPtr p) const; - - void AddFrameEnergy(RRMRotamerGroupPtr p) const; - - void AddFrameEnergy(FRMRotamerGroupPtr p) const; + const std::vector<uint>& GetResidueIndices() const { + return residue_indices_; + } const promod3::core::TetrahedralPolytopeTree& GetCollisionTree() const { return collision_tree_; } + const std::vector<Particle*>& GetParticles() const { + return frame_particles_; + } + private: std::vector<Particle*> frame_particles_; std::vector<FrameResiduePtr> frame_residues_; std::vector<uint> residue_indices_; promod3::core::TetrahedralPolytopeTree collision_tree_; - mutable std::vector<std::pair<int,int> > overlapping_particles_; }; diff --git a/sidechain/src/rotamer_group.cc b/sidechain/src/rotamer_group.cc index e4038ba2698b495ff668a64fba174f3f7714f871..a106aa5b63efe533df6c896bddb86f75ca06d6d5 100644 --- a/sidechain/src/rotamer_group.cc +++ b/sidechain/src/rotamer_group.cc @@ -72,6 +72,41 @@ void RRMRotamerGroup::ApplyOnResidue(uint rotamer_index, rotamers_[rotamer_index]->ApplyOnResidue(all_atom, res_idx); } +void RRMRotamerGroup::SetFrameEnergy(FramePtr frame) { + + for(iterator i = this->begin(); i != this->end(); ++i){ + (*i)->SetFrameEnergy(0.0); + } + this->AddFrameEnergy(frame); +} + +void RRMRotamerGroup::AddFrameEnergy(FramePtr frame) { + + core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( + "RRMRotamerGroup::AddFrameEnergy", 2); + + std::vector<std::pair<int, int> > overlapping_particles; + collision_tree_.OverlappingPolytopes(frame->GetCollisionTree(), + overlapping_particles); + + const std::vector<uint>& frame_residue_indices = frame->GetResidueIndices(); + const std::vector<Particle*>& frame_particles = frame->GetParticles(); + std::vector<Real> frame_energies(this->size(), 0.0); + Real energy; + + for(std::vector<std::pair<int,int> >::iterator i = overlapping_particles.begin(), + e = overlapping_particles.end(); i != e; ++i){ + if(frame_residue_indices[i->second] != residue_index_){ + energy = particles_[i->first]->PairwiseEnergy(frame_particles[i->second]); + frame_energies[rotamer_indices_[i->first]] += energy; + } + } + + for(uint i = 0; i < this->size(); ++i){ + rotamers_[i]->AddFrameEnergy(frame_energies[i]); + } +} + Real* RRMRotamerGroup::CalculatePairwiseEnergies(RRMRotamerGroupPtr other, Real threshold) { @@ -91,19 +126,21 @@ Real* RRMRotamerGroup::CalculatePairwiseEnergies(RRMRotamerGroupPtr other, int this_size = this->size(); int other_size = other->size(); Real* data = new Real[this_size * other_size]; - memset(data,0,this_size*other_size*sizeof(Real)); + memset(data, 0, this_size*other_size*sizeof(Real)); Real energy; for(std::vector<std::pair<int,int> >::iterator i = overlapping_particles.begin(), e = overlapping_particles.end(); i!= e; ++i){ energy = particles_[i->first]->PairwiseEnergy(other->particles_[i->second]); - if(energy == 0) continue; - data[rotamer_indices_[i->first]*other_size+other->rotamer_indices_[i->second]] += energy; + if(energy != Real(0.0)){ + data[rotamer_indices_[i->first]*other_size + + other->rotamer_indices_[i->second]] += energy; + } } Real max_e = 0.0; for(int i = 0; i < this_size*other_size; ++i){ - max_e = (max_e<data[i])?data[i]:max_e; + max_e = (max_e < data[i]) ? data[i] : max_e; } if(max_e < threshold){ @@ -234,6 +271,92 @@ void FRMRotamerGroup::Merge(const FRMRotamerGroupPtr other){ this->Init(); } +void FRMRotamerGroup::SetFrameEnergy(FramePtr frame){ + + for(iterator i = this->begin(); i != this->end(); ++i){ + (*i)->SetFrameEnergy(0); + for(uint j = 0, e = (*i)->subrotamer_size(); j < e; ++j){ + (*i)->SetFrameEnergy(0,j); + } + } + this->AddFrameEnergy(frame); +} + +void FRMRotamerGroup::AddFrameEnergy(FramePtr frame){ + + core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( + "FRMRotamerGroup::AddFrameEnergy", 2); + + std::vector<std::pair<int, int> > overlapping_particles; + collision_tree_.OverlappingPolytopes(frame->GetCollisionTree(), + overlapping_particles); + + const std::vector<uint>& frame_residue_indices = frame->GetResidueIndices(); + const std::vector<Particle*>& frame_particles = frame->GetParticles(); + + //reserve the memory + Real** frame_energies = new Real*[this->size()]; + int counter = 0; + int subrotamer_size; + int max_subrotamer_size = 0; + for(iterator i = this->begin(); i != this->end(); ++i){ + subrotamer_size = (*i)->subrotamer_size(); + frame_energies[counter] = new Real[subrotamer_size]; + //set the subrotamer frame energies to the already existing + //energies + for(int j = 0; j < subrotamer_size; ++j){ + frame_energies[counter][j] = (*i)->GetFrameEnergy(j); + } + max_subrotamer_size = std::max(max_subrotamer_size, subrotamer_size); + ++counter; + } + + Real energy; + Real* e_ptr; + for(std::vector<std::pair<int,int> >::iterator i = overlapping_particles.begin(), + e = overlapping_particles.end(); i != e ; ++i){ + + if(frame_residue_indices[i->second] != residue_index_){ + energy = particles_[i->first]->PairwiseEnergy(frame_particles[i->second]); + + if(energy != Real(0.0)){ + + e_ptr = frame_energies[rotamer_indices_[i->first]]; + + for(FRMRotamer::subrotamer_association_iterator j = ass_ind_b_[i->first]; + j != ass_ind_e_[i->first]; ++j){ + e_ptr[*j] += energy; + } + } + } + } + + std::vector<double> exp_buffer(max_subrotamer_size, 0.0); + counter = 0; + Real T; + Real one_over_T; + int actual_subrotamer_size; + + for(iterator i = this->begin(), e = this->end(); i != e; ++i){ + actual_subrotamer_size = (*i)->subrotamer_size(); + e_ptr = frame_energies[counter]; + T = (*i)->GetTemperature(); + one_over_T = 1.0/T; + for(int j = 0; j < actual_subrotamer_size; ++j){ + (*i)->SetFrameEnergy(e_ptr[j],j); + exp_buffer[j] = -e_ptr[j]*one_over_T; + } + (*i)->SetFrameEnergy(-T*LogSumExp(&exp_buffer[0], actual_subrotamer_size)); + ++counter; + } + + for(uint i = 0; i < this->size(); ++i){ + delete [] frame_energies[i]; + } + delete [] frame_energies; + +} + Real* FRMRotamerGroup::CalculatePairwiseEnergies(FRMRotamerGroupPtr other, Real threshold) { @@ -294,24 +417,25 @@ Real* FRMRotamerGroup::CalculatePairwiseEnergies(FRMRotamerGroupPtr other, energy = particles_[i->first]->PairwiseEnergy(other->particles_[i->second]); - if(energy == 0.0) continue; + if(energy != Real(0.0)){ - r_index_j = other->rotamer_indices_[i->second]; + r_index_j = other->rotamer_indices_[i->second]; - other_subrotamer_size = other_subrotamer_sizes[r_index_j]; + other_subrotamer_size = other_subrotamer_sizes[r_index_j]; - j_b = ass_ind_b_[i->first]; - j_e = ass_ind_e_[i->first]; - temp_it = other->ass_ind_b_[i->second]; - k_e = other->ass_ind_e_[i->second]; + j_b = ass_ind_b_[i->first]; + j_e = ass_ind_e_[i->first]; + temp_it = other->ass_ind_b_[i->second]; + k_e = other->ass_ind_e_[i->second]; - e_ptr = subrotamer_energies[rotamer_indices_[i->first]*other_size+r_index_j]; + e_ptr = subrotamer_energies[rotamer_indices_[i->first]*other_size+r_index_j]; - for(;j_b != j_e; ++j_b){ - k_b = temp_it; - temp = (*j_b)*other_subrotamer_size; - for(;k_b != k_e; ++k_b){ - e_ptr[temp+(*k_b)] += energy; + for(;j_b != j_e; ++j_b){ + k_b = temp_it; + temp = (*j_b)*other_subrotamer_size; + for(;k_b != k_e; ++k_b){ + e_ptr[temp+(*k_b)] += energy; + } } } } @@ -332,7 +456,7 @@ Real* FRMRotamerGroup::CalculatePairwiseEnergies(FRMRotamerGroupPtr other, return NULL; } - double exponent_buffer[max_subrotamer_combinations]; + std::vector<double> exp_buffer(max_subrotamer_combinations, 0.0); Real* energy_ptr; Real* data_ptr = &data[0]; double* double_ptr; @@ -370,7 +494,7 @@ Real* FRMRotamerGroup::CalculatePairwiseEnergies(FRMRotamerGroupPtr other, energy_ptr = subrotamer_energies[temp+j]; T = 0.5 * (T_i + j_temperatures[j]); one_over_T = 1.0/T; - double_ptr = &exponent_buffer[0]; + double_ptr = &exp_buffer[0]; b = other_subrotamer_sizes[j]; for(int k = 0; k < a; ++k){ @@ -382,7 +506,7 @@ Real* FRMRotamerGroup::CalculatePairwiseEnergies(FRMRotamerGroupPtr other, } } - *data_ptr = -T*LogSumExp(exponent_buffer,a*b) - i_frame_energy - j_frame_energies[j]; + *data_ptr = -T*LogSumExp(&exp_buffer[0], a*b) - i_frame_energy - j_frame_energies[j]; ++data_ptr; Real_ptr+=b; } diff --git a/sidechain/src/rotamer_group.hh b/sidechain/src/rotamer_group.hh index 2a48f9b56fecf7c552d73fc6fdcae3f7fe3a9cd8..5163b4dd09f907d93d3485579a22c35014dd1abe 100644 --- a/sidechain/src/rotamer_group.hh +++ b/sidechain/src/rotamer_group.hh @@ -9,6 +9,7 @@ #include <promod3/core/tetrahedral_polytope.hh> #include <promod3/sidechain/rotamer.hh> +#include <promod3/sidechain/frame.hh> #include <promod3/core/message.hh> #include <ost/mol/residue_handle.hh> @@ -37,6 +38,10 @@ public: void ApplyOnResidue(uint rotamer_index, loop::AllAtomPositions& all_atom, uint res_idx) const; + void SetFrameEnergy(FramePtr p); + + void AddFrameEnergy(FramePtr p); + Real* CalculatePairwiseEnergies(RRMRotamerGroupPtr other, Real threshold); Real GetMaxP() const; @@ -70,8 +75,8 @@ public: const_iterator begin() const { return rotamers_.begin(); } const_iterator end() const { return rotamers_.end(); } - const std::vector<int>& GetRotamerIndices() { return rotamer_indices_; } - const std::vector<Particle*>& GetParticles() { return particles_; } + //const std::vector<int>& GetRotamerIndices() { return rotamer_indices_; } + //const std::vector<Particle*>& GetParticles() { return particles_; } const std::vector<RRMRotamerPtr>& GetRotamers() { return rotamers_; } private: @@ -103,6 +108,10 @@ public: void ApplyOnResidue(uint rotamer_index, loop::AllAtomPositions& all_atom, uint res_idx) const; + void SetFrameEnergy(FramePtr p); + + void AddFrameEnergy(FramePtr p); + Real* CalculatePairwiseEnergies(FRMRotamerGroupPtr other, Real threshold); Real GetMaxP() const; @@ -136,13 +145,7 @@ public: const_iterator begin() const { return rotamers_.begin(); } const_iterator end() const { return rotamers_.end(); } - const std::vector<int>& GetRotamerIndices() { return rotamer_indices_; } - const std::vector<Particle*>& GetParticles() { return particles_; } const std::vector<FRMRotamerPtr>& GetRotamers() { return rotamers_; } - const std::vector<FRMRotamer::subrotamer_association_iterator>& - GetSubrotamerAssociationIndB() { return ass_ind_b_; } - const std::vector<FRMRotamer::subrotamer_association_iterator>& - GetSubrotamerAssociationIndE() { return ass_ind_e_; } private: diff --git a/sidechain/src/sidechain_reconstructor.cc b/sidechain/src/sidechain_reconstructor.cc index 7a8f72579506eb8511142178d1449e656c94e6f0..eb91586567f8e638d7c9baff82a2b708921ea5a2 100644 --- a/sidechain/src/sidechain_reconstructor.cc +++ b/sidechain/src/sidechain_reconstructor.cc @@ -155,12 +155,12 @@ template<typename RotamerGroup> void SidechainReconstructor::BuildDisulfids_( loop::AllAtomPositions& out_pos = *(res->env_pos->all_pos); // collect rotamers and assign frame energies with what we know so far - Frame frame(frame_residues); + FramePtr frame(new Frame(frame_residues)); std::vector<RotamerGroupPtr> cys_rotamer_groups(cys_indices.size()); for (uint i = 0; i < cys_indices.size(); ++i) { const uint res_idx = res_indices[cys_indices[i]]; env_->GetCydRotamerGroup(cys_rotamer_groups[i], res_idx); - if (cys_rotamer_groups[i]) frame.SetFrameEnergy(cys_rotamer_groups[i]); + if (cys_rotamer_groups[i]) cys_rotamer_groups[i]->SetFrameEnergy(frame); } // look for bridges (cys_rotamers[i] only non-null if bridge found) @@ -307,9 +307,9 @@ void SidechainReconstructor::SolveSystem_( CollectRotamerGroups_(res, rotamer_groups, has_sidechain); // set frame energies in rotamers (note: internal energies are precomputed) - Frame frame(frame_residues); + FramePtr frame(new Frame(frame_residues)); for (uint i = 0; i < rotamer_groups.size(); ++i) { - frame.SetFrameEnergy(rotamer_groups[i]); + rotamer_groups[i]->SetFrameEnergy(frame); rotamer_groups[i]->ApplySelfEnergyThresh(); }