saveContours_addExpt.C
Go to the documentation of this file.
1 #ifdef __CINT__
3  std::cout << "Sorry, you must run in compiled mode" << std::endl;
4 }
5 #else
6 
7 #include "../WorkshopIncludes.h"
8 
10 {
11  TFile* fPred_epoch3b = new TFile("../Predictions/Prediction_epoch3b_simple_kCCE.root");
12  PredictionNoExtrap* prediction_epoch3b = LoadFrom<PredictionNoExtrap>(fPred_epoch3b->GetDirectory("pred_3b")).release();
13  fPred_epoch3b->Close();
14 
15  TFile* fPred_addExpt = new TFile("../Predictions/Prediction_epoch3b_addExpt_kCCE.root");
16  PredictionNoExtrap* prediction_ugly = LoadFrom<PredictionNoExtrap>(fPred_addExpt->GetDirectory("pred_ugly")).release();
17  PredictionNoExtrap* prediction_bad = LoadFrom<PredictionNoExtrap>(fPred_addExpt->GetDirectory("pred_bad")).release();
18  PredictionNoExtrap* prediction_good = LoadFrom<PredictionNoExtrap>(fPred_addExpt->GetDirectory("pred_good")).release();
19 
20  fPred_addExpt->Close();
21 
22  osc::OscCalculatorPMNSOpt calc;
23  ana::ResetOscCalcToDefault(&calc); // This seeds to default (in particular max-mixing)
24 
25  calc.SetTh23(TMath::ASin(sqrt(0.403))); // Second analysis best fit
26  calc.SetDmsq32(2.67e-3); // Same
27 
28  Spectrum fakePred = prediction_epoch3b->Predict(&calc);
29  Spectrum fakeData = fakePred.FakeData(kSecondAnaPOT);
30 
31  // Split experiments data
32  Spectrum fakePred_ugly = prediction_ugly->Predict(&calc);
33  Spectrum fakeData_ugly = fakePred_ugly.FakeData(kSecondAnaPOT);
34 
35  Spectrum fakePred_bad = prediction_bad->Predict(&calc);
36  Spectrum fakeData_bad = fakePred_bad.FakeData(kSecondAnaPOT);
37 
38  Spectrum fakePred_good = prediction_good->Predict(&calc);
39  Spectrum fakeData_good = fakePred_good.FakeData(kSecondAnaPOT);
40 
41  // Experiments
42 
43  SingleSampleExperiment fakeExpt(prediction_epoch3b, fakeData);
44 
45  // Multiple exprts
46  SingleSampleExperiment fakeExpt_ugly(prediction_ugly, fakeData_ugly);
47  SingleSampleExperiment fakeExpt_bad(prediction_bad, fakeData_bad);
48  SingleSampleExperiment fakeExpt_good(prediction_good, fakeData_good);
49  /// Combine experiments
50  std::vector<const IExperiment*> vExpts;
51  vExpts.push_back(&fakeExpt_ugly);
52  vExpts.push_back(&fakeExpt_bad);
53  vExpts.push_back(&fakeExpt_good);
54 
55  MultiExperiment multiExpt(vExpts);
56 
57  // 2D Fits
58 
59  Surface* surface = new Surface( &fakeExpt, &calc,
60  &kFitSinSqTheta23, 40, 0.3, 0.7,
61  &kFitDmSq32Scaled, 35, 2.3, 3.0);
62 
63  Surface* surfaceExtra = new Surface( &multiExpt, &calc,
64  &kFitSinSqTheta23, 40, 0.3, 0.7,
65  &kFitDmSq32Scaled, 35, 2.3, 3.0);
66 
67  // Plot
68  TCanvas* c = new TCanvas();
69  c->SetLeftMargin(0.12);
70  c->SetBottomMargin(0.15);
71 
72  surface->DrawBestFit(kBlue);
73  surface->DrawContour( Gaussian90Percent2D(*surface), kSolid, kBlue);
74 
75  surfaceExtra->DrawBestFit(kRed);
76  surfaceExtra->DrawContour( Gaussian90Percent2D(*surfaceExtra), kDashed, kRed);
77 
78  c->SaveAs("../Contours/ExtraExpt.root");
79  c->SaveAs("../Contours/ExtraExpt.pdf");
80 
81  c->Clear();
82 
83  // Get 1D marginalised plots
84  TH1* profileSinsq = Profile(&fakeExpt, &calc,
85  &kFitSinSqTheta23, 40, 0.3, 0.7, -1,
86  {&kFitDmSq32Scaled} );
87 
88  TH1* profileSinsqMulti = Profile(&multiExpt, &calc,
89  &kFitSinSqTheta23, 40, 0.3, 0.7, -1,
90  {&kFitDmSq32Scaled} );
91 
92  profileSinsq->SetLineColor(kBlue);
93  profileSinsq->Draw("hist");
94  profileSinsqMulti->SetLineColor(kRed);
95  profileSinsqMulti->SetLineStyle(7);
96  profileSinsqMulti->Draw("hist same");
97 
98  std::cout << "Max-mix rej 1 : " << sqrt(profileSinsq->Interpolate(0.5)) << std::endl;
99  std::cout << "Max-mix rej Multi: " << sqrt(profileSinsqMulti->Interpolate(0.5)) << std::endl;
100 
101  c->SaveAs("../Contours/1D_addExpt.root");
102  c->SaveAs("../Contours/1D_addExpt.pdf");
103 
104  c->Close();
105 
106 } // End of function
107 
108 
109 #endif
enum BeamMode kRed
void ResetOscCalcToDefault(osc::IOscCalcAdjustable *calc)
Reset calculator to default assumptions for all parameters.
Definition: Calcs.cxx:23
T sqrt(T number)
Definition: d0nt_math.hpp:156
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
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
Spectrum Predict(osc::IOscCalc *calc) const override
Spectrum FakeData(double pot) const
Synonymous with AsimovData(). Retained for compatibility.
Definition: Spectrum.cxx:349
const ConstrainedFitVarWithPrior fitDmSq32Scaled_UniformPrior & kFitDmSq32Scaled
OStream cout
Definition: OStream.cxx:6
Combine multiple component experiments.
TH2 * Gaussian90Percent2D(const FrequentistSurface &s)
Up-value surface for 90% confidence in 2D in gaussian approximation.
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
Prediction that just uses FD MC, with no extrapolation.
void saveContours_addExpt()
Float_t e
Definition: plot.C:35
const double kSecondAnaPOT
Definition: Exposures.h:87
enum BeamMode kBlue
Compare a single data spectrum to the MC + cosmics expectation.