FCBin.cxx
Go to the documentation of this file.
1 #include "CAFAna/FC/FCBin.h"
2 
3 #include "CAFAna/FC/FCPoint.h"
4 
5 #include <algorithm>
6 #include <cassert>
7 
8 #include "TH1.h"
9 #include "TTree.h"
10 
11 #include <iostream>
12 // Anonymous namespace, no-one else needs to see these
13 namespace{
14  int gNVals;
15  std::vector<float> gVals;
16 }
17 
18 namespace ana
19 {
20  //----------------------------------------------------------------------
21  void FCBin::Add(const FCPoint& pt)
22  {
23  fVals.push_back(pt.Up());
24  }
25 
26  //----------------------------------------------------------------------
27  double FCBin::UpValueForSignificance(double sig)
28  {
29  assert(sig > 0 && sig < 1);
30 
31  std::sort(fVals.begin(), fVals.end());
32 
33  const int N = fVals.size();
34 
35  const int idx = sig*(N+1)-1;
36  // Bin must have enough points for this significance
37  assert(idx >= 0);
38  assert(idx < N-1);
39 
40  // Up values for the points either side of this significance
41  const double up0 = fVals[idx];
42  const double up1 = fVals[idx+1];
43  // What significance taking the points themselves would be worth. Consider
44  // some small-N cases. If there is only one point it represents the 50%
45  // mark. If there are two, 50% is halfway between them, and they are 33%
46  // and 67%.
47  const double sig0 = double(idx+1)/(N+1);
48  const double sig1 = double(idx+2)/(N+1);
49  const double dsig = 1./(N+1); //sig1-sig0;
50 
51  // Make sure we picked the right points
52  assert(sig >= sig0 && sig <= sig1);
53 
54  // Linear interpolation
55  return (up0*(sig1-sig)+up1*(sig-sig0))/dsig;
56  }
57 
58  //----------------------------------------------------------------------
59  double FCBin::SignificanceForUpValue(double up) const
60  {
61  std::sort(fVals.begin(), fVals.end());
62  auto it = std::lower_bound(fVals.begin(), fVals.end(), up);
63 
64  auto v1 = *it;
65  auto v0 = *(it-1);
66 
67  // up value must be covered in range
68  assert(up >= v0 && up <= v1);
69  // should never happen with lower_bound
70  if (v0==v1) return double(it-fVals.begin()-1);
71 
72  return double((up-v0)/(v1-v0) + it-fVals.begin()-1) / fVals.size();
73  }
74 
75  //----------------------------------------------------------------------
76  void FCBin::Draw(double maxx) const
77  {
78  TH1* h = new TH1F("", ";#Delta#chi^{2};Experiments", 200, 0, maxx);
79  for(float v: fVals) h->Fill(v);
80  h->Draw();
81 
82  /*
83  double* integ = h->GetIntegral();
84  for(int n = 0; n < 102; ++n) h->SetBinContent(n, integ[n]);
85  */
86  }
87 
88  //----------------------------------------------------------------------
89  void FCBin::InitToTree(TTree* tr)
90  {
91  tr->Branch("nvals", &gNVals);
92  tr->Branch("vals", &gVals.front(), "vals[nvals]/F");
93  }
94 
95  //----------------------------------------------------------------------
96  void FCBin::ToTree(TTree* tr) const
97  {
98  gNVals = fVals.size();
99 
100  gVals = fVals;
101  tr->SetBranchAddress("vals", &gVals.front());
102 
103  tr->Fill();
104  }
105 
106  //----------------------------------------------------------------------
107  void FCBin::InitFromTree(TTree* tr)
108  {
109  tr->SetBranchAddress("nvals", &gNVals);
110 
111  // It's a shape we have to run through the tree upfront to determine this,
112  // but TTree doesn't seem to like having the branch address changed
113  // mid-entry as you'd have to do otherwise.
114  int maxN = 0;
115  for(int i = 0; i < tr->GetEntries(); ++i){
116  tr->GetEntry(i);
117  maxN = std::max(maxN, gNVals);
118  }
119 
120  gVals.resize(maxN);
121  tr->SetBranchAddress("vals", &gVals.front());
122  }
123 
124  //----------------------------------------------------------------------
126  {
127  tr->GetEntry(idx);
128 
129  assert(gNVals <= int(gVals.size()));
130 
131  FCBin ret;
132  ret.fVals = std::vector<float>(gVals.begin(), gVals.begin()+gNVals);
133  return ret;
134  }
135 } // namespace
T max(const caf::Proxy< T > &a, T b)
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
set< int >::iterator it
void Draw(double maxx=10) const
Definition: FCBin.cxx:76
double UpValueForSignificance(double sig)
Definition: FCBin.cxx:27
double SignificanceForUpValue(double up) const
Definition: FCBin.cxx:59
Represents the results of a single Feldman-Cousins pseudo-experiment.
Definition: FCPoint.h:8
void ToTree(TTree *tr) const
Definition: FCBin.cxx:96
A collection of Feldman-Cousins experiments at the same oscillation parameters.
Definition: FCBin.h:13
static FCBin FromTree(TTree *tr, int idx)
Definition: FCBin.cxx:125
static void InitFromTree(TTree *tr)
Definition: FCBin.cxx:107
std::vector< float > fVals
Definition: FCBin.h:31
static void InitToTree(TTree *tr)
Definition: FCBin.cxx:89
assert(nhit_max >=nhit_nbins)
double Up() const
Definition: FCPoint.h:38
void Add(const FCPoint &pt)
Definition: FCBin.cxx:21