genie_contours.C
Go to the documentation of this file.
1 #include "CAFAna/Core/Registry.h"
3 #include "CAFAna/Fit/Fit.h"
6 #include "CAFAna/Core/Spectrum.h"
7 #include "CAFAna/FC/FCSurface.h"
13 #include "CAFAna/Systs/Systs.h"
19 #include "CAFAna/Systs/BeamSysts.h"
21 #include "CAFAna/Vars/FitVars.h"
22 #include "OscLib/IOscCalc.h"
23 #include "TFile.h"
24 #include "TCanvas.h"
25 #include "TH2.h"
26 #include "TLegend.h"
27 #include "TLatex.h"
28 
29 using namespace ana;
30 void make_contours(MultiExperiment* expt, osc::IOscCalcAdjustable* calc, std::vector<const ISyst*> pc_syst,
31  FrequentistSurface surf_numu, TH1* prof_numu, int npc)
32 {
33  std::vector<const ISyst*> numuSysts = pc_syst;
34 
35  // ResetOscCalcToDefault(calc);
36  FrequentistSurface surf_numu_syst(expt, calc, &kFitSinSqTheta23, 50, 0, 1,
37  &kFitDmSq32Scaled, 50, 2, 3.5,
39  {{&kFitDeltaInPiUnits,{0.,0.5,1.,1.5}}});
40 
41  TCanvas* cnumu = new TCanvas("numu", "numu", 200, 10, 1600, 600);
42  cnumu->Divide(2,1);
43  cnumu->cd(1);
44  surf_numu.DrawContour(Gaussian90Percent2D(surf_numu), kSolid, kBlue);
45  surf_numu.DrawBestFit(kBlue);
46  surf_numu_syst.DrawContour(Gaussian90Percent2D(surf_numu_syst), kDashed, kRed);
47  surf_numu_syst.DrawBestFit(kRed);
48 
49  TH1* prof_numu_syst = Profile(expt, calc, &kFitSinSqTheta23, 50, 0.3, 0.7, -1,
50  {&kFitDmSq32Scaled, &kFitSinSq2Theta13}, numuSysts);
51  cnumu->cd(2);
52  prof_numu->SetLineColor(kBlue);
53  prof_numu_syst->SetLineColor(kRed);
54  prof_numu_syst->SetLineStyle(kDashed);
55  TLegend* leg_numu = new TLegend(0.5, 0.4, 0.7, 0.8);
56  std::string cv_leg_numu = "CV : "+std::to_string(sqrt(prof_numu->Interpolate(0.5)));
57  std::string pc_leg_numu = "PC"+std::to_string(npc)+" : "+std::to_string(sqrt(prof_numu_syst->Interpolate(0.5)));
58  leg_numu->AddEntry(prof_numu, cv_leg_numu.c_str(), "l");
59  leg_numu->AddEntry(prof_numu_syst, pc_leg_numu.c_str(), "l");
60  prof_numu->Draw("hist");
61  prof_numu->SetTitle("#theta_{23} Sensitivity");
62  prof_numu_syst->Draw("hist same");
63  leg_numu->Draw();
64  std::string numufile = "contours/numu_pc"+std::to_string(npc)+".root";
65  std::string numufile2 = "contours/numu_pc"+std::to_string(npc)+".png";
66  cnumu->Print(numufile.c_str());
67  cnumu->Print(numufile2.c_str());
68 
69  // ResetOscCalcToDefault(calc);
70  // FrequentistSurface surf_nue_syst(expt, calc,
71  // &kFitDeltaInPiUnits,100,0,2,
72  // &kFitSinSqTheta23, 50, 0, 1,
73  // {&kFitDmSq32Scaled, &kFitSinSq2Theta13}, nueSysts);
74  //
75  // TCanvas* cnue = new TCanvas("nue", "nue", 200, 10, 1600, 600);
76  // cnue->Divide(2,1);
77  // cnue->cd(1);
78  // surf_nue.DrawContour(Gaussian90Percent2D(surf_nue), kSolid, kBlue);
79  // surf_nue.DrawBestFit(kBlue);
80  // surf_nue_syst.DrawContour(Gaussian90Percent2D(surf_nue_syst), kDashed, kRed);
81  // surf_nue_syst.DrawBestFit(kRed);
82  //
83  //
84  // TH1* prof_nue_syst = Profile(expt, calc, &kFitDeltaInPiUnits, 100, 0, 2, -1,
85  // {&kFitDmSq32Scaled, &kFitSinSq2Theta13, &kFitSinSqTheta23}, numuSysts);
86  // cnue->cd(2);
87  // prof_nue->SetLineColor(kBlue);
88  // prof_nue_syst->SetLineColor(kRed);
89  // prof_nue_syst->SetLineStyle(kDashed);
90  // TLegend* leg_nue = new TLegend(0.3, 0.4, 0.5, 0.8);
91  // std::string cv_leg_nue = "CV : "+std::to_string(sqrt(prof_nue->Interpolate(0)));
92  // std::string pc_leg_nue = "PC"+std::to_string(npc)+" : "+std::to_string(sqrt(prof_nue_syst->Interpolate(0)));
93  // leg_nue->AddEntry(prof_nue, cv_leg_nue.c_str(), "l");
94  // leg_nue->AddEntry(prof_nue_syst, pc_leg_nue.c_str(), "l");
95  // prof_nue->Draw("hist");
96  // prof_nue->SetTitle("CP Sensitivity");
97  // prof_nue_syst->Draw("hist same");
98  // leg_nue->Draw();
99  // std::string nuefile = "contours/nue_pc"+std::to_string(npc)+".root";
100  // std::string nuefile2 = "contours/nue_pc"+std::to_string(npc)+".png";
101  // cnue->Print(nuefile.c_str());
102  // cnue->Print(nuefile2.c_str());
103 
104 }
105 
106 void make_nom_expt(TFile* predFile, std::string nue_folder, std::string numu_folder, const Spectrum& nue_data, const Spectrum& numu_data)
107 {
108  auto nuePred = ana::LoadFrom<IPrediction>(predFile, nue_folder.c_str());
109  auto numuPred = ana::LoadFrom<IPrediction>(predFile, numu_folder.c_str());
110 
112  const double potFD = 9E20;
113 
114  ResetOscCalcToDefault(calc);
115  calc->SetTh23(asin(sqrt(0.5)));
116  calc->SetDmsq32(2.67e-3);
117  calc->SetdCP(1.5*M_PI);
118  auto nueExpt = new SingleSampleExperiment(nuePred.get(), nue_data);
119  auto numuExpt = new SingleSampleExperiment(numuPred.get(), numu_data);
120  auto reactor = WorldReactorConstraint2016();
121  auto expt = new MultiExperiment({nueExpt, numuExpt, reactor});
122 
123  // ResetOscCalcToDefault(calc);
124  FrequentistSurface surf_numu(expt, calc, &kFitSinSqTheta23, 50, 0, 1,
125  &kFitDmSq32Scaled, 50, 2, 3.5,
127  {{&kFitDeltaInPiUnits, {0.,0.5,1.,1.5}}});
128 
129  surf_numu.DrawContour(Gaussian90Percent2D(surf_numu), kSolid, kGray);
130  surf_numu.DrawBestFit(kGray);
131 
132  gPad->Update();
133 
134 }
135 
136 TH1* make_expt_profile(TFile* predFile, std::string nue_folder, std::string numu_folder, const Spectrum& nue_data, const Spectrum& numu_data)
137 {
138  auto nuePred = ana::LoadFrom<IPrediction>(predFile, nue_folder.c_str());
139  auto numuPred = ana::LoadFrom<IPrediction>(predFile, numu_folder.c_str());
140 
142  const double potFD = 9E20;
143 
144  ResetOscCalcToDefault(calc);
145  calc->SetTh23(asin(sqrt(0.4)));
146  calc->SetDmsq32(2.67e-3);
147  calc->SetdCP(1.5*M_PI);
148  auto nueExpt = new SingleSampleExperiment(nuePred.get(), nue_data);
149  auto numuExpt = new SingleSampleExperiment(numuPred.get(), numu_data);
150  auto reactor = WorldReactorConstraint2016();
151  auto expt = new MultiExperiment({nueExpt, numuExpt, reactor});
152 
153  // ResetOscCalcToDefault(calc);
154  TH1* prof_numu = Profile(expt, calc, &kFitSinSqTheta23, 50, 0.3, 0.7, -1,
156  prof_numu->SetLineColor(kGray+1);
157  return prof_numu;
158 }
159 
160 double get_univ_chisq(TFile* predFile, std::string nue_fhc_folder, std::string numu_fhc_folder, std::string nue_rhc_folder, std::string numu_rhc_folder, const Spectrum& nue_fhc_data, const Spectrum& numu_fhc_data, const Spectrum& nue_rhc_data, const Spectrum& numu_rhc_data, const SystShifts shift, osc::IOscCalcAdjustable* calc)
161 {
162  // auto nuePred_fhc = ana::LoadFrom<PredictionInterp>(predFile, nue_fhc_folder.c_str());
163  // auto numuPred_fhc = ana::LoadFrom<PredictionInterp>(predFile, numu_fhc_folder.c_str());
164  // auto nuePred_rhc = ana::LoadFrom<PredictionInterp>(predFile, nue_rhc_folder.c_str());
165  // auto numuPred_rhc = ana::LoadFrom<PredictionInterp>(predFile, numu_rhc_folder.c_str());
166  auto nuePred_fhc = PredictionInterp::LoadFrom(predFile, nue_fhc_folder.c_str()).release();
167  auto numuPred_fhc = PredictionInterp::LoadFrom(predFile, numu_fhc_folder.c_str()).release();
168  auto nuePred_rhc = PredictionInterp::LoadFrom(predFile, nue_rhc_folder.c_str()).release();
169  auto numuPred_rhc = PredictionInterp::LoadFrom(predFile, numu_rhc_folder.c_str()).release();
170 
171  // osc::IOscCalcAdjustable* calc = DefaultOscCalc();
172  const double potFD_fhc = 8.85E20;
173  const double potFD_rhc = 8.1E20;
174 
175  auto nueExpt_fhc = new SingleSampleExperiment(nuePred_fhc, nue_fhc_data);
176  auto numuExpt_fhc = new SingleSampleExperiment(numuPred_fhc, numu_fhc_data);
177  auto nueExpt_rhc = new SingleSampleExperiment(nuePred_rhc, nue_rhc_data);
178  auto numuExpt_rhc = new SingleSampleExperiment(numuPred_rhc, numu_rhc_data);
179  auto reactor = WorldReactorConstraint2017();
180  auto expt = new MultiExperiment({nueExpt_fhc, numuExpt_fhc, nueExpt_rhc, numuExpt_rhc, reactor});
181 
182  // ResetOscCalcToDefault(calc);
183  return sqrt(expt->ChiSq(calc, shift));
184 }
185 
187 {
188 
189  // osc::IOscCalcAdjustable* calc = DefaultOscCalc();
190  // ResetOscCalcToDefault(calc);
191  return sqrt(expt->ChiSq(calc, shift));
192 }
193 
195 {
196 
197  TString filename1 = "genie_preds_000_499.root";
198  TString filename2 = "genie_preds_500_999.root";
199  std::vector<const ISyst*> pc_no_sup_shifts;
200  std::vector<const ISyst*> pc_fhc_sup_shifts;
201  std::vector<const ISyst*> pc_fhc_ws_sup_shifts;
202  // for(int i = 0; i < 10; i++)
203  // pc_fhc_sup_shifts.push_back(GetPPFXPrincipalsFHCSup(i));
204  // for(int i = 0; i < 10; i++)
205  // pc_fhc_ws_sup_shifts.push_back(GetPPFXPrincipalsFHCWSSup(i)); for(int i = 0; i < 10; i++)
206  for(int i = 0; i < 10; i++)
207  pc_no_sup_shifts.push_back(GetGeniePrincipals2018Small(i));
208  TFile* predFile = new TFile(filename1, "read");
209  TFile* predFile2 = new TFile(filename2, "read");
210 
211  auto nuePred_fhc = PredictionInterp::LoadFrom(predFile, "nue_prop_fhc").release();
212  auto numuPred_fhc = PredictionInterp::LoadFrom(predFile, "numu_fhc").release();
213  auto nuePred_rhc = PredictionInterp::LoadFrom(predFile, "nue_prop_rhc").release();
214  auto numuPred_rhc = PredictionInterp::LoadFrom(predFile, "numu_rhc").release();
215 
217  const double potFD_fhc = 8.85E20;
218  const double potFD_rhc = 8.1E20;
219 
220  ResetOscCalcToDefault(calc);
221  calc->SetTh23(asin(sqrt(0.4)));
222  calc->SetDmsq32(2.44e-3);
223  calc->SetdCP(1.5*M_PI);
224  auto nueFakeData_fhc = nuePred_fhc->Predict(calc).FakeData(potFD_fhc);
225  auto numuFakeData_fhc = numuPred_fhc->Predict(calc).FakeData(potFD_fhc);
226  auto nueExpt_fhc = new SingleSampleExperiment(nuePred_fhc, nueFakeData_fhc);
227  auto numuExpt_fhc = new SingleSampleExperiment(numuPred_fhc, numuFakeData_fhc);
228  auto nueFakeData_rhc = nuePred_rhc->Predict(calc).FakeData(potFD_rhc);
229  auto numuFakeData_rhc = numuPred_rhc->Predict(calc).FakeData(potFD_rhc);
230  auto nueExpt_rhc = new SingleSampleExperiment(nuePred_rhc, nueFakeData_rhc);
231  auto numuExpt_rhc = new SingleSampleExperiment(numuPred_rhc, numuFakeData_rhc);
232  auto reactor = WorldReactorConstraint2017();
233  auto expt = new MultiExperiment({nueExpt_fhc, numuExpt_fhc, nueExpt_rhc, numuExpt_rhc, reactor});
234 
235  /*
236  std::string title = "CP Violation";
237  std::string label = "cp";
238  calc->SetdCP(0);
239  */
240  /*
241  std::string title = "Max-Mixing Rejection";
242  std::string label = "mmix";
243  calc->SetTh23(asin(sqrt(0.5)));
244  */
245 
246  std::string title = "IH Rejection";
247  std::string label = "ih";
248  calc->SetDmsq32(-2.44e-3);
249 
250  // ResetOscCalcToDefault(calc);
251  // for(int npc_tot = 6; npc_tot > 0; npc_tot--){
252 
253  int npc_tot = 5;
254  double rlow = 1.5;
255  double rhigh = 2.5;
256  TH1D* sens_syst_no_sup = new TH1D("sens_syst_no_sup", "sens_syst_no_sup", 200, rlow, rhigh);
257  // TH1D* sens_syst_fhc_sup = new TH1D("sens_syst_fhc_sup", "sens_syst_fhc_sup", 200, 0.75, 1.15);
258  // TH1D* sens_syst_fhc_ws_sup = new TH1D("sens_syst_fhc_ws_sup", "sens_syst_fhc_ws_sup", 200, 0.75, 1.15);
259  gRandom->SetSeed(0);
260  for(int i = 0; i < 10000; i++){
261  SystShifts shift_no_sup;
262  // SystShifts shift_fhc_sup;
263  // SystShifts shift_fhc_ws_sup;
264  // int npc = 0;
265  // for(auto syst: pc_no_sup_shifts){
266  for(int npc = 0; npc < npc_tot; npc++){
267  // if(npc >= npc_tot) continue;
268  double shift = 2.0*gRandom->Gaus();
269  shift_no_sup.SetShift(pc_no_sup_shifts[npc], shift);
270  // shift_fhc_sup.SetShift(pc_fhc_sup_shifts[npc], (1.25/1.)*gRandom->Gaus());
271  // shift_fhc_ws_sup.SetShift(pc_fhc_ws_sup_shifts[npc], (1.25/1.)*gRandom->Gaus());
272  //npc++;
273  }
274  //std::cout << "-----" << std::endl;
275  double chisq_no_sup = get_chisq(expt, shift_no_sup, calc);
276  // double chisq_fhc_sup = get_chisq(expt, shift_fhc_sup, calc);
277  // double chisq_fhc_ws_sup = get_chisq(expt, shift_fhc_ws_sup, calc);
278  sens_syst_no_sup->Fill(chisq_no_sup);
279  // sens_syst_fhc_sup->Fill(chisq_fhc_sup);
280  // sens_syst_fhc_ws_sup->Fill(chisq_fhc_ws_sup);
281  //std::cout << i << " " << chisq_no_sup << std::endl;
282  }
283  TH1D* sensitivity = new TH1D("sensitivity", "sensitivity", 20, rlow, rhigh);
284  for(int i = 0; i < 1000; i++){
285  if( i==499 ) continue;
286  if( i==999 ) continue;
287  double chisq = i < 500 ?
288  get_univ_chisq(predFile, "nue_prop_univ_fhc_"+std::to_string(i), "numu_univ_fhc_"+std::to_string(i), "nue_prop_univ_rhc_"+std::to_string(i), "numu_univ_rhc_"+std::to_string(i), nueFakeData_fhc, numuFakeData_fhc, nueFakeData_rhc, numuFakeData_rhc, kNoShift, calc) :
289  get_univ_chisq(predFile2, "nue_prop_univ_fhc_"+std::to_string(i), "numu_univ_fhc_"+std::to_string(i), "nue_prop_univ_rhc_"+std::to_string(i), "numu_univ_rhc_"+std::to_string(i), nueFakeData_fhc, numuFakeData_fhc, nueFakeData_rhc, numuFakeData_rhc, kNoShift, calc) ;
290  sensitivity->Fill(chisq);
291  //std::cout << i << " " << chisq << std::endl;
292  }
293 
294  TCanvas *c_no_sup = new TCanvas("c_no_sup", "c_no_sup", 200, 10, 800, 600);
295  c_no_sup->cd();
296  sensitivity->Scale(1./sensitivity->Integral());
297  sensitivity->SetMaximum(2.0*sensitivity->GetMaximum());
298  sensitivity->Draw("hist same");
299  sens_syst_no_sup->Scale(sens_syst_no_sup->GetNbinsX()/(sensitivity->GetNbinsX()*sens_syst_no_sup->Integral()));
300  sens_syst_no_sup->Draw("hist same");
301  sensitivity->SetTitle(Form("%s w/ %d PC", title.c_str(), npc_tot));
302  sensitivity->GetXaxis()->SetTitle("#sigma");
303  std::cout << "PC Mean : " << sens_syst_no_sup->GetMean() << std::endl;
304  std::cout << "PC RMS : " << sens_syst_no_sup->GetRMS() << std::endl;
305  std::cout << "PPFX Univ Mean : " << sensitivity->GetMean() << std::endl;
306  std::cout << "PPFX Univ RMS : " << sensitivity->GetRMS() << std::endl;
307  sens_syst_no_sup->SetLineColor(kRed);
308  c_no_sup->Update();
309  c_no_sup->Print(Form("plots/%s_joint_no_sup_%dpc.pdf", label.c_str(), npc_tot));
310 
311  // TCanvas *c_fhc_sup = new TCanvas("c_fhc_sup", "c_fhc_sup", 200, 10, 800, 600);
312  // c_fhc_sup->cd();
313  // sensitivity->Scale(1./sensitivity->Integral());
314  // sensitivity->Draw("hist same");
315  // sens_syst_fhc_sup->Scale(sens_syst_fhc_sup->GetNbinsX()/(sensitivity->GetNbinsX()*sens_syst_fhc_sup->Integral()));
316  // sens_syst_fhc_sup->Draw("hist same");
317  // sensitivity->SetTitle(Form("CP violation w/ %d PC", npc_tot));
318  // sensitivity->GetXaxis()->SetTitle("#sigma");
319  // std::cout << "PC Mean : " << sens_syst_fhc_sup->GetMean() << std::endl;
320  // std::cout << "PC RMS : " << sens_syst_fhc_sup->GetRMS() << std::endl;
321  // std::cout << "PPFX Univ Mean : " << sensitivity->GetMean() << std::endl;
322  // std::cout << "PPFX Univ RMS : " << sensitivity->GetRMS() << std::endl;
323  // sens_syst_fhc_sup->SetLineColor(kRed);
324  // c_fhc_sup->Update();
325  // c_fhc_sup->Print(Form("deltacp_fhc_only_fhc_sup_%dpc.pdf", npc_tot));
326  //
327  // TCanvas *c_fhc_ws_sup = new TCanvas("c_fhc_ws_sup", "c_fhc_ws_sup", 200, 10, 800, 600);
328  // c_fhc_ws_sup->cd();
329  // sensitivity->Scale(1./sensitivity->Integral());
330  // sensitivity->Draw("hist same");
331  // sens_syst_fhc_ws_sup->Scale(sens_syst_fhc_ws_sup->GetNbinsX()/(sensitivity->GetNbinsX()*sens_syst_fhc_ws_sup->Integral()));
332  // sens_syst_fhc_ws_sup->Draw("hist same");
333  // sensitivity->SetTitle(Form("CP violation w/ %d PC", npc_tot));
334  // sensitivity->GetXaxis()->SetTitle("#sigma");
335  // std::cout << "PC Mean : " << sens_syst_fhc_ws_sup->GetMean() << std::endl;
336  // std::cout << "PC RMS : " << sens_syst_fhc_ws_sup->GetRMS() << std::endl;
337  // std::cout << "PPFX Univ Mean : " << sensitivity->GetMean() << std::endl;
338  // std::cout << "PPFX Univ RMS : " << sensitivity->GetRMS() << std::endl;
339  // sens_syst_fhc_ws_sup->SetLineColor(kRed);
340  // c_fhc_ws_sup->Update();
341  // c_fhc_ws_sup->Print(Form("deltacp_fhc_only_fhc_ws_sup_%dpc.pdf", npc_tot));
342  // }
343 
344 
345 }
virtual double ChiSq(osc::IOscCalcAdjustable *osc, const SystShifts &syst=SystShifts::Nominal()) const
Definition: IExperiment.h:18
enum BeamMode kRed
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
void ResetOscCalcToDefault(osc::IOscCalcAdjustable *calc)
Reset calculator to default assumptions for all parameters.
Definition: Calcs.cxx:23
Simple record of shifts applied to systematic parameters.
Definition: SystShifts.h:20
T sqrt(T number)
Definition: d0nt_math.hpp:156
double get_univ_chisq(TFile *predFileNue, TFile *predFileNumu, std::string nue_fhc_folder, std::string numu_fhc_folder, std::string nue_rhc_folder, std::string numu_rhc_folder, const Spectrum &nue_fhc_data, const Spectrum &numu_fhc_data, const Spectrum &nue_rhc_data, const Spectrum &numu_rhc_data, const SystShifts shift, osc::IOscCalcAdjustable *calc)
TGraph * Profile(const IExperiment *expt, osc::IOscCalcAdjustable *calc, const IFitVar *v, int nbinsx, double minx, double maxx, double input_minchi, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts, const SeedList &seedPts, const std::vector< SystShifts > &systSeedPts, std::map< const IFitVar *, TGraph * > &profVarsMap, std::map< const ISyst *, TGraph * > &profSystsMap, MinuitFitter::FitOpts opts)
scan in one variable, profiling over all others
Definition: Fit.cxx:48
osc::OscCalcDumb calc
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
Definition: Calcs.cxx:49
virtual void SetDmsq32(const T &dmsq32)=0
#define M_PI
Definition: SbMath.h:34
void genie_contours(int npc=5)
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
void make_contours(MultiExperiment *expt, osc::IOscCalcAdjustable *calc, std::vector< const ISyst * > pc_syst, FrequentistSurface surf_numu, TH1 *prof_numu, int npc)
const char * label
Log-likelihood scan across two parameters.
void DrawBestFit(Color_t color, Int_t marker=kFullCircle) const
Draw the best fit point.
Definition: ISurface.cxx:32
expt
Definition: demo5.py:34
void make_nom_expt(TFile *predFile, std::string nue_folder, std::string numu_folder, const Spectrum &nue_data, const Spectrum &numu_data)
void DrawContour(TH2 *fc, Style_t style, Color_t color, double minchi=-1)
Definition: ISurface.cxx:44
const ReactorExperiment * WorldReactorConstraint2016()
Updated value for SecondAna based on the latest PDG.
const ReactorExperiment * WorldReactorConstraint2017()
Reactor constraint from PDG 2017.
const ConstrainedFitVarWithPrior fitDmSq32Scaled_UniformPrior & kFitDmSq32Scaled
std::unique_ptr< IPrediction > LoadFrom< IPrediction >(TDirectory *dir, const std::string &label)
Definition: IPrediction.cxx:18
GeniePCASyst * GetGeniePrincipals2018Small(int PCIdx)
const SystShifts kNoShift
Definition: SystShifts.cxx:21
OStream cout
Definition: OStream.cxx:6
double get_chisq(IExperiment *expt, const SystShifts shift, osc::IOscCalcAdjustable *calc)
TH1 * make_expt_profile(TFile *predFile, std::string nue_folder, std::string numu_folder, const Spectrum &nue_data, const Spectrum &numu_data)
Combine multiple component experiments.
static std::unique_ptr< PredictionInterp > LoadFrom(TDirectory *dir, const std::string &name)
TH2 * Gaussian90Percent2D(const FrequentistSurface &s)
Up-value surface for 90% confidence in 2D in gaussian approximation.
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
virtual void SetTh23(const T &th23)=0
Base class defining interface for experiments.
Definition: IExperiment.h:14
const FitSinSq2Theta13 kFitSinSq2Theta13
Definition: FitVars.cxx:13
std::string to_string(ModuleType mt)
Definition: ModuleType.h:32
Float_t e
Definition: plot.C:35
const FitVarWithPrior fitDeltaInPiUnits_UniformPriordCP & kFitDeltaInPiUnits
enum BeamMode kBlue
void SetShift(const ISyst *syst, double shift, bool force=false)
Definition: SystShifts.cxx:80
virtual void SetdCP(const T &dCP)=0
T asin(T number)
Definition: d0nt_math.hpp:60
Compare a single data spectrum to the MC + cosmics expectation.
enum BeamMode string