52 outfilename =
"histograms_forFNEX_CVN_prop_nueOnly_stats_no_rock_muons_realdata_IH.root";
55 TFile *
infile =
new TFile(
"/nova/ana/nu_e_ana/Ana2017/Predictions/provisional_singles/latest/other/pred_xp_prop_EnergyCVN2D_full_nueconcat_Nominal.root",
"read");
57 TFile *
outfile =
new TFile(outfilename,
"recreate");
61 info +=
"2017 analysis";
62 info +=
"Nue prediction uses CVN and proportional decomposition";
63 info +=
"ND MC scaled to data POT";
64 info +=
"Cosmic spectrum is scaled to data livetime \n";
65 info +=
"FD predicted components are oscillated and scaled to FD POT \n";
66 info +=
"Reactor constraint 0.085 \\pm 0.05\n";
67 info +=
"PDG dmsq32 constraint 2.44 \\pm 0 .06 e-3 and -2.49 \\pm 0.06 e-3 \n";
68 info +=
"Fit delta, ssth23, dmsq32, ss2th13, no systematics \n";
69 info +=
"NH IH surfaces wrt best fit point overall \n";
70 info +=
"POT and Best fit values stored in TVectorD's \n";
71 TNamed readme (
"README",info.Data());
75 auto data =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/data_comp");
76 double ndpot =
data->POT();
80 data->ToTH1(ndpot)->Write(
"nd_data");
82 outfile->mkdir(
"nd_mc");
86 auto nd_mc_nue_comp =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/nue_comp" );
87 auto nd_mc_antinue_comp =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/antinue_comp" );
88 auto nd_mc_numu_comp =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/numu_comp" );
89 auto nd_mc_antinumu_comp =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/antinumu_comp" );
90 auto nd_mc_nc_comp =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/nc_comp" );
91 auto nd_mc_total_comp =
Spectrum::LoadFrom(infile,
"pred_xp_prop_Nominal/noShift/nue_pred_EnergyCVN2D/extrap/EEextrap/Decomp/total_comp" );
93 TH1D* h_nd_mc_nue = nd_mc_nue_comp ->ToTH1(ndpot);
94 TH1D* h_nd_mc_anue = nd_mc_antinue_comp ->ToTH1(ndpot);
95 TH1D* h_nd_mc_numu = nd_mc_numu_comp ->ToTH1(ndpot);
96 TH1D* h_nd_mc_anumu = nd_mc_antinumu_comp ->ToTH1(ndpot);
97 TH1D* h_nd_mc_nc = nd_mc_nc_comp ->ToTH1(ndpot);
98 TH1D* h_nd_mc_total = (TH1D*)h_nd_mc_nue ->
Clone(
"h_nd_mc_total");
99 h_nd_mc_total ->Add(h_nd_mc_anue);
100 h_nd_mc_total ->Add(h_nd_mc_numu);
101 h_nd_mc_total ->Add(h_nd_mc_anumu);
102 h_nd_mc_total ->Add(h_nd_mc_nc);
104 h_nd_mc_nue ->Write(
"nue" );
105 h_nd_mc_anue ->Write(
"anue" );
106 h_nd_mc_numu ->Write(
"numu" );
107 h_nd_mc_anumu ->Write(
"anumu" );
108 h_nd_mc_nc ->Write(
"nc" );
109 h_nd_mc_total ->Write(
"total" );
113 Spectrum*
intime = LoadFromFile<Spectrum>(
"/nova/ana/nu_e_ana/Ana2017/data/data_spectra_fd.root",
"All/spect_merged_4bin").
release();
114 Spectrum*
outtime = LoadFromFile<Spectrum>(
"/nova/ana/nu_e_ana/Ana2017/Predictions/cosmic/cosmic_spectra_prediction.root",
"All/cosm_merged_4bin").
release();
118 double cosmicscale = innN/outN;
120 double fdpot = intime->
POT();
124 intime ->
ToTH1(fdpot) ->Write(
"fd_data");
130 <<
"Calling GetNuePrediction2017 function in the getHists_FNEX.C CAF macro for - " 146 double minchi23 = 1E20;
148 auto bestcalc = calc->
Copy();
149 for(
double ssth = 0.4; ssth < 0.6; ssth += 0.05){
150 for(
int hie : {+1, -1}){
155 double thisminchi = fit23.Fit(calc)->EvalMetricVal();
156 if(thisminchi < minchi23){
157 minchi23 = thisminchi;
158 bestcalc = calc->
Copy();
160 const double dcp = calc->
GetdCP();
161 for(
int n = 1;
n <= 3; ++
n){
163 thisminchi = fit23.Fit(calc)->EvalMetricVal();
164 if(thisminchi<minchi23){
165 minchi23 = thisminchi;
166 bestcalc = calc->
Copy();
174 for(
int hie: {+1, -1}){
184 auto surf = surf23.ToTH2(minchi23);
186 surf->Write((TString)
"dchisq_"+(hie > 0 ?
"NH" :
"IH"));
196 outfile->mkdir(
"fd_extrap_oscil_BF");
197 outfile->cd(
"fd_extrap_oscil_BF");
209 TH1D* h_fd_mc_nue_app = fd_mc_nue_app . ToTH1(fdpot);
210 TH1D* h_fd_mc_nue_surv = fd_mc_nue_surv . ToTH1(fdpot);
211 TH1D* h_fd_mc_numu_surv = fd_mc_numu_surv . ToTH1(fdpot);
212 TH1D* h_fd_mc_numu_app = fd_mc_numu_app . ToTH1(fdpot);
213 TH1D* h_fd_mc_anue_app = fd_mc_anue_app . ToTH1(fdpot);
214 TH1D* h_fd_mc_anue_surv = fd_mc_anue_surv . ToTH1(fdpot);
215 TH1D* h_fd_mc_anumu_surv = fd_mc_anumu_surv . ToTH1(fdpot);
216 TH1D* h_fd_mc_anumu_app = fd_mc_anumu_app . ToTH1(fdpot);
217 TH1D* h_fd_mc_tau_cc_all = fd_mc_tau_cc_all . ToTH1(fdpot);
218 TH1D* h_fd_mc_nc = fd_mc_nc . ToTH1(fdpot);
220 TH1D* h_fd_mc_total = (TH1D*)h_fd_mc_nue_app ->
Clone(
"h_fd_mc_total");
221 h_fd_mc_total ->Add(h_fd_mc_nue_surv );
222 h_fd_mc_total ->Add(h_fd_mc_numu_surv );
223 h_fd_mc_total ->Add(h_fd_mc_numu_app );
224 h_fd_mc_total ->Add(h_fd_mc_anue_app );
225 h_fd_mc_total ->Add(h_fd_mc_anue_surv );
226 h_fd_mc_total ->Add(h_fd_mc_anumu_surv);
227 h_fd_mc_total ->Add(h_fd_mc_anumu_app );
228 h_fd_mc_total ->Add(h_fd_mc_tau_cc_all);
229 h_fd_mc_total ->Add(h_fd_mc_nc );
231 h_fd_mc_nue_app ->Write(
"nue_app" );
232 h_fd_mc_nue_surv ->Write(
"nue_surv" );
233 h_fd_mc_numu_surv ->Write(
"numu_surv" );
234 h_fd_mc_numu_app ->Write(
"numu_app" );
235 h_fd_mc_anue_app ->Write(
"anue_app" );
236 h_fd_mc_anue_surv ->Write(
"anue_surv" );
237 h_fd_mc_anumu_surv ->Write(
"anumu_surv" );
238 h_fd_mc_anumu_app ->Write(
"anumu_app" );
239 h_fd_mc_tau_cc_all ->Write(
"tau_cc_all" );
240 h_fd_mc_nc ->Write(
"nc" );
241 h_fd_mc_total ->Write(
"total" );
248 outfile->mkdir(
"fd_noextrap_oscil_BF");
249 outfile->cd(
"fd_noextrap_oscil_BF");
261 TH1D* noextrap_h_fd_mc_nue_app = noextrap_fd_mc_nue_app . ToTH1(fdpot);
262 TH1D* noextrap_h_fd_mc_nue_surv = noextrap_fd_mc_nue_surv . ToTH1(fdpot);
263 TH1D* noextrap_h_fd_mc_numu_surv = noextrap_fd_mc_numu_surv . ToTH1(fdpot);
264 TH1D* noextrap_h_fd_mc_numu_app = noextrap_fd_mc_numu_app . ToTH1(fdpot);
265 TH1D* noextrap_h_fd_mc_anue_app = noextrap_fd_mc_anue_app . ToTH1(fdpot);
266 TH1D* noextrap_h_fd_mc_anue_surv = noextrap_fd_mc_anue_surv . ToTH1(fdpot);
267 TH1D* noextrap_h_fd_mc_anumu_surv = noextrap_fd_mc_anumu_surv . ToTH1(fdpot);
268 TH1D* noextrap_h_fd_mc_anumu_app = noextrap_fd_mc_anumu_app . ToTH1(fdpot);
269 TH1D* noextrap_h_fd_mc_tau_cc_all = noextrap_fd_mc_tau_cc_all . ToTH1(fdpot);
270 TH1D* noextrap_h_fd_mc_nc = noextrap_fd_mc_nc . ToTH1(fdpot);
272 TH1D* noextrap_h_fd_mc_total = (TH1D*)noextrap_h_fd_mc_nue_app ->
Clone(
"noextrap_h_fd_mc_total");
273 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_nue_surv );
274 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_numu_surv );
275 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_numu_app );
276 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_anue_app );
277 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_anue_surv );
278 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_anumu_surv);
279 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_anumu_app );
280 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_tau_cc_all);
281 noextrap_h_fd_mc_total ->Add(noextrap_h_fd_mc_nc );
283 noextrap_h_fd_mc_nue_app ->Write(
"nue_app" );
284 noextrap_h_fd_mc_nue_surv ->Write(
"nue_surv" );
285 noextrap_h_fd_mc_numu_surv ->Write(
"numu_surv" );
286 noextrap_h_fd_mc_numu_app ->Write(
"numu_app" );
287 noextrap_h_fd_mc_anue_app ->Write(
"anue_app" );
288 noextrap_h_fd_mc_anue_surv ->Write(
"anue_surv" );
289 noextrap_h_fd_mc_anumu_surv ->Write(
"anumu_surv" );
290 noextrap_h_fd_mc_anumu_app ->Write(
"anumu_app" );
291 noextrap_h_fd_mc_tau_cc_all ->Write(
"tau_cc_all" );
292 noextrap_h_fd_mc_nc ->Write(
"nc" );
293 noextrap_h_fd_mc_total ->Write(
"total" );
300 pots.Write(
"pot_nd_fd");
304 bests[0] = bestcalc->GetdCP();
305 bests[1] = bestcalc->GetTh13();
306 bests[2] = bestcalc->GetTh23();
307 bests[3] = bestcalc->GetDmsq32();
308 bests.Write(
"best_dcp_th13_th23_dmsq32");
312 bestchi[0] = minchi23;
313 bestchi.Write(
"best_chisq");
const Dmsq32Constraint kDmsq32ConstraintPDG2015(2.44e-3, 0.06e-3, 2.49e-3, 0.06e-3)
const XML_Char XML_Encoding * info
TCut intime("tave>=217.0 && tave<=227.8")
virtual _IOscCalcAdjustable< T > * Copy() const =0
void ResetOscCalcToDefault(osc::IOscCalcAdjustable *calc)
Reset calculator to default assumptions for all parameters.
fvar< T > fabs(const fvar< T > &x)
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.
const FitDmSq32 kFitDmSq32
double Integral(double exposure, double *err=0, EExposureType expotype=kPOT) const
Return total number of events scaled to pot.
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
std::unique_ptr< T > LoadFrom(TDirectory *dir, const std::string &label)
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
Charged-current interactions.
void SetValue(osc::IOscCalcAdjustable *osc, double val) const override
virtual T GetDmsq32() const
Very simple model allowing inclusion of reactor constraints.
Combine multiple component experiments.
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
Base class defining interface for experiments.
TCut outtime("(tave>0.0&&tave<200.0) || (tave>250.0&&tave<450.0)")
const FitSinSq2Theta13 kFitSinSq2Theta13
All neutrinos, any flavor.
const FitVarWithPrior fitDeltaInPiUnits_UniformPriordCP & kFitDeltaInPiUnits
double Livetime() const
Seconds. For informational purposes only. No calculations use this.
virtual void SetdCP(const T &dCP)=0
Compare a single data spectrum to the MC + cosmics expectation.
Perform MINUIT fits in one or two dimensions.