DAS  3.0
Das Analysis System
Namespaces | Classes | Functions
DAS::Prefiring Namespace Reference

Namespaces

 Prefiring
 

Classes

struct  Functor
 

Functions

void applyPrefiringWeights (const vector< fs::path > &inputs, const fs::path &output, const pt::ptree &config, const int steering, const DT::Slice slice={1, 0})
 
void fitPrefiringCorrections (const fs::path &input, const fs::path &output, const int steering)
 

Function Documentation

◆ applyPrefiringWeights()

void DAS::Prefiring::applyPrefiringWeights ( const vector< fs::path > &  inputs,
const fs::path &  output,
const pt::ptree &  config,
const int  steering,
const DT::Slice  slice = {1,0} 
)

Application of prefiring weights in data or simulation.

Calls a functor

Parameters
inputsinput ROOT files (n-tuple)
outputoutput ROOT file (n-tuple)
configconfig handled with `Darwin::Tools::Options`
steeringparameters obtained from explicit options
slicenumber and index of slice
38  {1,0}
39  )
40 {
41  cout << __func__ << ' ' << slice << " start" << endl;
42 
43  unique_ptr<TChain> tIn = DT::GetChain(inputs);
44  unique_ptr<TFile> fOut(DT_GetOutput(output));
45  auto tOut = unique_ptr<TTree>(tIn->CloneTree(0));
46 
47  DT::MetaInfo metainfo(tOut);
48  metainfo.Check(config);
49 
50  auto isMC = metainfo.Get<bool>("flags", "isMC");
51  if (isMC)
52  BOOST_THROW_EXCEPTION( DE::BadInput("Only real data may be used as input.", metainfo) );
53 
54  auto year = metainfo.Get<int>("flags", "year");
55  if (year != 2016 && year != 2017)
56  BOOST_THROW_EXCEPTION( DE::BadInput("Only 2016 and 2017 data may be used as input.", metainfo) );
57 
58  GenEvent * genEvt = nullptr;
59  RecEvent * recEvt = nullptr;
60  vector<RecJet> * recJets = nullptr;
61  if (isMC)
62  tIn->SetBranchAddress("genEvent", &genEvt);
63  tIn->SetBranchAddress("recEvent", &recEvt);
64  tIn->SetBranchAddress("recJets", &recJets);
65 
66  ControlPlots raw("raw");
67  vector<ControlPlots> calib { ControlPlots("nominal") };
68  bool applySyst = (steering & DT::syst) == DT::syst;
69  if (applySyst) {
70  metainfo.Set<string>("variations", RecEvent::WeightVar, "Pref" + SysDown);
71  metainfo.Set<string>("variations", RecEvent::WeightVar, "Pref" + SysUp);
72 
73  calib.push_back(ControlPlots("Pref" + SysDown));
74  calib.push_back(ControlPlots("Pref" + SysUp));
75  }
76 
77  auto file = config.get<fs::path>("corrections.prefiring.file");
78  auto name = config.get<string> ("corrections.prefiring.name");
79  auto type = config.get<string> ("corrections.prefiring.type");
80  Prefiring::Functor apply(file, name, type, applySyst);
81  metainfo.Set<fs::path>("corrections", "prefiring", "file", file);
82  metainfo.Set<string> ("corrections", "prefiring", "name", name);
83  metainfo.Set<string> ("corrections", "prefiring", "type", type);
84 
85  auto totRecWgt = [&](size_t i) {
86  return (isMC ? genEvt->weights.front().v : 1.) * recEvt->weights.at(i).v;
87  };
88 
89  for (DT::Looper looper(tIn, slice); looper(); ++looper) {
90  [[ maybe_unused ]]
91  static auto& cout = (steering & DT::verbose) == DT::verbose ? ::cout : DT::dev_null;
92 
93  raw(*recJets, totRecWgt(0));
94  apply(*recEvt, *recJets);
95  for (size_t i = 0; i < calib.size(); ++i)
96  calib.at(i)(*recJets, totRecWgt(i));
97 
98  if ((steering & DT::fill) == DT::fill) tOut->Fill();
99  }
100 
101  raw.Write(fOut.get());
102  for (size_t i = 0; i < calib.size(); ++i)
103  calib.at(i).Write(fOut.get());
104 
105  TDirectory * controlplots = fOut->mkdir("Prefiring");
106  apply.Write(controlplots);
107 
108  metainfo.Set<bool>("git", "complete", true);
109  tOut->Write();
110 
111  cout << __func__ << ' ' << slice << " stop" << endl;
112 }

◆ fitPrefiringCorrections()

void DAS::Prefiring::fitPrefiringCorrections ( const fs::path &  input,
const fs::path &  output,
const int  steering 
)

(i) symmetrize the JME-provided tables between plus and minus sides (ii) fit the turn-on with an error function (iii) use the smooth efficiency to correct data

Todo:
proper propagation of uncertainties
Todo:
smooth fit also for the average map
Todo:
check if files exist
Parameters
inputinput ROOT files (histograms)
outputoutput ROOT file (histograms)
steeringparameters obtained from explicit options
77 {
78  cout << __func__ << " start" << endl;
79 
81  auto fIn = make_unique<TFile>(input.c_str(), "READ");
82  auto fOut = make_unique<TFile>(output.c_str(), "RECREATE");
83 
84  const vector<TString> eras {"16BCD", "16EF", "16FGH",
85  "17B", "17C", "17D", "17E", "17F"};
86  auto N = eras.size();
87 
88  auto hA = make_unique<TH1D>("a" , ";;a" , N, 0.5, 0.5 + N);
89  auto hMu = make_unique<TH1D>("mu" , ";;#mu" , N, 0.5, 0.5 + N);
90  auto hSigma = make_unique<TH1D>("sigma" , ";;#sigma", N, 0.5, 0.5 + N);
91  auto hA2 = make_unique<TH1D>("a2" , ";;a" , N, 0.5, 0.5 + N);
92  auto hMu2 = make_unique<TH1D>("mu2" , ";;#mu" , N, 0.5, 0.5 + N);
93  auto hSigma2 = make_unique<TH1D>("sigma2", ";;#sigma", N, 0.5, 0.5 + N);
94  for (TH1 * h: { hA .get(), hMu .get(), hSigma .get(),
95  hA2.get(), hMu2.get(), hSigma2.get() })
96  h->SetDirectory(fOut.get());
97 
98  auto maxEta = 2.0;
99  auto maxPt = 6500/cosh(maxEta);
100  auto maxPt2 = 6500/cosh(maxEta+0.5);
101 
102  auto aMin = 0.0, aMax = 1.3;
103  auto muMin = 100., muMax = 200.;
104  auto sigmaMin = 0.1, sigmaMax = 2.;
105 
106  hA->SetMinimum(aMin);
107  hA->SetMaximum(aMax);
108  hMu->SetMinimum(muMin);
109  hMu->SetMaximum(muMax);
110  hSigma->SetMinimum(sigmaMin);
111  hSigma->SetMaximum(sigmaMax);
112  for (size_t i = 0; i < N; ++i) {
113  TString era = eras.at(i);
114 
115  hA ->GetXaxis()->SetBinLabel(i+1,era);
116  hMu ->GetXaxis()->SetBinLabel(i+1,era);
117  hSigma->GetXaxis()->SetBinLabel(i+1,era);
118 
119  TString hname = "L1prefiring_jetpt_20" + era;
120  auto h2 = unique_ptr<TH2>(fIn->Get<TH2>(hname));
121  if (h2.get() == nullptr) {
122  TString what = "Map " + hname + " can't be found";
123  BOOST_THROW_EXCEPTION( DE::BadInput(what.Data(), h2) );
124  }
125  h2->SetDirectory(0);
126  int yn = h2->GetXaxis()->FindBin(-maxEta-0.1),
127  yp = h2->GetXaxis()->FindBin( maxEta+0.1),
128  yn2 = h2->GetXaxis()->FindBin(-maxEta-0.6),
129  yp2 = h2->GetXaxis()->FindBin( maxEta+0.6);
130  auto h = unique_ptr<TH1>(h2->ProjectionY(Form("%d", rand()), yn, yn));
131  h->Add( h2->ProjectionY(Form("%d", rand()), yp, yp));
132  h->Scale(0.5);
133  auto k = unique_ptr<TH1>(h2->ProjectionY(Form("%d", rand()), yn2, yn2));
134  k->Add( h2->ProjectionY(Form("%d", rand()), yp2, yp2));
135  k->Scale(0.5);
136 
137  auto g = Prefiring::Copy(h, maxPt);
138  auto g2 = Prefiring::Copy(k, maxPt2);
139  g->SetLineColor(kBlue);
140  g->SetMarkerColor(kBlue);
141  g2->SetLineColor(kRed);
142  g2->SetMarkerColor(kRed);
143 
144  auto f = new TF1("L1prefiring_jetpt_20" + era + "_ybin5", Prefiring::ansatz, 30, maxPt, 3);
145  auto f2 = new TF1("L1prefiring_jetpt_20" + era + "_ybin6", Prefiring::ansatz, 30, maxPt2, 3);
146  f->SetTitle("2.0 < |y| < 2.5");
147  f->SetLineColor(g->GetLineColor());
148  f->SetLineWidth(1);
149  f2->SetTitle("2.5 < |y| < 3.0");
150  f2->SetLineColor(g2->GetLineColor());
151  f2->SetLineWidth(1);
152 
153  f->FixParameter(0,0.5*(aMin+aMax)); f2->FixParameter(0,0.5*(aMin+aMax));
154  f->SetParLimits(0,aMin,aMax); f2->SetParLimits(0,aMin,aMax);
155 
156  f->SetParameter(1,0.5*(muMin+muMax)); f2->SetParameter(1,0.5*(muMin+muMax));
157  f->SetParLimits(1,muMin,muMax); f2->SetParLimits(1,muMin,muMax);
158 
159  f->SetParameter(2,0.5*(sigmaMin+sigmaMax)); f2->SetParameter(2,0.5*(sigmaMin+sigmaMax));
160  f->SetParLimits(2,sigmaMin,sigmaMax); f2->SetParLimits(2,sigmaMin,sigmaMax);
161 
162  auto r = g->Fit(f, "RS");
163  auto r2 = g2->Fit(f2, "RS");
164 
165  fOut->cd();
166  f ->Write();
167  f2->Write();
168 
169  auto a = r->Parameter(0), Ea = r->Error(0),
170  mu = r->Parameter(1), Emu = r->Error(1),
171  sigma = r->Parameter(2), Esigma = r->Error(2);
172 
173  auto fUp = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt, 3);
174  auto fDn = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt, 3);
175  fUp->SetLineStyle(kDotted);
176  fUp->SetLineColor(kBlue);
177  fUp->SetLineWidth(1);
178  fUp->SetParameter(0,a*1.2);
179  fUp->SetParameter(1,mu);
180  fUp->SetParameter(2,sigma);
181  fDn->SetLineStyle(kDotted);
182  fDn->SetLineColor(kBlue);
183  fDn->SetLineWidth(1);
184  fDn->SetParameter(0,a*0.8);
185  fDn->SetParameter(1,mu);
186  fDn->SetParameter(2,sigma);
187 
188  hA ->SetBinContent(i+1,a); hA ->SetBinError(i+1,Ea);
189  hMu ->SetBinContent(i+1,mu); hMu ->SetBinError(i+1,Emu);
190  hSigma->SetBinContent(i+1,sigma); hSigma->SetBinError(i+1,Esigma);
191 
192  a = r2->Parameter(0); Ea = r2->Error(0);
193  mu = r2->Parameter(1); Emu = r2->Error(1);
194  sigma = r2->Parameter(2); Esigma = r2->Error(2);
195 
196  auto f2Up = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt2, 3);
197  auto f2Dn = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt2, 3);
198  f2Up->SetLineStyle(kDotted);
199  f2Up->SetLineColor(kRed);
200  f2Up->SetLineWidth(1);
201  f2Up->SetParameter(0,a*1.2);
202  f2Up->SetParameter(1,mu);
203  f2Up->SetParameter(2,sigma);
204  f2Dn->SetLineStyle(kDotted);
205  f2Dn->SetLineColor(kRed);
206  f2Dn->SetLineWidth(1);
207  f2Dn->SetParameter(0,a*0.8);
208  f2Dn->SetParameter(1,mu);
209  f2Dn->SetParameter(2,sigma);
210 
211  hA2 ->SetBinContent(i+1,a); hA2 ->SetBinError(i+1,Ea);
212  hMu2 ->SetBinContent(i+1,mu); hMu2 ->SetBinError(i+1,Emu);
213  hSigma2->SetBinContent(i+1,sigma); hSigma2->SetBinError(i+1,Esigma);
214  }
215 
216  cout << __func__ << " stop" << endl;
217 }
DYToLL_M-50_13TeV_pythia8_cff_GEN_SIM_RECOBEFMIX_DIGI_L1_DIGI2RAW_L1Reco_RECO.name
name
Definition: DYToLL_M-50_13TeV_pythia8_cff_GEN_SIM_RECOBEFMIX_DIGI_L1_DIGI2RAW_L1Reco_RECO.py:48
Darwin::Tools::syst
@ syst
activate -s to systematic uncertainties
Definition: Options.h:30
metPhiCorrectionExample.eras
list eras
Definition: metPhiCorrectionExample.py:21
Step::verbose
static bool verbose
Definition: Step.h:40
DAS::SysUp
const std::string SysUp
Suffix used for "up" uncertainties. Follows the Combine convention.
Definition: Format.h:8
Ntupliser_cfg.year
int year
Definition: Ntupliser_cfg.py:63
DAS::Prefiring::Prefiring::Copy
unique_ptr< TH1 > Copy(const unique_ptr< TH1 > &h, double maxPt)
Definition: fitPrefiringCorrections.cc:39
Darwin::Tools::Looper
Facility to loop over a n-tuple, including parallelisation and printing.
Definition: Looper.h:33
Darwin::Tools::GetChain
std::unique_ptr< TChain > GetChain(std::vector< std::filesystem::path > inputs, const char *name="events")
Load chain from a list of files.
Definition: FileUtils.cc:67
Ntupliser_cfg.f
f
Definition: Ntupliser_cfg.py:255
generate_html.era
era
Definition: generate_html.py:85
Darwin::Tools::MetaInfo
Generic meta-information for n-tuple (including speficities to Darwin).
Definition: MetaInfo.h:65
DAS::Unfolding::DrellYan::maxEta
static const double maxEta
Definition: ZPtY.h:36
Darwin::Tools::fill
@ fill
activate -f to fill the tree
Definition: Options.h:28
Ntupliser_cfg.config
config
Definition: Ntupliser_cfg.py:263
DAS::Prefiring::Prefiring::ansatz
double ansatz(double *x, double *p)
NOTE: this analytical is partly arbitrary.
Definition: fitPrefiringCorrections.cc:34
DYToLL_M-50_13TeV_pythia8_cff_GEN_SIM_RECOBEFMIX_DIGI_L1_DIGI2RAW_L1Reco_RECO.input
input
Definition: DYToLL_M-50_13TeV_pythia8_cff_GEN_SIM_RECOBEFMIX_DIGI_L1_DIGI2RAW_L1Reco_RECO.py:35
DAS::SysDown
const std::string SysDown
Suffix used for "down" uncertainties. Follows the Combine convention.
Definition: Format.h:10
Ntupliser_cfg.isMC
string isMC
Definition: Ntupliser_cfg.py:59
jercExample.inputs
def inputs
Definition: jercExample.py:118
Darwin::Tools::dev_null
static std::ostream dev_null(nullptr)
to redirect output stream to nowhere
DT_GetOutput
#define DT_GetOutput(output)
Definition: FileUtils.h:96
Darwin::Exceptions::BadInput
Generic exception for ill-defined input (before the event loop).
Definition: exceptions.h:68