Jpp  18.4.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JTestSignificance.hh
Go to the documentation of this file.
1 #ifndef __JCOMPAREHISTOGRAMS__JTESTSIGNIFICANCE__
2 #define __JCOMPAREHISTOGRAMS__JTESTSIGNIFICANCE__
3 
4 #include <istream>
5 #include <ostream>
6 
7 #include "JLang/JException.hh"
8 
10 
11 #include "TH1.h"
12 
13 
14 /**
15  * \author rgruiz, bofearreigh, bjung
16  */
17 namespace JCOMPAREHISTOGRAMS {}
18 namespace JPP { using namespace JCOMPAREHISTOGRAMS; }
19 
20 namespace JCOMPAREHISTOGRAMS {
21 
22  /**
23  * Significance test
24  *
25  * Compares two histograms \f$ H_a \f$ and \f$ H_b \f$ by calculating the average bin significance \n
26  *
27  * \f[
28  * 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)}}
29  * \f]
30  *
31  * where:
32  * - \f$ N \f$ is the number of bins in \f$ H_a \f$ and \f$ H_b \f$.
33  * - \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$.
34  * - \f$ K \f$ is a normalisation parameter which can be specified by the user.
35  * - It should be interpreted as the integral of \f$ H_a \f$ divided by the integral of \f$ H_b \f$.
36  * - Note that if the specified value of K is negative, the value will be reset to the integral of \f$ H_a \f$ divided by the integral of \f$ H_b \f$.
37  * - \f$ \sigma \f$ is calculated assuming that the bin contents follow Poisson distributions.
38  *
39  * A threshold value for \f$ S \f$ is used to tag the test as passed or failed.
40  */
42  public JTest_t
43  {
44  public:
45 
46  /**
47  * Default constructor.
48  */
50  JTest_t("Significance", "Significance")
51  {}
52 
53 
54  /**
55  * Applies Significance test for two ROOT TH1 histograms.
56  *
57  * \param o1 First histogram
58  * \param o2 Second histogram
59  */
60  void test(const TObject* o1, const TObject* o2) override
61  {
62  using namespace std;
63  using namespace JPP;
64 
65  const TH1* h1 = dynamic_cast<const TH1*>(o1);
66  const TH1* h2 = dynamic_cast<const TH1*>(o2);
67 
68  if (h1 == NULL || h2 == NULL) {
69  THROW(JValueOutOfRange, "JTestSignificance::test(): Could not cast given TObjects to TH1.");
70  }
71 
72  if(h1->GetNbinsX() != h2->GetNbinsX() ||
73  h1->GetNbinsY() != h2->GetNbinsY() ||
74  h1->GetNbinsZ() != h2->GetNbinsZ()) {
75  THROW(JValueOutOfRange, "JTestSignificance::test(): Histograms with different bining. The objects: " <<
76  h1->GetName() << " and " << h2->GetName() << " can not be compared." << endl);
77  }
78 
79  if (K < 0) {
80  K = h2->GetEntries() / h1->GetEntries();
81  }
82 
83  TH1* h3 = (TH1*) h1->Clone(h1->GetName() == h2->GetName() ?
84  MAKE_CSTRING(h1->GetName() << "_" << testName) :
85  MAKE_CSTRING(h1->GetName() << "_VS_" << h2->GetName() << "_" << testName));
86 
87  h3->Add(h2,-1*K);
88 
89  double S = 0;
90 
91  for (int i=1 ; i <= h1->GetNbinsX() ; ++i){
92  for (int j=1 ; j <= h1->GetNbinsY() ; ++j){
93  for (int k=1 ; k <= h1->GetNbinsZ() ; ++k){
94 
95  const double a = h1->GetBinContent(i,j,k);
96  const double b = h2->GetBinContent(i,j,k);
97 
98  if (a!=0 || b!=0) {
99  S += fabs((a - K*b)/sqrt(a + K*K*b));
100  }
101  }
102  }
103  }
104 
105  S /= (h1->GetNbinsX()*h1->GetNbinsY()*h1->GetNbinsZ());
106 
107  const bool passed = (S > threshold);
108 
109  const JResultTitle title(testName, resultType, passed, S);
110 
111  h3->SetTitle(title.getTitle().c_str());
112 
113  const int Ndims = h3->GetDimension();
114 
115  if (Ndims == 1) {
116  h3->GetYaxis()->SetTitle(resultType.c_str());
117  } else if (Ndims == 2) {
118  h3->GetZaxis()->SetTitle(resultType.c_str());
119  }
120 
121  const JTestResult r(testName,
122  JRootObjectID(MAKE_STRING(h1->GetDirectory()->GetPath() << h1->GetName())),
123  JRootObjectID(MAKE_STRING(h2->GetDirectory()->GetPath() << h1->GetName())),
124  resultType, S, threshold, h3, passed);
125 
126  this->push_back(r);
127  }
128 
129 
130  /**
131  * Read test parameters from input.
132  *
133  * \param in input stream
134  * \return input stream
135  */
136  std::istream& read(std::istream& in) override
137  {
138  using namespace JPP;
139 
140  in >> threshold >> K;
141 
142  if (threshold < 0.0) {
143  THROW(JValueOutOfRange, "JTestSignificance::read(): Invalid threshold value " << threshold);
144  }
145 
146  return in;
147  }
148 
149  private:
150 
151  double threshold; //!< threshold p-value to decide if test is passed.
152  double K; //!< normalization factor between histograms.
153  };
154 }
155 
156 #endif
Exceptions.
then fatal No hydrophone data file $HYDROPHONE_TXT fi sort gr k
Interface to read input and write output for TObject tests.
Definition: JTest_t.hh:40
std::string getTitle() const
Returns a standard string to be used as title of a graphical root object.
Definition: JResultTitle.hh:57
Class dedicated to standardize the title of the graphical objects produced by the JTest_t() derived c...
Definition: JResultTitle.hh:25
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:712
Definition: JRoot.hh:19
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:136
Auxiliary class to handle file name, ROOT directory and object name.
void test(const TObject *o1, const TObject *o2) override
Applies Significance test for two ROOT TH1 histograms.
data_type r[M+1]
Definition: JPolint.hh:868
std::istream & read(std::istream &in) override
Read test parameters from input.
double K
normalization factor between histograms.
#define MAKE_STRING(A)
Make string.
Definition: JPrint.hh:127
then JCalibrateToT a
Definition: JTuneHV.sh:113
const std::string resultType
test result type
Definition: JTest_t.hh:181
const std::string testName
test name
Definition: JTest_t.hh:180
double threshold
threshold p-value to decide if test is passed.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:48
Structure containing the result of the test.
Definition: JTestResult.hh:28
int j
Definition: JPolint.hh:792