drawVsPOT.C
Go to the documentation of this file.
1 #ifdef __CINT__
2 void drawVsPOT(std::string samdef, unsigned int timestamp, unsigned int nDays, std::string naming){
3  std::cout << "Sorry, you must run in compiled mode" << std::endl;
4 }
5 #else
6 
7 #include "CAFAna/Cuts/Cuts.h"
10 #include "CAFAna/Cuts/NumuCuts2018.h"
11 #include "CAFAna/Core/Loaders.h"
13 #include "CAFAna/Core/Spectrum.h"
15 #include "CAFAna/Vars/Vars.h"
16 #include "CAFAna/Analysis/Calcs.h"
17 #include "CAFAna/Analysis/Plots.h"
18 
19 #include "TStyle.h"
20 #include "TCanvas.h"
21 #include "TGraph.h"
22 #include "TGraphErrors.h"
23 #include "TGraphAsymmErrors.h"
24 #include "TFile.h"
25 #include "TH1.h"
26 #include "TH2.h"
27 #include "TLegend.h"
28 #include "TLatex.h"
29 #include "TProfile.h"
30 #include "TTimeStamp.h"
31 #include "TPaveLabel.h"
32 
33 #include "VarsAndCuts.h"
34 #include "Periods.h"
35 
36 using namespace ana;
37 
38 void drawVsPOT(std::string samdef, unsigned int timestamp, unsigned int nDays, std::string naming)
39 {
40  std::cout << "Drawing plots..." << std::endl;
41 
42  gStyle->SetOptFit(1111);
43 
44  const std::string outDir = "/nova/ana/users/novadq/BeamMonitoring/Plots/stability/";
45  const std::string fFiles = "dataset_def_name_newest_snapshot " + samdef + Form(" and online.subrunendtime >= '%u' and online.subrunendtime <= '%u'", timestamp - 86400*nDays, timestamp);
46 
47  const Binning kSingleBin = Binning::Simple(1,0.5,1.5);
48 
49  const Cut kNumuNDCut = kNumuCutND2018;
50 
51  SpectrumLoader loader(fFiles);
53 
54  TH1D* hPOT[kNumPeriods];
55  TH1D* hNumPOT[kNumPeriods];
56 
57  for(int j = 0; j < kNumPeriods; j++)
58  {
59  hPOT[j] = new TH1D(Form("pot_%i",j+1),";POT (10^{15});POT (10^{15})",1000,0,100);
60  hNumPOT[j] = new TH1D(Form("npot_%i",j+1),";POT (10^{15});POT (10^{15})",1,0.5,1.5);
61  loader.AddSpillHistogram(hPOT[j], kPOTe12, getYearAndMonth(datasets[j].year, datasets[j].month));
62  loader.AddSpillHistogram(hNumPOT[j], kSpillUnweighted, getYearAndMonth(datasets[j].year, datasets[j].month), kPOTe12);
63  }
64 
65  const int kNumPlots = 1;
66 
67  Plot plots[kNumPlots] = {
68  {";Central Time;Number of slices / 10^{15} POT", "slices", kSingleBin, kUnweighted, kNumuNDCut, 0, 2.0}
69  };
70 
71  /// ====================================================================== ///
72 
74 
75  for(int i = 0; i < kNumPlots; ++i){
76  for(int j = 0; j < kNumPeriods; ++j){
77  hData[i][j] = new Spectrum(plots[i].label, loader, kSingleBin, kUnweighted, plots[i].bins, plots[i].var, kYear == datasets[j].year && kMonth == datasets[j].month && plots[i].cut);
78  } // end loop over j
79  } // end loop over i
80 
81  // GO!
82  loader.Go();
83 
84  for(int i = 0; i < kNumPlots; ++i){
85 
86  std::string fullName = plots[i].fname+"_"+naming;
87  std::string fullNameTotal = plots[i].fname+"_Total_"+naming;
88  TCanvas* c = new TCanvas("c",(fullName).c_str());
89  c->SetBottomMargin(0.15);
90  c->SetLeftMargin(0.15);
91  c->SetRightMargin(0.15);
92 
93  TLegend* leg = new TLegend(0.17,0.17,0.83,0.4);
94  leg->SetFillColor(0);
95  leg->SetNColumns(5);
96 
97  TGraphErrors* grTotal = new TGraphErrors;
98 
99  for(int j = 0; j < kNumPeriods; ++j){
100 
101  std::string monthName = plots[i].fname+Form("_period%i_",j+1);
102 
103  TH2* Histogram = hData[i][j]->ToTH2(hData[i][j]->POT());
104 
105  TProfile* Profile = Histogram->ProfileX((monthName+"Profile").c_str());
106  Profile->SetMarkerStyle(datasets[j].style);
107  Profile->SetMarkerColor(datasets[j].color);
108  Profile->SetLineColor(datasets[j].color);
109  Profile->SetMarkerSize(1.4);
110 
111  Profile->GetXaxis()->CenterTitle();
112  Profile->GetYaxis()->CenterTitle();
113  Profile->GetYaxis()->SetTitle(Histogram->GetYaxis()->GetTitle());
114 
115  TGraphErrors* gr = new TGraphErrors;
116  gr->SetName((monthName+"_1D").c_str());
117 
118  double y = Histogram->Integral();
119  double dy = sqrt(y);
120  double x = hPOT[j]->GetMean();
121  double npot = hNumPOT[j]->Integral()/1000.;
122 
123  if(x > 0)
124  {
125  std::cout << x << " " << y/npot << " " << dy/npot << std::endl;
126  gr->SetPoint(0, x, y/npot);
127  gr->SetPointError(0, 0, dy/npot);
128 
129  grTotal->SetPoint(grTotal->GetN(), x, y/npot);
130  grTotal->SetPointError(grTotal->GetN()-1, 0, dy/npot);
131  }
132 
133  gr->SetTitle("");
134  gr->SetMarkerStyle(datasets[j].style);
135  gr->SetMarkerColor(datasets[j].color);
136  gr->SetLineColor(datasets[j].color);
137 
138  if(j == 0)
139  {
140  TH2D* Hist0 = new TH2D("Hist0","",100, 20, 70, 100, plots[i].yRangeLow, plots[i].yRangeHigh);
141  Hist0->GetXaxis()->CenterTitle();
142  Hist0->GetYaxis()->CenterTitle();
143  Hist0->GetYaxis()->SetTitle("Slices / 10^{9} POT");
144  Hist0->GetXaxis()->SetTitle("Average POT / 10^{12}");
145  Hist0->SetLineColor(0);
146  Hist0->SetLineWidth(0);
147  Hist0->Draw("hist");
148 
149  gr->Draw("P");
150  }
151 
152  else gr->Draw("P");
153 
154  leg->AddEntry(Profile, datasets[j].label.c_str(), "p");
155 
156  } // End loop over j
157 
158  leg->Draw();
159 
160  //c->SaveAs(("stability/"+fullName+".pdf").c_str());
161  //c->SaveAs(("stability/"+fullName+".png").c_str());
162  //c->SaveAs(("stability/"+fullName+".root").c_str());
163  c->SaveAs((outDir+fullName+".pdf").c_str());
164  c->SaveAs((outDir+fullName+".png").c_str());
165  c->SaveAs((outDir+fullName+".root").c_str());
166 
167  c->Clear();
168 
169  grTotal->SetTitle("");
170  grTotal->SetMarkerStyle(20);
171  grTotal->SetMarkerColor(kBlack);
172  grTotal->SetLineColor(kBlack);
173  grTotal->GetXaxis()->CenterTitle();
174  grTotal->GetYaxis()->CenterTitle();
175  grTotal->GetYaxis()->SetTitle("Slices / 10^{9} POT");
176  grTotal->GetXaxis()->SetTitle("Average POT / 10^{12}");
177 
178  grTotal->Draw("AP");
179  grTotal->Fit("pol1");
180 
181  //c->SaveAs(("stability/"+fullNameTotal+".pdf").c_str());
182  //c->SaveAs(("stability/"+fullNameTotal+".png").c_str());
183  //c->SaveAs(("stability/"+fullNameTotal+".root").c_str());
184  c->SaveAs((outDir+fullNameTotal+".pdf").c_str());
185  c->SaveAs((outDir+fullNameTotal+".png").c_str());
186  c->SaveAs((outDir+fullNameTotal+".root").c_str());
187 
188 
189  } // End loop over i
190 
191  std::cout << "All plots finished!" << std::endl;
192 
193 } // End of function
194 
195 #endif
const SpillVar kPOTe12([](const caf::SRSpillProxy *spill){return spill->spillpot/1e12;})
TH2 * ToTH2(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 2D to obtain TH2.
Definition: Spectrum.cxx:165
Represent the binning of a Spectrum&#39;s x-axis.
Definition: Binning.h:16
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Cut kNumuCutND2018
Definition: NumuCuts2018.h:41
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
void SetSpillCut(const SpillCut &cut)
std::string outDir
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
const char * label
const SpillVar kSpillUnweighted
Definition: Var.h:98
double dy[NP][NC]
====================================================================== ///
Definition: CutFlow_Data.C:28
virtual void AddSpillHistogram(TH1 *h, const SpillVar &var, const SpillCut &cut, const SpillVar &wei=kSpillUnweighted)
Uses include counting the total POT or spills in a run.
const int kNumPeriods
Definition: Periods.h:8
const int kNumPlots
Definition: GetSpectra.h:1
const std::vector< Plot > plots
SpillCut getYearAndMonth(unsigned int iYear, unsigned int iMonth)
Definition: VarsAndCuts.h:71
virtual void Go() override
Load all the registered spectra.
const double j
Definition: BetheBloch.cxx:29
loader
Definition: demo0.py:10
std::vector< float > Spectrum
Definition: Constants.h:610
OStream cout
Definition: OStream.cxx:6
std::vector< double > POT
const Binning bins
Definition: NumuCC_CPiBin.h:8
std::string fname
Definition: CutFlow_Data.C:30
const Cut cut
Definition: exporter_fd.C:30
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
void drawVsPOT(std::string samdef, unsigned int timestamp, unsigned int nDays, std::string naming)
Definition: drawVsPOT.C:38
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
const Var kYear
Definition: VarsAndCuts.h:34
std::vector< Dataset > datasets
Definition: Periods.h:3
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
Definition: Binning.cxx:107
const Var kUnweighted
The simplest possible Var, always 1. Used as a default weight.
Definition: Var.h:96
const Var kMonth
Definition: VarsAndCuts.h:35
static constexpr Double_t year
Definition: Munits.h:185
enum BeamMode string