46 #include "CAFAna/Core/Binning.h" 52 #include "CAFAna/Core/Var.h" 67 #include "TGraphAsymmErrors.h" 73 #include "G4RwgtPlots/Geant4WeightVars.h" 82 #include "/nova/app/users/csweeney/xsec-TR/Utilities/misc.h" 99 TH1 *
ret = (TH1*) h_up->Clone();
100 for(
auto ibinx = 1; ibinx <= ret->GetNbinsX(); ibinx++) {
101 for(
auto ibiny = 1; ibiny <= ret->GetNbinsY(); ibiny++) {
102 for(
auto ibinz = 1; ibinz <= ret->GetNbinsZ(); ibinz++) {
103 double up = h_up ->GetBinContent(ibinx, ibiny, ibinz);
104 double down = h_down ->GetBinContent(ibinx, ibiny, ibinz);
105 double nominal = h_nominal ->GetBinContent(ibinx, ibiny, ibinz);
107 ret->SetBinContent(ibinx, ibiny, ibinz, abs_shift + nominal);
121 TDirectory *
tmp = gDirectory;
122 TCanvas *
c =
new TCanvas();
123 c->SetLeftMargin(0.15);
124 c->SetRightMargin(0.01);
126 double max_val = -1e9;
127 for(
auto ihist = 0
u; ihist < universes.size(); ihist++) {
130 h_nominal->GetYaxis()->SetRangeUser(0, max_val * 1.2);
131 h_nominal->SetTitle(title.c_str());
132 h_nominal->Draw(
"hist");
133 for(
auto ihist = 0
u; ihist < universes.size(); ihist++) {
134 universes[ihist]->Draw(
"hist same");
136 h_nominal->Draw(
"hist same");
138 c->Print(name.c_str());
189 std::string inName =
"/pnfs/nova/scratch/users/csweeney/final_grid_prod5_g4rwgt_miniprod/*";
197 const Var * common_weight =
new Var(kUnity);
202 std::vector<std::vector<Var>> weights_vec{
203 GetkGeantPiplusWeights(), GetkGeantPiplusWeights(), GetkGeantPiplusWeights(),
204 GetkGeantPiminusWeights(), GetkGeantPiminusWeights(), GetkGeantPiminusWeights(),
205 GetkGeantAllPionsWeights(), GetkGeantAllPionsWeights(), GetkGeantAllPionsWeights()
207 std::vector<std::string> myLabel_vec{
208 "Vertex in X Direction (cm)",
"Vertex in Y Direction (cm)",
"Vertex in Z Direction (cm)",
209 "Vertex in X Direction (cm)",
"Vertex in Y Direction (cm)",
"Vertex in Z Direction (cm)",
210 "Vertex in X Direction (cm)",
"Vertex in Y Direction (cm)",
"Vertex in Z Direction (cm)" 212 std::vector<std::string> myTitle_vec{
213 "Piplus 100 universes",
"Piplus 100 universes",
"Piplus 100 universes",
214 "Piminus 100 universes",
"Piminus 100 universes",
"Piminus 100 universes",
215 "Total 100 universes",
"Total 100 universes",
"Total 100 universes" 217 std::vector<std::string> myName_vec = {
218 "eff_piplus_vtxX",
"eff_piplus_vtxY",
"eff_piplus_vtxZ",
219 "eff_piminus_vtxX",
"eff_piminus_vtxY",
"eff_piminus_vtxZ",
220 "eff_total_vtxX",
"eff_total_vtxY",
"eff_total_vtxZ" 222 std::vector<Binning> myBins_vec{
227 std::vector<const Var *> kMyVar_vec{
238 std::vector<Multiverse *> myMV_vec;
239 std::vector<Spectrum *> mySpec_vec;
243 for(
uint i=0;
i<myName_vec.size(); ++
i){
247 const Var * kMyVar = kMyVar_vec[
i];
248 std::vector<Var>
weights = weights_vec[
i];
270 myMV_vec.push_back(myMV);
271 mySpec_vec.push_back(sSig);
282 for(
uint i=0;
i<myMV_vec.size(); ++
i){
306 TH1 * hSig = sSig->
ToTH1(pot);
315 TH1 * hUp = sUp->
ToTH1(pot);
316 std::vector<TH1*> hUp_vec{hUp};
320 TH1 * hDown = sDown->
ToTH1(pot);
321 std::vector<TH1*> hDown_vec{hDown};
325 int nBins = hSig->GetNbinsX();
326 double xLow = hSig->GetXaxis()->GetXmin();
327 double xHigh = hSig->GetXaxis()->GetXmax();
329 TH1D * fracUncert =
new TH1D(
"frac",
"", nBins, xLow, xHigh);
331 double up = hUp->GetBinContent(
i);
332 double down = hDown->GetBinContent(
i);
333 double nom = hSig->GetBinContent(
i);
336 double frac = abs_shift / nom;
339 fracUncert->SetBinContent(
i, frac);
342 fracUncert->GetXaxis()->SetTitle(myLabel.c_str());
351 graph->SetTitle(myTitle.c_str());
355 hSig->SetLineColor(1);
356 hSig->SetLineWidth(2);
365 TString rCname =
"rC";
366 TCanvas* rC =
new TCanvas(rCname, rCname);
368 rC -> SetBottomMargin(0.);
370 TPad* P1 =
new TPad(
"Temp_1",
"", 0.0, Spl, 1.0, 1.0, 0 );
371 TPad* P2 =
new TPad(
"Temp_2",
"", 0.0, 0.0, 1.0, Spl, 0 );
372 P2 -> SetRightMargin (.03);
373 P2 -> SetTopMargin (.00);
374 P2 -> SetBottomMargin(.3);
375 P2 -> SetLeftMargin (.1);
377 P1 -> SetRightMargin (.03);
378 P1 -> SetLeftMargin (.1);
379 P1 -> SetTopMargin (.10);
380 P1 -> SetBottomMargin(.00);
389 fracUncert->GetYaxis()->SetTitleSize( Lb2 );
390 fracUncert->GetYaxis()->SetTitleOffset(0.4);
391 fracUncert->GetYaxis()->SetLabelSize( Lb2 );
392 fracUncert->GetXaxis()->SetTitleSize( Lb2 );
393 fracUncert->GetXaxis()->SetLabelSize( Lb2 );
394 fracUncert->Draw(
"");
399 hSig->SetTitle(myTitle.c_str());
400 hSig->GetYaxis()->SetTitleSize( Lb1 );
401 hSig->GetYaxis()->SetLabelSize( Lb1 );
402 hSig->GetYaxis()->SetTitleOffset( 0.7 );
404 hSig->GetXaxis()->SetLabelSize (0 );
405 hSig->GetXaxis()->SetLabelOffset(99);
406 hSig->SetLineColor(4);
407 hSig->DrawCopy(
"hist");
409 hSig->Draw(
"hist ][");
410 graph->Draw(
"e2 same");
411 hSig->Draw(
"hist ][ same");
413 rC->SaveAs((
"./plots/multiverse/" + DateTime() + myName +
".png").c_str());
T max(const caf::Proxy< T > &a, T b)
_HistAxis< Var > HistAxis
Represent the binning of a Spectrum's x-axis.
Cuts and Vars for the 2020 FD DiF Study.
Spectrum * GetPlusOneSigmaShift(const Spectrum *)
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.
Simple record of shifts applied to systematic parameters.
TGraphAsymmErrors * PlotWithSystErrorBand(IPrediction *pred, const std::vector< const ISyst * > &systs, osc::IOscCalc *calc, double pot, int col, int errCol, float headroom, bool newaxis, EBinType bintype, double alpha)
Plot prediction with +/-1sigma error band.
Proxy for caf::StandardRecord.
Spectrum * GetMinusOneSigmaShift(const Spectrum *)
_Var< caf::SRProxy > Var
Representation of a variable to be retrieved from a caf::StandardRecord object.
Representation of a spectrum in any variable, with associated POT.
const NuTruthCut kTrueDetector_NT([](const caf::SRNeutrinoProxy *nu){bool fid=(nu->vtx.X()< 192 && nu->vtx.X() >-191 && nu->vtx.Y()< 194 && nu->vtx.Y() >-187 && nu->vtx.Z()< 1270 && nu->vtx.Z() > 0);return fid; })
_Cut< caf::SRProxy > Cut
Representation of a cut (selection) to be applied to a caf::StandardRecord object.
virtual void Go() override
Load all the registered spectra.
void PlotDebug(std::vector< TH1 * > universes, TH1 *h_nominal, std::string title, std::string name)
double frac(double x)
Fractional part.
UpDownPair< TH1 > ToUpDownHist(Multiverse *mv, const TH1 *h_nominal)
std::vector< float > Spectrum
const NuTruthCut kIsNumuCC1Pi_NT([](const caf::SRNeutrinoProxy *nu){int nPions=0;int nMuons=0;int nPrims=nu->prim.size();for(int prim_idx=0;prim_idx< nPrims;prim_idx++){auto &prim=nu->prim[prim_idx];int pdg=prim.pdg;if(abs(pdg)==211)++nPions;else if(pdg==13)++nMuons;}bool is_numucc=nu->iscc &&nu->pdg==14;bool ret=is_numucc &&(nPions==1)&&(nMuons==1);return ret;})
const SystShifts kNoShift
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
Var kUnity([](const caf::SRProxy *sr){return 1.f;})
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
Cut CutFromNuTruthCut(const NuTruthCut &stc)