#ifndef __JCOMPAREHISTOGRAMS__JTESTCHI2_SLICE2D__ #define __JCOMPAREHISTOGRAMS__JTESTCHI2_SLICE2D__ #include #include #include "JLang/JException.hh" #include "JCompareHistograms/JTest_t.hh" #include "JCompareHistograms/JResultTitle.hh" #include "TH1.h" #include "TH2.h" /** * \author rgruiz, bjung */ namespace JCOMPAREHISTOGRAMS {} namespace JPP { using namespace JCOMPAREHISTOGRAMS; } namespace JCOMPAREHISTOGRAMS { /** * Implementation of the Chi2 test for 2D histograms.\n * This class is derived from the abstract class JTest_t(). For a general description of the implementation of this and other tests derived from JTest_t(), see its documentation.\n * The parameter slice() can have the values x, X, y or Y. The histograms are sliced along the corresponding axis, and the JChi2Test() test is applied to each slice.\n * The input parameter `threshold`, is used to evaluate whether the test is passed or failed for each slice.\n * The evaluation is done by comparing the `threshold` value with the result produced by the Chi2 test for each slice. The output of JChi2Test() is a p-value.\n * The parameter `threshold` should therefore be a real value between 0 and 1. */ class JTestChi2_Slice2D: public JTest_t { public: /** * Default constructor. */ JTestChi2_Slice2D() : JTest_t("Chi2_Slice2D", "p-Value(chi2)") {} /** * Applies Chi2 test for two ROOT TH2 histograms. * * \param o1 First histogram * \param o2 Second histogram */ void test(const TObject* o1, const TObject* o2) override { using namespace std; using namespace JPP; const TH2* h1 = dynamic_cast(o1); const TH2* h2 = dynamic_cast(o2); if (h1 == NULL || h2 == NULL) { THROW(JValueOutOfRange, "JTestChi2_Slice2D::test(): Could not cast given TObjects to TH2"); } TH1* h3 = NULL; const char* const h3name = (h1->GetName() == h2->GetName() ? MAKE_CSTRING(h1->GetName() << "_" << testName << "_" << slice) : MAKE_CSTRING(h1->GetName() << "_VS_" << h2->GetName() << "_" << testName << "_" << slice)); int nFailures = 0; int nSlices = 0; if (slice == 'x' || slice == 'X') { const int nSlices1 = h1->GetNbinsX(); const int nSlices2 = h2->GetNbinsX(); if (nSlices1 != nSlices2) { THROW(JValueOutOfRange, "JTestChi2_Slice2D::test(): Histograms with different binning. The objects: " << h1->GetName() << " and " << h2->GetName() << " can not be compared." << endl); } nSlices = nSlices1; h3 = h1->ProjectionX(h3name); for (int i=1 ; i<=nSlices1 ; ++i){ const std::string sliceName = MAKE_STRING(h3->GetName() << "_" << to_string(i)); const TH1* s1 = h1->ProjectionY(sliceName.c_str(),i,i); const TH1* s2 = h2->ProjectionY(sliceName.c_str(),i,i); if (!(s1->GetSumOfWeights() > 0 && s2->GetSumOfWeights() > 0)) { continue; } const double chi2 = s1->Chi2Test(s2 , options.c_str()); const bool failed = (options.find("CHI2") != std::string::npos ? (chi2 > threshold ? true : false) : (chi2 < threshold ? true : false)); if (failed) nFailures++; h3->SetBinContent(i,chi2); } } else if (slice == 'y' || slice == 'Y') { const int nSlices1 = h1->GetNbinsY(); const int nSlices2 = h2->GetNbinsY(); if (nSlices1 != nSlices2) { THROW(JValueOutOfRange, "JTestChi2_Slice2D::test(): Histograms with different binning. The objects: " << h1->GetName() << " and " << h2->GetName() << " can not be compared." << endl); } nSlices = nSlices1; h3 = h1->ProjectionY(h3name); for (int i=1 ; i<=nSlices1 ; ++i){ const std::string sliceName = MAKE_STRING(h3->GetName() << "_" << to_string(i)); const TH1* s1 = h1->ProjectionX(sliceName.c_str(),i,i); const TH1* s2 = h2->ProjectionX(sliceName.c_str(),i,i); if (!(s1->GetSumOfWeights() > 0 && s2->GetSumOfWeights() > 0)) { continue; } const double chi2 = s1->Chi2Test (s2 , options.c_str()); const bool failed = (options.find("CHI2") != std::string::npos ? (chi2 > threshold ? true : false) : (chi2 < threshold ? true : false));; if (failed) nFailures++; h3->SetBinContent(i,chi2); } } const bool passed = (nFailures/nSlices < failuresThreshold); const JResultTitle title(testName, resultType, passed , nFailures); h3->SetTitle(title.getTitle().c_str()); h3->GetYaxis()->SetTitle(resultType.c_str()); const JTestResult r(testName, JRootObjectID(MAKE_STRING(h1->GetDirectory()->GetPath() << h1->GetName())), JRootObjectID(MAKE_STRING(h2->GetDirectory()->GetPath() << h1->GetName())), resultType, nFailures, failuresThreshold, h3, passed); this->push_back(r); } /** * Read test parameters from input. * * \param in input stream * \return input stream */ std::istream& read(std::istream& in) override { using namespace JPP; in >> threshold >> failuresThreshold >> slice >> options; if (threshold < 0.0) { THROW(JValueOutOfRange, "JTestChi2_Slice2D::read(): Invalid threshold value " << threshold); } if (failuresThreshold < 0.0 || failuresThreshold > 1.0) { THROW(JValueOutOfRange, "JTestChi2_Slice2D::read(): Invalid failuresThreshold value " << failuresThreshold); } if (slice != 'x' && slice != 'X' && slice != 'y' && slice != 'Y') { THROW(JValueOutOfRange, "JTestChi2_Slice2D::read(): Invalid slice option " << slice); } return in; } private: double threshold; //!< threshold p-value to decide if test is passed. double failuresThreshold; //!< threshold p-value to decide if test is passed. char slice; //!< axis to slice. std::string options; //!< options for the ROOT chi2 test. }; } #endif