35 SystName =
"WithSysts";
38 if(usesysts ==
false){
40 SystName =
"StatsOnly";
50 if(calculator ==
"sa"){
52 CalcName =
"SABestFit";
56 if(calculator ==
"3a"){
58 CalcName =
"3ABestFit";
73 pot = kPeriod1POT + kPeriod2POT + kPeriod3POT + kPeriod5POT;
80 std::string inDirData =
"/nova/ana/nu_mu_ana/Ana2017/Data";
82 std::string inName =
"Prediction_StatsOnly_3ACut_3AEnergy_OptBinning_";
84 TFile* fDataQ0 =
new TFile( (inDirData +
"/fd_data_3ACut_3AEnergy_OptBinning_full_" + SampleName +
"_Quant0.root").c_str() );
85 TFile* fDataQ1 =
new TFile( (inDirData +
"/fd_data_3ACut_3AEnergy_OptBinning_full_" + SampleName +
"_Quant1.root").c_str() );
86 TFile* fDataQ2 =
new TFile( (inDirData +
"/fd_data_3ACut_3AEnergy_OptBinning_full_" + SampleName +
"_Quant2.root").c_str() );
87 TFile* fDataQ3 =
new TFile( (inDirData +
"/fd_data_3ACut_3AEnergy_OptBinning_full_" + SampleName +
"_Quant3.root").c_str() );
88 TFile* fDataQ4 =
new TFile( (inDirData +
"/fd_data_3ACut_3AEnergy_OptBinning_full_" + SampleName +
"_Quant4.root").c_str() );
89 Spectrum* FD_Quant0 = LoadFrom< Spectrum >(fDataQ0,
"fd_data").
release();
90 Spectrum* FD_Quant1 = LoadFrom< Spectrum >(fDataQ1,
"fd_data").
release();
91 Spectrum* FD_Quant2 = LoadFrom< Spectrum >(fDataQ2,
"fd_data").
release();
92 Spectrum* FD_Quant3 = LoadFrom< Spectrum >(fDataQ3,
"fd_data").
release();
93 Spectrum* FD_Quant4 = LoadFrom< Spectrum >(fDataQ4,
"fd_data").
release();
94 TH1* hData = FD_Quant0->
ToTH1(kHereSAPOT);
117 CosmicsSpecFile =
"Cosmics/cosmicsScaled_NewBinning_" + SampleName +
".root";
118 CosmicsHistName =
"cosmics3A";
119 cosmicsName =
"WithCosmics";
125 TFile inCosmicsFile(CosmicsSpecFile.c_str());
126 auto *CosQuant1 = (TH1*)inCosmicsFile.Get(
"q1all");
127 auto *CosQuant2 = (TH1*)inCosmicsFile.Get(
"q2all");
128 auto *CosQuant3 = (TH1*)inCosmicsFile.Get(
"q3all");
129 auto *CosQuant4 = (TH1*)inCosmicsFile.Get(
"q4all");
130 auto *CosQuant0 = (TH1*)CosQuant1->Clone(
"CosQuant0");
131 CosQuant0->Add(CosQuant2);
132 CosQuant0->Add(CosQuant3);
133 CosQuant0->Add(CosQuant4);
135 Spectrum cosmics1(CosQuant1, pot, livetime);
136 Spectrum cosmics2(CosQuant2, pot, livetime);
137 Spectrum cosmics3(CosQuant3, pot, livetime);
138 Spectrum cosmics4(CosQuant4, pot, livetime);
165 std::vector<PredictionInterp*> predictionVec;
170 predictionVec.push_back(&predQ1);
171 predictionVec.push_back(&predQ2);
172 predictionVec.push_back(&predQ3);
173 predictionVec.push_back(&predQ4);
177 for(
int epochId = 1; epochId <
nquant; epochId++){
181 hPred->Add(CosQuant0);
186 Spectrum fakePred_quant1 = predictionVec[0]->Predict(calc);
187 Spectrum fakePred_quant2 = predictionVec[1]->Predict(calc);
188 Spectrum fakePred_quant3 = predictionVec[2]->Predict(calc);
189 Spectrum fakePred_quant4 = predictionVec[3]->Predict(calc);
195 fakeData_quant1 += cosmics1;
196 fakeData_quant2 += cosmics2;
197 fakeData_quant3 += cosmics3;
198 fakeData_quant4 += cosmics4;
204 predictionVec.clear();
206 std::vector <const IExperiment*> vExpts;
207 vExpts.push_back(&fakeExpt_quant1);
208 vExpts.push_back(&fakeExpt_quant2);
209 vExpts.push_back(&fakeExpt_quant3);
210 vExpts.push_back(&fakeExpt_quant4);
215 TFile*
fout =
new TFile(outDir +
"ContProfSens_" + SystName +
"_" + CalcName +
"_3ACut_3AEnergy_OptBinning_" + SampleName +
"_Cosmics.root",
"RECREATE");
217 FrequentistSurface* surface =
new FrequentistSurface( &multExpt, calc, &
kFitSinSqTheta23, 40, 0.3, 0.7, &
kFitDmSq32Scaled, 50, 2.0, 3.0, {&
kFitDeltaInPiUnits,&
kFitSinSq2Theta13},
systs);
221 surface->
SaveTo(fout,
"surface");
222 hProfileTheta23->Write(
"profiletheta23");
223 hProfileDelta32->Write(
"profiledelta32");
224 hData->Write(
"Data");
225 hPred->Write(
"Pred");
Cuts and Vars for the 2020 FD DiF Study.
TH1D * ToTH1(double exposure, Color_t col=kBlack, Style_t style=kSolid, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Histogram made from this Spectrum, scaled to some exposure.
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
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.
Log-likelihood scan across two parameters.
Charged-current interactions.
const double kAna2017Livetime
Spectrum FakeData(double pot) const
Synonymous with AsimovData(). Retained for compatibility.
Loads shifted spectra from files.
const ReactorExperiment * WorldReactorConstraint2017()
Reactor constraint from PDG 2017.
const ConstrainedFitVarWithPrior fitDmSq32Scaled_UniformPrior & kFitDmSq32Scaled
Combine multiple component experiments.
const SolarConstraints kSolarConstraintsPDG2017(7.53e-5, 0.18e-5, 0.851, 0.020)
void SaveTo(TDirectory *dir, const std::string &name) const
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
virtual void SetTh23(const T &th23)=0
Neutral-current interactions.
Both neutrinos and antineutrinos.
const SpillCut kStandardSpillCuts
Apply this unless you're doing something special.
const FitSinSq2Theta13 kFitSinSq2Theta13
Template for Cut and SpillCut.
void getContProf_Sensitivity(int sampleCut=9, std::string calculator="3a", bool usesysts=true, bool cosmics=true)
std::vector< const ISyst * > getAllAna2017Systs(const EAnaType2017 ana, const bool smallgenie)
All neutrinos, any flavor.
const FitVarWithPrior fitDeltaInPiUnits_UniformPriordCP & kFitDeltaInPiUnits
Compare a single data spectrum to the MC + cosmics expectation.