diff --git a/modules/mol/alg/doc/molalg.rst b/modules/mol/alg/doc/molalg.rst index 795d4f96e2613fcb50c66abcde51ffe91075f724..51de4cca63a8b7f794d4a671ca5cff1b779d7cdb 100644 --- a/modules/mol/alg/doc/molalg.rst +++ b/modules/mol/alg/doc/molalg.rst @@ -685,18 +685,20 @@ Local Distance Test scores (lDDT, DRMSD) Dictionary-like object containing the list of interatomic distances that originate from a single residue to be checked during a run of the Local Distance Difference Test algorithm - (key = pair of :class:`UniqueAtomIdentifier`, value = pair of floats) + (key = pair of :class:`UniqueAtomIdentifier`, value = pair of floats + representing min and max distance observed in the structures used to build + the map). .. class:: GlobalRDMap - Dictionary-like object containing all the :class:`~ost.mol.alg.ResidueRDMap` objects related to residues of a single structure - (key = :class:`~ost.mol.ResNum`, value = :class:`ResidueRDMap`) + Dictionary-like object containing all the :class:`~ost.mol.alg.ResidueRDMap` objects related to all the residues + (key = :class:`~ost.mol.ResNum`, value = :class:`ResidueRDMap`). .. function:: PrintResidueRDMap(residue_distance_list) Prints to standard output all the distances contained in a - :class:`~ost.mol.alg.ResidueRDMap` object + :class:`~ost.mol.alg.ResidueRDMap` object. .. function:: PrintGlobalRDMap(global_distance_list) diff --git a/modules/mol/alg/src/filter_clashes.cc b/modules/mol/alg/src/filter_clashes.cc index e0d2bee05a80d0959b3528776d501151f2b388c0..8623d8779342b089366d6d43bf835cdc1c51377b 100644 --- a/modules/mol/alg/src/filter_clashes.cc +++ b/modules/mol/alg/src/filter_clashes.cc @@ -404,19 +404,28 @@ std::pair<EntityView,StereoChemistryInfo> CheckStereoChemistry(const EntityView& } if (other_atom.GetHashCode() > atom.GetHandle().GetHashCode()) { Real blength = bond.GetLength(); - String bond_str = bond_string(atom,other_atom); - std::pair<Real,Real> length_stddev = bond_table.GetParam(bond_str,residue_str); + String bond_str = bond_string(atom, other_atom); + std::pair<Real,Real> length_stddev = bond_table.GetParam(bond_str, residue_str); Real ref_length = length_stddev.first; Real ref_stddev = length_stddev.second; Real min_length = ref_length - bond_tolerance*ref_stddev; Real max_length = ref_length + bond_tolerance*ref_stddev; Real zscore = (blength - ref_length)/ref_stddev; if (blength < min_length || blength > max_length) { - LOG_INFO("BOND:" << " " << res.GetChain() << " " << res.GetName() << " " << res.GetNumber() << " " << bond_str << " " << min_length << " " << max_length << " " << blength << " " << zscore << " " << "FAIL") + LOG_INFO("BOND:" << " " << res.GetChain() << " " + << res.GetName() << " " << res.GetNumber() << " " + << bond_str << " " << min_length << " " << max_length + << " " << blength << " " << zscore << " " << "FAIL") bad_bond_count++; - UniqueAtomIdentifier atom_ui(atom.GetResidue().GetChain().GetName(),atom.GetResidue().GetNumber(),atom.GetResidue().GetName(),atom.GetName()); - UniqueAtomIdentifier other_atom_ui(other_atom.GetResidue().GetChain().GetName(),other_atom.GetResidue().GetNumber(),other_atom.GetResidue().GetName(),other_atom.GetName()); - StereoChemicalBondViolation bond_v(atom_ui,other_atom_ui,blength,std::make_pair(min_length,max_length)); + UniqueAtomIdentifier atom_ui(atom.GetResidue().GetChain().GetName(), + atom.GetResidue().GetNumber(), + atom.GetResidue().GetName(), atom.GetName()); + UniqueAtomIdentifier other_atom_ui(other_atom.GetResidue().GetChain().GetName(), + other_atom.GetResidue().GetNumber(), + other_atom.GetResidue().GetName(), + other_atom.GetName()); + StereoChemicalBondViolation bond_v(atom_ui, other_atom_ui, blength, + std::make_pair(min_length, max_length)); bond_violation_list.push_back(bond_v); remove_sc=true; if (always_remove_bb==true) { @@ -432,7 +441,10 @@ std::pair<EntityView,StereoChemistryInfo> CheckStereoChemistry(const EntityView& } } } else { - LOG_VERBOSE("BOND:" << " " << res.GetChain() << " " << res.GetName() << " " << res.GetNumber() << " " << bond_str << " " << min_length << " " << max_length << " " << blength << " " << zscore << " " << "PASS") + LOG_VERBOSE("BOND:" << " " << res.GetChain() << " " + << res.GetName() << " " << res.GetNumber() << " " + << bond_str << " " << min_length << " " << max_length + << " " << blength << " " << zscore << " " << "PASS") } bond_count++; running_sum_zscore_bonds+=zscore; @@ -486,12 +498,25 @@ std::pair<EntityView,StereoChemistryInfo> CheckStereoChemistry(const EntityView& Real max_width = ref_width + angle_tolerance*ref_stddev; Real zscore = (awidth - ref_width)/ref_stddev; if (awidth < min_width || awidth > max_width) { - LOG_INFO("ANGLE:" << " " << res.GetChain() << " " << res.GetName() << " " << res.GetNumber() << " " << angle_str << " " << min_width << " " << max_width << " " << awidth << " " << zscore << " " << "FAIL") + LOG_INFO("ANGLE:" << " " << res.GetChain() << " " + << res.GetName() << " " << res.GetNumber() << " " + << angle_str << " " << min_width << " " << max_width + << " " << awidth << " " << zscore << " " << "FAIL") bad_angle_count++; - UniqueAtomIdentifier atom1_ui(atom1.GetResidue().GetChain().GetName(),atom1.GetResidue().GetNumber(),atom1.GetResidue().GetName(),atom1.GetName()); - UniqueAtomIdentifier atom_ui(atom.GetResidue().GetChain().GetName(),atom.GetResidue().GetNumber(),atom.GetResidue().GetName(),atom.GetName()); - UniqueAtomIdentifier atom2_ui(atom2.GetResidue().GetChain().GetName(),atom2.GetResidue().GetNumber(),atom2.GetResidue().GetName(),atom2.GetName()); - StereoChemicalAngleViolation angle_v(atom1_ui,atom_ui,atom2_ui,awidth,std::make_pair(min_width,max_width)); + UniqueAtomIdentifier atom1_ui(atom1.GetResidue().GetChain().GetName(), + atom1.GetResidue().GetNumber(), + atom1.GetResidue().GetName(), + atom1.GetName()); + UniqueAtomIdentifier atom_ui(atom.GetResidue().GetChain().GetName(), + atom.GetResidue().GetNumber(), + atom.GetResidue().GetName(), + atom.GetName()); + UniqueAtomIdentifier atom2_ui(atom2.GetResidue().GetChain().GetName(), + atom2.GetResidue().GetNumber(), + atom2.GetResidue().GetName(), + atom2.GetName()); + StereoChemicalAngleViolation angle_v(atom1_ui,atom_ui,atom2_ui,awidth, + std::make_pair(min_width,max_width)); angle_violation_list.push_back(angle_v); remove_sc=true; if (always_remove_bb==true) { @@ -502,7 +527,10 @@ std::pair<EntityView,StereoChemistryInfo> CheckStereoChemistry(const EntityView& remove_bb=true; } } else { - LOG_VERBOSE("ANGLE:" << " " << res.GetChain() << " " << res.GetName() << " " << res.GetNumber() << " " << angle_str << " " << min_width << " " << max_width << " " << awidth << " " << zscore << " " << "PASS") + LOG_VERBOSE("ANGLE:" << " " << res.GetChain() << " " + << res.GetName() << " " << res.GetNumber() << " " + << angle_str << " " << min_width << " " << max_width + << " " << awidth << " " << zscore << " " << "PASS") } angle_count++; running_sum_zscore_angles+=zscore; @@ -610,7 +638,14 @@ std::pair<EntityView,ClashingInfo> FilterClashes(const EntityView& ent, const Cl Real threshold=distance-tolerance; if (d<threshold*threshold) { - LOG_INFO(atom.GetResidue().GetChain() << " " << atom.GetResidue().GetName() << " " << atom.GetResidue().GetNumber() << " " << atom.GetName() << " " << atom2.GetResidue().GetChain() << " " << atom2.GetResidue().GetName() << " " << atom2.GetResidue().GetNumber() << " " << atom2.GetName() << " " << threshold << " " << sqrt(d) << " " << sqrt(d)-threshold << " " << "FAIL") + LOG_INFO("CLASH:" << " " << atom.GetResidue().GetChain() << " " + << atom.GetResidue().GetName() << " " + << atom.GetResidue().GetNumber() << " " << atom.GetName() + << " " << atom2.GetResidue().GetChain() << " " + << atom2.GetResidue().GetName() << " " + << atom2.GetResidue().GetNumber() << " " << atom2.GetName() + << " " << threshold << " " << sqrt(d) << " " + << sqrt(d)-threshold << " " << "FAIL") bad_distance_count++; UniqueAtomIdentifier atom_ui(atom.GetResidue().GetChain().GetName(),atom.GetResidue().GetNumber(),atom.GetResidue().GetName(),atom.GetName()); UniqueAtomIdentifier atom2_ui(atom2.GetResidue().GetChain().GetName(),atom2.GetResidue().GetNumber(),atom2.GetResidue().GetName(),atom2.GetName()); @@ -626,7 +661,14 @@ std::pair<EntityView,ClashingInfo> FilterClashes(const EntityView& ent, const Cl remove_bb=true; } } else { - LOG_VERBOSE("CLASH:" << " " << atom.GetResidue().GetChain() << " " << atom.GetResidue().GetName() << " " << atom.GetResidue().GetNumber() << " " << atom.GetName() << " " << atom2.GetResidue().GetChain() << " " << atom2.GetResidue().GetNumber() << " " << atom2.GetResidue().GetName() << " " << atom2.GetName() << " " << threshold << " " << sqrt(d) << " " << sqrt(d)-threshold << " " << "PASS") + LOG_VERBOSE("CLASH:" << " " << atom.GetResidue().GetChain() << " " + << atom.GetResidue().GetName() << " " + << atom.GetResidue().GetNumber() << " " << atom.GetName() + << " " << atom2.GetResidue().GetChain() << " " + << atom2.GetResidue().GetName() << " " + << atom2.GetResidue().GetNumber() << " " << atom2.GetName() + << " " << threshold << " " << sqrt(d) << " " + << sqrt(d)-threshold << " " << "PASS") } distance_count++; } diff --git a/modules/mol/alg/src/lddt.cc b/modules/mol/alg/src/lddt.cc index a64f5bb192d5abab6777b750786858d1fd477d16..dbe122d80d0b13e250edaec00291ea268bfa2bc6 100644 --- a/modules/mol/alg/src/lddt.cc +++ b/modules/mol/alg/src/lddt.cc @@ -311,11 +311,14 @@ int main (int argc, char **argv) settings.PrintParameters(); if (structural_checks) { LOG_INFO("Log entries format:"); + // verbosity level/format must match filter_clashes::CheckStereoChemistry LOG_INFO("BOND INFO FORMAT: Chain Residue ResNum Bond Min Max Observed Z-score Status"); LOG_INFO("ANGLE INFO FORMAT: Chain Residue ResNum Angle Min Max Observed Z-score Status"); - LOG_INFO("CLASH INFO FORMAT: Chain1 Residue1 ResNum1 Atom1 Chain2 Residue2 ResNum2 Atom2 Observed Difference Status"); + // verbosity level/format must match filter_clashes::FilterClashes + LOG_INFO("CLASH INFO FORMAT: Chain1 Residue1 ResNum1 Atom1 Chain2 Residue2 ResNum2 Atom2 Threshold Observed Difference Status"); } - LOG_INFO("LDDT INFO FORMAT: Chain1 Residue1 ResNum1 Atom1 Chain2 Residue2 ResNum2 Atom2 ModelDist TargetDist Difference Tolerance Status"); + // verbosity level/format must match local_dist_diff_test::calc_overlap1 + LOG_VERBOSE("LDDT INFO FORMAT: Chain1 Residue1 ResNum1 Atom1 Chain2 Residue2 ResNum2 Atom2 ModelDist TargetDistMin TargetDistMax Tolerance Status"); // error if the reference structure is empty if (glob_dist_list.size()==0) { diff --git a/modules/mol/alg/src/local_dist_diff_test.cc b/modules/mol/alg/src/local_dist_diff_test.cc index c52bf115666962310aef3de44b126eb5dd300e9b..1182dc375e8c44c01685cc55341ea35269c776fd 100644 --- a/modules/mol/alg/src/local_dist_diff_test.cc +++ b/modules/mol/alg/src/local_dist_diff_test.cc @@ -27,7 +27,8 @@ String vector_to_string(std::vector<Real> vec) { // helper function bool within_tolerance(Real mdl_dist, const std::pair<Real,Real>& values, Real tol) { - + // this considers min & max distance observed in all reference structures + // -> for single reference case values.first and values.second are the same return (values.first-tol)<=mdl_dist && (values.second+tol)>=mdl_dist; } @@ -82,10 +83,15 @@ std::pair<long int, long int> calc_overlap1(const ResidueRDMap& res_distance_lis Real tol = * tol_list_it; if (within_tolerance(mdl_dist,values,tol)) { if (log) { - LOG_VERBOSE("lddt:" << " " << av1.GetResidue().GetChain() << " " << av1.GetResidue().GetName() << " " << av1.GetResidue().GetNumber() << " " << av1.GetName() - << " " << av2.GetResidue().GetChain() << " " << av2.GetResidue().GetName() << " " << av2.GetResidue().GetNumber() << " " << av2.GetName() << " " - << mdl_dist << " " << values.first << " " << values.second << " " << tol << " " << "PASS") - } + LOG_VERBOSE("LDDT:" << " " << av1.GetResidue().GetChain() << " " + << av1.GetResidue().GetName() << " " + << av1.GetResidue().GetNumber() << " " << av1.GetName() + << " " << av2.GetResidue().GetChain() << " " + << av2.GetResidue().GetName() << " " + << av2.GetResidue().GetNumber() << " " << av2.GetName() + << " " << mdl_dist << " " << values.first << " " + << values.second << " " << tol << " " << "PASS") + } overlap.first+=1; if (!only_fixed) { overlap_list[rindex1].first+=1; @@ -93,10 +99,15 @@ std::pair<long int, long int> calc_overlap1(const ResidueRDMap& res_distance_lis } } else { if (log) { - LOG_VERBOSE("lddt:" << " " << av1.GetResidue().GetChain() << " " << av1.GetResidue().GetName() << " " << av1.GetResidue().GetNumber() << " " << av1.GetName() - << " " << av2.GetResidue().GetChain() << " " << av2.GetResidue().GetName() << " " << av2.GetResidue().GetNumber() << " " << av2.GetName() << " " - << mdl_dist << " " << values.first << " " << values.second << " " << tol << " " << "FAIL"); - } + LOG_VERBOSE("LDDT:" << " " << av1.GetResidue().GetChain() << " " + << av1.GetResidue().GetName() << " " + << av1.GetResidue().GetNumber() << " " << av1.GetName() + << " " << av2.GetResidue().GetChain() << " " + << av2.GetResidue().GetName() << " " + << av2.GetResidue().GetNumber() << " " << av2.GetName() + << " " << mdl_dist << " " << values.first << " " + << values.second << " " << tol << " " << "FAIL"); + } break; } }