Jpp test-rotations-old
the software that should make you happy
Loading...
Searching...
No Matches
JTestPoissonLogLikelihoodRatio.hh
Go to the documentation of this file.
1#ifndef __JCOMPAREHISTOGRAMS__JTESTPOISSONLOGLIKELIHOODRATIO__
2#define __JCOMPAREHISTOGRAMS__JTESTPOISSONLOGLIKELIHOODRATIO__
3
4#include <istream>
5#include <ostream>
6
7#include "JLang/JException.hh"
8
11
12#include "TH1.h"
13
14
15/**
16 * \author bjung
17 */
18namespace JCOMPAREHISTOGRAMS {
19
20 /**
21 * Implementation of the Poisson log-likelihood ratio test.\n
22 * The first histogram is treated as an Asimov dataset corresponding to a given null hypothesis,\n
23 * which is compared to an alternative hypothesis given by the second histogram.\n\n
24 *
25 * 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
26 */
28 public JTest_t
29 {
30 public:
31
32 /**
33 * Default constructor.
34 */
36 JTest_t("Poissonian Log-Likelihood Ratio", "NLLR")
37 {}
38
39
40 /**
41 * Applies a Poissonian log-likelihood ratio test to the two given histograms.\n
42 * The first histogram is treated as an Asimov dataset corresponding to a given null hypothesis.\n
43 * The second histogram is treated as the expectation for the alternative hypothesis, to which the null hypothesis is compared.
44 *
45 * \param o1 First histogram
46 * \param o2 Second histogram
47 */
48 void test(const TObject* o1, const TObject* o2) override
49 {
50 using namespace std;
51 using namespace JPP;
52
53 const TH1* h1 = dynamic_cast<const TH1*>(o1);
54 const TH1* h2 = dynamic_cast<const TH1*>(o2);
55
56 if (h1 == NULL || h2 == NULL) {
57 THROW(JValueOutOfRange, "JTestPoissonLogLikelihood::test(): Could not cast given TObjects to TH1.");
58 }
59
60 if(h1->GetNbinsX() != h2->GetNbinsX() ||
61 h1->GetNbinsY() != h2->GetNbinsY() ||
62 h1->GetNbinsZ() != h2->GetNbinsZ()) {
63 THROW(JValueOutOfRange, "JTestPoissonLogLikelihood::test(): Histograms with different binning. The objects: " <<
64 h1->GetName() << " and " << h2->GetName() << " can not be compared." << endl);
65 }
66
67 TH1* h3 = (TH1*) h1->Clone(h1->GetName() == h2->GetName() ?
68 MAKE_CSTRING(h1->GetName() << "_" << testName) :
69 MAKE_CSTRING(h1->GetName() << "_VS_" << h2->GetName() << "_" << testName));
70
71 h3->Reset();
72
73 double LnLratio = 0.0;
74
75 for (int i=1 ; i <= h1->GetNbinsX() ; ++i){
76 for (int j=1 ; j <= h1->GetNbinsY() ; ++j){
77 for (int k=1 ; k <= h1->GetNbinsZ() ; ++k){
78
79 double contribution = 0.0;
80
81 const double n0 = h1->GetBinContent(i,j,k);
82 const double n1 = h2->GetBinContent(i,j,k);
83
84 if (n0 > 0.0 && n1 > 0.0) {
85 contribution = 2 * (n1 - n0 + n0 * log(n0/n1));
86 } else if (n1 > 0.0) {
87 contribution = 2 * n1;
88 }
89
90 h3->SetBinContent(i,j,k, contribution);
91
92 LnLratio += contribution;
93 }
94 }
95 }
96
97 const bool passed = (LnLratio < threshold);
98
99 const JResultTitle title(testName, resultType, passed, LnLratio);
100
101 h3->SetTitle(title.getTitle().c_str());
102 h3->GetYaxis()->SetTitle(resultType.c_str());
103
104 const JTestResult r (testName,
105 JRootObjectID(MAKE_STRING(h1->GetDirectory()->GetPath() << h1->GetName())),
106 JRootObjectID(MAKE_STRING(h2->GetDirectory()->GetPath() << h1->GetName())),
107 resultType, LnLratio, threshold, h3, passed);
108
109 this->push_back(r);
110 }
111
112
113 /**
114 * Read test parameters from input.
115 *
116 * \param in input stream
117 * \return input stream
118 */
119 std::istream& read(std::istream& in) override
120 {
121 using namespace JPP;
122
123 in >> threshold;
124
125 if (!(threshold > 0.0)) {
126 THROW(JValueOutOfRange, "JTestPoissonLogLikelihoodRatio::read(): Given chi-square threshold " << threshold << " is invalid");
127 }
128
129 return in;
130 }
131
132 private:
133
134 double threshold; //!< threshold chi-square to decide if test is passed.
135 };
136}
137
138#endif
Exceptions.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
#define MAKE_CSTRING(A)
Make C-string.
Definition JPrint.hh:72
#define MAKE_STRING(A)
Make string.
Definition JPrint.hh:63
Class dedicated to standardize the title of the graphical objects produced by the JTest_t() derived c...
std::string getTitle() const
Returns a standard string to be used as title of a graphical root object.
Implementation of the Poisson log-likelihood ratio test.
double threshold
threshold chi-square to decide if test is passed.
std::istream & read(std::istream &in) override
Read test parameters from input.
void test(const TObject *o1, const TObject *o2) override
Applies a Poissonian log-likelihood ratio test to the two given histograms.
Interface to read input and write output for TObject tests.
Definition JTest_t.hh:42
const std::string resultType
test result type
Definition JTest_t.hh:181
const std::string testName
test name
Definition JTest_t.hh:180
Auxiliary class to handle file name, ROOT directory and object name.
Exception for accessing a value in a collection that is outside of its range.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Structure containing the result of the test.