From 2cfc42c14539fe80ae30ca452e481fd662335953 Mon Sep 17 00:00:00 2001
From: Valerio Mariani <valerio.mariani@unibas.ch>
Date: Fri, 23 Mar 2012 11:33:18 +0100
Subject: [PATCH] Added coverage resport

---
 modules/mol/alg/src/ldt.cc             | 22 ++++++++++++++++++++++
 modules/mol/alg/src/local_dist_test.cc |  4 ++--
 2 files changed, 24 insertions(+), 2 deletions(-)

diff --git a/modules/mol/alg/src/ldt.cc b/modules/mol/alg/src/ldt.cc
index 3b37839c6..d11ce8f22 100644
--- a/modules/mol/alg/src/ldt.cc
+++ b/modules/mol/alg/src/ldt.cc
@@ -27,6 +27,7 @@
 #include <ost/io/mol/pdb_reader.hh>
 #include <ost/io/io_exception.hh>
 #include <ost/conop/conop.hh>
+#include <ost/conop/amino_acids.hh>
 #include <ost/mol/iterator.hh>
 #include <ost/platform.hh>
 #include <ost/log.hh>
@@ -75,6 +76,25 @@ void usage()
   std::cerr << "   -e         print version" << std::endl;
 }
 
+std::pair<int,int> compute_coverage (const EntityView& v,const GlobalDistanceList& glob_dist_list)
+{
+  int second=0;
+  int first=0;  
+  ChainView vchain=v.GetChainList()[0]; 
+  for (std::vector<ResidueDistanceList>::const_iterator i=glob_dist_list.begin();i!=glob_dist_list.end();++i)
+  {
+    ResNum rnum = (*i)[0].GetFirstAtom().GetResNum();
+    String rname = (*i)[0].GetFirstAtom().GetResidueName();
+    if (ost::conop::ResidueNameToOneLetterCode(rname)!='X') {
+      second++;
+      if (vchain.FindResidue(rnum)) {
+        first++;     
+      }
+    }    
+  }
+  return std::make_pair<int,int>(first,second);  
+}
+
 int main (int argc, char **argv)
 {
   String version = "Beta - 2012-01-17"; 
@@ -236,9 +256,11 @@ int main (int argc, char **argv)
     }
     
     GlobalDistanceList glob_dist_list = CreateDistanceList(ref_view,radius);
+    std::pair<int,int> cov = compute_coverage(v,glob_dist_list);
     Real ldt=LDTHA(v, glob_dist_list);
     
     std::cout << "File: " << files[i] << std::endl; 
+    std::cout << "Coverage: " << (float(cov.first)/float(cov.second)) << " (" << cov.first << " out of " << cov.second << " residues)" << std::endl;
     std::cout << "Global LDT score: " << ldt << std::endl;
     std::cout << "Local LDT Score:" << std::endl;
     std::cout << "Chain\tResName\tResNum\tScore" << std::endl;
diff --git a/modules/mol/alg/src/local_dist_test.cc b/modules/mol/alg/src/local_dist_test.cc
index 13f37e2e0..d65205128 100644
--- a/modules/mol/alg/src/local_dist_test.cc
+++ b/modules/mol/alg/src/local_dist_test.cc
@@ -1,5 +1,6 @@
 #include <ost/log.hh>
 #include <ost/mol/mol.hh>
+#include <ost/conop/amino_acids.hh>
 #include "local_dist_test.hh"
 #include <boost/concept_check.hpp>
 
@@ -259,8 +260,7 @@ GlobalDistanceList CreateDistanceList(const EntityView& ref,Real max_dist)
  ResidueViewList ref_residues=ref.GetChainList()[0].GetResidueList();
  for (ResidueViewList::iterator i=ref_residues.begin(), e=ref_residues.end(); i!=e; ++i) {
    ResidueView rview = (*i);
-//#   if (rview.IsPeptideLinking()) {
-   if (true) {
+   if (ost::conop::ResidueNameToOneLetterCode(rview.GetName())!='X') {
      ResidueDistanceList res_dist_list;
      AtomViewList ref_atoms=(*i).GetAtomList();
      AtomViewList within;
-- 
GitLab