56 bool corrSysts =
false;
58 bool th13Surf =
false;
59 bool dmsqSurf =
false;
60 bool numuOnly =
false;
65 bool fakedata =
false;
70 if(fakedata && joint && onlyNH) savefile =
"sensitivity_numu_nue_joint_NH";
71 if(fakedata && joint && onlyIH) savefile =
"sensitivity_numu_nue_joint_IH";
72 if(fakedata && nueOnly && onlyNH) savefile =
"sensitivity_nueOnly_NH";
73 if(fakedata && nueOnly && onlyIH) savefile =
"sensitivity_nueOnly_IH";
74 if(fakedata && numuOnly && onlyNH) savefile =
"sensitivity_numuOnly_NH";
75 if(fakedata && numuOnly && onlyIH) savefile =
"sensitivity_numuOnly_IH";
76 if(realdata && joint && onlyNH) savefile =
"contours_numu_nue_joint_NH";
77 if(realdata && joint && onlyIH) savefile =
"contours_numu_nue_joint_IH";
78 if(realdata && nueOnly && onlyNH) savefile =
"contours_nueOnly_NH";
79 if(realdata && nueOnly && onlyIH) savefile =
"contours_nueOnly_IH";
80 if(realdata && numuOnly && onlyNH) savefile =
"contours_numuOnly_NH";
81 if(realdata && numuOnly && onlyIH) savefile =
"contours_numuOnly_IH";
85 TString
outdir =
"/nova/app/users/prabhjot/ReadCAFMakeEventList/";
86 TString
outfilename (outdir +
"hist_contours_2017_" + decomp +
"_" + savefile +
"_2017");
89 TFile *
outfile =
new TFile(outfilename+
".root",
"recreate");
95 std::vector <const IPrediction * > preds;
96 std::vector <std::pair <TH1D*, double > >
cosmics;
97 std::vector <Spectrum * >
data;
98 std::vector <const IExperiment * > expts;
102 if(fakedata && onlyNH)
SetFakeCalc(calc_fake, 0.477, 2.43
e-3, 1.432);
103 else if(fakedata && onlyIH)
SetFakeCalc(calc_fake, 0.564422, -2.498
e-3, 1.45913);
110 data.insert(data.end(), nue_data);
118 preds.insert(preds.end(), numu_preds.begin(), numu_preds.end());
120 cosmics.insert(cosmics.end(), numu_cosmics.begin(), numu_cosmics.end());
123 data.insert(data.end(), numu_data.begin(), numu_data.end());
129 for(
int i = 0;
i <
int(preds.size()); ++
i){
132 if(realdata) thisdata = data[
i];
141 cosmics[
i].first->SetMarkerStyle(34);
143 cosmics[
i].first->Draw(
"hist p same");
144 gPad->Print(outdir +
"debug_predictions.pdf");
154 std::cout <<
"Adding WorldReactorConstraint2017\n";
159 std::cout <<
"Adding Dmsq32ConstraintPDG2017\n";
163 std::cout <<
"Creating Multiexperiment with a total of " << expts.size() <<
" experiments\n\n" <<
std::endl;
169 std::vector<const ISyst*>
systs = {};
173 std::vector <SystShifts> seedShifts = {};
186 double minchi23 = 1E20;
188 for(
double thseed : {0.3, 0.7}){
190 for(
int hie:{-1, 1}){
192 std::cout <<
"\n\nFinding best fit " << (thseed < 0.5 ?
"LO " :
"UO ") << (hie > 0 ?
"NH " :
"IH ") <<
"\n\n";
195 auto thisminchi = fit23.Fit(calc, auxShifts ,{
202 minchi23 =
min(minchi23, thisminchi);
216 v.Write(
"overall_minchi");
219 bests[0]= calc_fake->GetdCP();
220 bests[1]= calc_fake->GetTh13();
221 bests[2]= calc_fake->GetTh23();
222 bests[3]= calc_fake->GetDmsq32();
223 bests.Write(
"best_dcp_th13_th23_dmsq32");
227 std::cout <<
"Using 30 bins for dmsq surface\n\n";
230 for(
int hie: {-1, +1}){
232 if((onlyNH && hie < 0) || (onlyIH && hie > 0))
continue;
234 std::cout <<
"Starting surface " << (hie > 0?
"NH ":
"IH") <<
"\n\n";
237 if(!th13Surf && ! dmsqSurf){
248 auto surf1 = surf23.
ToTH2(minchi23);
253 surf1->Write((TString)
"delta_th23_"+hieStr);
254 surf23.SaveTo(outfile, (TString)
"surf_delta_th23_"+hieStr);
256 std::vector<TH2*> profhists = surf23.GetProfiledHists();
259 profhists[0]->Write((hieStr+
"_DmSq32").c_str());
260 profhists[1]->Write((hieStr+
"_SinSq2Theta13").c_str());
264 if(!th13Surf && dmsqSurf){
271 &
kFitDmSq32, steps, (hie > 0 ? 2
e-3 : -3.0
e-3 ), (hie > 0 ? 3.0
e-3 : -2
e-3),
275 auto surf6=surf23m.
ToTH2(minchi23);
279 surf6->Write((TString)
"th23_dm32_"+(hie > 0 ?
"NH" :
"IH"));
281 surf23m.SaveTo(outfile, (TString)
"surf_th23_dm32_"+(hie > 0 ?
"NH" :
"IH"));
291 &kFitSinSq2Theta13, steps, 0, 0.35,
296 auto surf4 = surf13.
ToTH2(-1);
300 surf4->Write((TString)
"th13_delta_"+(hie > 0 ?
"NH":
"IH"));
301 surf13.SaveTo(outfile, (TString)
"surf_th13_delta_"+(hie > 0 ?
"NH" :
"IH"));
Cuts and Vars for the 2020 FD DiF Study.
void ResetOscCalcToDefault(osc::IOscCalcAdjustable *calc)
Reset calculator to default assumptions for all parameters.
fvar< T > fabs(const fvar< T > &x)
const Color_t kCosmicBackgroundColor
const FitDmSq32 kFitDmSq32
Simple record of shifts applied to systematic parameters.
const Dmsq32Constraint kDmsq32ConstraintPDG2017(2.45e-3, 0.05e-3, 2.52e-3, 0.05e-3)
static SystShifts Nominal()
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
virtual void SetDmsq32(const T &dmsq32)=0
Representation of a spectrum in any variable, with associated POT.
const XML_Char const XML_Char * data
Log-likelihood scan across two parameters.
string outfilename
knobs that need extra care
const ReactorExperiment * WorldReactorConstraint2017()
Reactor constraint from PDG 2017.
virtual T GetDmsq32() const
static float min(const float a, const float b, const float c)
TH2 * ToTH2(double minchi=-1) const
Combine multiple component experiments.
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
const FitSinSq2Theta13 kFitSinSq2Theta13
TH1 * DataMCComparison(const Spectrum &data, const Spectrum &mc, EBinType bintype)
const FitVarWithPrior fitDeltaInPiUnits_UniformPriordCP & kFitDeltaInPiUnits
void CAF_makeCAFSensitivities_for_FNEX()
Compare a single data spectrum to the MC + cosmics expectation.
Perform MINUIT fits in one or two dimensions.