monoprob.C
Go to the documentation of this file.
7 #include "CAFAna/Vars/FitVars.h"
11 #include "OscLib/IOscCalc.h"
12 #include "Utilities/rootlogon.C"
13 
14 #include "TString.h"
15 #include "TFeldmanCousins.h"
16 #include "TGraph.h"
17 #include "TH2.h"
18 #include "TLatex.h"
19 #include "TLegend.h"
20 #include "TPad.h"
21 
23 using namespace ana;
24 
25 void monoprob(bool data = false)
26 {
28  double pot = kAna2017POT;
29  double dmsq32IH = -2.48e-3;//71e-3;
30  double dmsq32NH = 2.44e-3;//67e-3;
31  double ssth23Lo = 0.429;
32  double ssth23A = 0.468;
33  double ssth23B = 0.560;
34  double ssth23 = 0.5;
35 
36  double ssth23Hi = 0.597;
37  double ss2th13 = util::sqr(sin(2*calc->GetTh13()));
38 
39  std::string decomp = "combo";
40  auto pred = GetNuePrediction2017(decomp, DefaultOscCalc(),{});
41  auto cosm = GetNueCosmics2017();
42  double nCosmics = cosm.first->Integral();
43 
44  auto c1 = new TCanvas(ana::UniqueName().c_str());
45  c1->SetFillStyle(0);
46  c1->SetLeftMargin(0.14);
47  c1->SetBottomMargin(0.15);
48 
49  const double ymax = 95;
50  TH2* axes = new TH2F("", data? ";;Total events":
51  ";;Total events expected",
52  100, 0, 2, 100, 0, ymax);
53  CenterTitles(axes);
54  axes->Draw();
55  axes->GetXaxis()->SetTitleSize(kBlessedTitleFontSize);
56  axes->GetYaxis()->SetTitleSize(kBlessedTitleFontSize);
57  axes->GetXaxis()->SetLabelSize(kBlessedLabelFontSize);
58  axes->GetYaxis()->SetLabelSize(kBlessedLabelFontSize);
59  axes->GetXaxis()->SetTitleOffset(0.8);
60  axes->GetYaxis()->SetTitleOffset(0.75);
61  XAxisDeltaCPLabels(axes);
62  c1->RedrawAxis();
63 
64 
65  TGraph* g[2][2] = {{new TGraph, new TGraph}, {new TGraph, new TGraph}};
66  TGraph* gerr[2] = {new TGraph, new TGraph};
67 
68  for(int hie: {+1, -1})
69  {
70  if(hie < 0)
71  calc->SetDmsq32(dmsq32IH);
72  else
73  calc->SetDmsq32(dmsq32NH);
74 
75  const int hieIdx = (hie+1)/2;
76 
77  const int col = (hie > 0) ? kPrimColorNH : kPrimColorIH;
78  const int fillC = (hie > 0) ? kSecoColorNH : kSecoColorIH;
79 
80  for (int i =0; i< 2; ++i){
81  g[hieIdx][i]->SetLineWidth(3);
82  g[hieIdx][i]->SetLineColor(col);
83  gerr[hieIdx]->SetLineWidth(3);
84  gerr[hieIdx]->SetFillColorAlpha(fillC, .33);
85 
86  kFitSinSqTheta23.SetValue(calc, i==0?ssth23A:ssth23B);
87 // if(i==1) g[hieIdx][i]->SetLineStyle(7);
88 
89  for(double delta = 0; delta < 2.025; delta += .05){
90  calc->SetdCP(delta*M_PI);
91 
92  g[hieIdx][i]->SetPoint(g[hieIdx][i]->GetN(), delta,
93  pred->Predict(calc).Integral(pot) + nCosmics);
94  } // delta
95  }
96 
97  kFitSinSqTheta23.SetValue(calc, ssth23Lo);
98 
99  for(double delta = 0; delta < 2.025; delta += .05){
100  calc->SetdCP(delta*M_PI);
101 
102  gerr[hieIdx]->SetPoint(gerr[hieIdx]->GetN(), delta,
103  pred->Predict(calc).Integral(pot) + nCosmics);
104  } // delta
105 
106  kFitSinSqTheta23.SetValue(calc, ssth23Hi);
107 
108  for(double delta = 2.025; delta > -.1; delta -= .05){
109  calc->SetdCP(delta*M_PI);
110 
111  gerr[hieIdx]->SetPoint(gerr[hieIdx]->GetN(), delta,
112  pred->Predict(calc).Integral(pot) + nCosmics);
113  } // delta
114  } // hie
115 
116  for(int h = 0; h< 2; ++h){
117  gerr[h]->Draw("f same");
118  for (int i = 0; i<2; ++i){
119  g[h][i]->Draw("l same");
120  }
121  }
122 
123  TGraph* daterr = 0;
124  if(data){
125  TFeldmanCousins fc(0.6827); // 1 sigma
126  int ndata = 66;
127  fc.SetMuMax (100);
128 
129  daterr = new TGraph;
130 
131  daterr->SetPoint(0, -.1, fc.CalculateLowerLimit(ndata, 0));
132  daterr->SetPoint(1, 2.1, fc.CalculateLowerLimit(ndata, 0));
133  daterr->SetPoint(2, 2.1, fc.CalculateUpperLimit(ndata, 0));
134  daterr->SetPoint(3, -.1, fc.CalculateUpperLimit(ndata, 0));
135 
136  daterr->SetFillColorAlpha(kCentralValueColor, .5);//.33);
137  std::cout << "lower" << fc.CalculateLowerLimit(ndata, 0) <<"\n";
138  std::cout << "upper" << fc.CalculateUpperLimit(ndata, 0) <<"\n";
139 
140  std::cout << "lower" << fc.CalculateLowerLimit(33, 0) <<"\n";
141  std::cout << "upper" << fc.CalculateUpperLimit(33, 0) <<"\n";
142  daterr->Draw("f same");
143 
144  TGraph* gdat = new TGraph;
145  gdat->SetLineWidth(3);
146  gdat->SetLineColor(kCentralValueColor);
147  gdat->SetPoint(0, 0, ndata);
148  gdat->SetPoint(1, 2, ndata);
149  gdat->Draw("l same");
150  }
151 
152  TLatex* ltx = new TLatex(.16, .8, TString::Format(
153  "#splitline{NOvA FD}{%g#times10^{20} POT eq.}",
154  kAna2017FullDetEquivPOT/1.e20));
155 
156  ltx->SetNDC();
157  ltx->SetTextAlign(11);
158  ltx->SetTextSize(kBlessedLabelFontSize);
159  ltx->Draw();
160 
161  ltx = new TLatex(.875, .85, TString::Format("sin^{2}2#theta_{13}=%.3f",ss2th13));
162  ltx->SetNDC();
163  ltx->SetTextAlign(32);
164  ltx->SetTextSize(kBlessedLabelFontSize);
165  ltx->Draw();
166 
167  ltx = new TLatex(.875, .79, TString::Format(
168  "sin^{2}#theta_{23}=%.2f-%.2f",
169  ssth23Lo, ssth23Hi).Data());
170  ltx->SetNDC();
171  ltx->SetTextAlign(32);
172  ltx->SetTextSize(kBlessedLabelFontSize);
173  ltx->Draw();
174 
175  // Make legend
176  TLegend* leg = new TLegend(.35, .15, .85, .4);
177  leg->SetFillStyle(0);
178  leg->SetTextSize(kBlessedLabelFontSize);
179 
180  TGraph* dummy = new TGraph;
181  dummy->SetLineWidth(3);
182 
183  if(data){
184  dummy->SetLineColor(kCentralValueColor);
185  dummy->SetFillColorAlpha(kCentralValueColor, .5);
186  leg->AddEntry(dummy->Clone(), "Data (#pm1#sigma)", "lf");
187  }
188  else leg->AddEntry((TObject*)0, "", "");
189  dummy->SetLineColor(kPrimColorNH);
190  dummy->SetFillColorAlpha(kSecoColorNH, .33);
191  leg->AddEntry(dummy->Clone(), TString::Format(
192  "NH: #Deltam_{32}^{2}=#plus%g#times10^{-3}eV^{2}", dmsq32NH/1e-3).Data(), "lf");
193  dummy->SetLineColor(kPrimColorIH);
194  dummy->SetFillColorAlpha(kSecoColorIH, .33);
195  leg->AddEntry(dummy->Clone(), TString::Format(
196  "IH: #Deltam_{32}^{2}=#minus%g #times10^{-3}eV^{2}", fabs(dmsq32IH)/1e-3).Data(), "lf");
197 
198 
199  leg->Draw();
200 
201 
202  if(data) Preliminary(); else Simulation();
203 
204  gPad->SetFillStyle(0);
205  gPad->RedrawAxis();
206 
207  TBox* reax = new TBox(0, 0, 2, ymax);
208  reax->SetFillStyle(0);
209  reax->SetLineWidth(2);
210  reax->Draw("l");
211 
212  for(TString ext: {".png", ".pdf"}){
213  gPad->Print((TString)"monoprob_2017" + (data? "_data":"") + ext);
214  }
215 }
const Color_t kPrimColorIH
Definition: Style.h:64
void Simulation()
Definition: tools.h:16
const Color_t kPrimColorNH
Definition: Style.h:61
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
double ssth23
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
double delta
Definition: runWimpSim.h:98
void XAxisDeltaCPLabels(TH1 *axes, bool t2kunits)
Label the x-axis with fractions of pi.
void monoprob(bool data=false)
Definition: monoprob.C:25
const IPrediction * GetNuePrediction2017(std::string decomp, osc::IOscCalc *calc, bool corrSysts, bool GetFromUPS=false)
const double kAna2017POT
Definition: Exposures.h:174
void CenterTitles(TH1 *histo)
Definition: Plots.cxx:1481
osc::OscCalcDumb calc
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
Definition: Calcs.cxx:49
T sqr(T x)
More efficient square function than pow(x,2)
Definition: MathUtil.h:23
virtual void SetDmsq32(const T &dmsq32)=0
#define M_PI
Definition: SbMath.h:34
const Color_t kSecoColorIH
Definition: Style.h:65
const XML_Char const XML_Char * data
Definition: expat.h:268
Double_t ymax
Definition: plot.C:25
std::pair< TH1D *, double > GetNueCosmics2017(bool GetFromUPS=false)
const Color_t kCentralValueColor
Definition: Style.h:86
void SetValue(osc::IOscCalcAdjustable *osc, double val) const override
Int_t col[ntarg]
Definition: Style.C:29
#define pot
void Preliminary()
OStream cout
Definition: OStream.cxx:6
virtual T GetTh13() const
Definition: IOscCalc.h:93
T sin(T number)
Definition: d0nt_math.hpp:132
const double kAna2017FullDetEquivPOT
Definition: Exposures.h:195
double ss2th13
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
const Color_t kSecoColorNH
Definition: Style.h:62
c1
Definition: demo5.py:24
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
const Float_t kBlessedLabelFontSize
Definition: Style.h:90
Float_t e
Definition: plot.C:35
const Float_t kBlessedTitleFontSize
Definition: Style.h:89
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:29
virtual void SetdCP(const T &dCP)=0
enum BeamMode string