JTestSignificance_t.hh

Go to the documentation of this file.

#ifndef __JTESTSIGNIFICANCE_T__
#define __JTESTSIGNIFICANCE_T__

#include <istream>
#include <ostream>
#include "JCompareHistograms/JTest_t.hh"
#include "TMath.h"


/**
 *
 * \author bofearraigh
 * \author rgruiz
 *
 * Compares two histograms \f$ H_a \f$ and \f$ H_b \f$ by calculating the average bin significance \n
 *
 * \f[
 * S = \frac{1}{N} \cdot \sum_{i=0}^{N} \frac{|a_i - K \cdot b_i|}{\sqrt{\sigma ^2(a_i) + K^2 \cdot \sigma ^2(b_i)}}
 * \f]
 *
 * where:
 *  - \f$ N \f$ is the number of bins in \f$ H_a \f$ and \f$ H_b \f$.
 *  - \f$ a_i \f$ and \f$ b_i \f$ are the contents of bin \f$ i \f$ in \f$ H_a \f$ and \f$ H_b \f$.
 *  - \f$ K \f$ is a normalisation parameter. It should be interpreted as the integral of \f$ H_a \f$ divided by the integral of \f$ H_b \f$
 *  - \f$ \sigma \f$  is calculated assuming that the bin contents follow Poisson distributions.
 *
 * A threshold value for \f$ S \f$ is used to tag the test as passed or failed.
 */
class JTestSignificance_t 
{ 
public:  

  /**
   * Default constructor.
   */  
  JTestSignificance_t(){}

  /**
   * Significance test for histograms.\n
   *
   * \param h1 First object
   * \param h2 Second object
   * \param threshold 
   * \param K Normalisation parameter   
   * \param parameterName Name of the parameter used to test the histograms
   * \param testName Name of the test used to compare the histograms
   *
   * \return Test result.
   */
  template<class T>
  JTestResult JSignificanceTest(T* h1, T* h2, double threshold, double K, std::string testName, std::string parameterName) {    

    using namespace std;
    using namespace JPP;
    
    if(h1 -> GetNbinsX() != h2 -> GetNbinsX())
      ERROR("Histograms with different bining. The objects: " << h1 -> GetName() << " can not be compared." << endl);
        

    T* h3 = (T*)h1->Clone(h1->GetName()==h2->GetName() ?
                                MAKE_CSTRING(to_string(h1->GetName())) :
                                MAKE_CSTRING(to_string(h1->GetName()) + "_VS_" + to_string(h2->GetName())));
    
    h3->Add(h2,-1*K);
    
    double S = 0;

    for (int i=0 ; i < h1->GetNbinsX() ; ++i){
      for (int j=0 ; j< h1->GetNbinsY() ; ++j){
        for (int k=0 ; k< h1->GetNbinsY() ; ++k){
          
          double a = h1->GetBinContent(i+1,j+1,k+1);
          double b = h2->GetBinContent(i+1,j+1,k+1);

          if(a!=0 || b!=0){  
            S += fabs((a - K*b)/sqrt(a + K*K*b));
          }
        }       
      }
    }

    S /= (h1->GetNbinsX()*h1->GetNbinsY()*h1->GetNbinsZ());

    bool passed;
    
    (S < threshold ? passed = false : passed = true);
    
    JResultTitle title(testName, parameterName, passed, S);
    
    h3->SetTitle(title.getTitle().c_str());
  
    JTestResult r (testName,
                   string (h1->GetDirectory()->GetPath()).append(h1->GetName()),
                   string (h2->GetDirectory()->GetPath()).append(h2->GetName()),
                   h1->GetDirectory()->GetFile()->GetName(),
                   h2->GetDirectory()->GetFile()->GetName(),
                   parameterName, S , threshold, h3, passed);

    return r;
  };

}; 
  
#endif