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

inputs	input ROOT files (n-tuple)
output	output ROOT file (n-tuple)
config	config handled with `Darwin::Tools::Options`
steering	parameters obtained from explicit options
slice	number and index of slice

                               {1,0} 
        )
{
    cout << __func__ << ' ' << slice << " start" << endl;
 
    DT::Flow flow(steering);
    auto tIn = flow.GetInputTree(inputs, slice);
    auto [fOut, tOut] = flow.GetOutput(output);
 
    DT::MetaInfo metainfo(tOut);
    metainfo.Check(config);
 
    auto isMC = metainfo.Get<bool>("flags", "isMC");
    if (isMC)
        BOOST_THROW_EXCEPTION( DE::BadInput("Only real data may be used as input.", metainfo) );
 
    auto year = metainfo.Get<int>("flags", "year");
    if (year != 2016 && year != 2017)
        BOOST_THROW_EXCEPTION( DE::BadInput("Only 2016 and 2017 data may be used as input.", metainfo) );
 
    auto genEvt = flow.GetBranchReadOnly<GenEvent>("genEvent");
    auto recEvt = flow.GetBranchReadWrite<RecEvent>("recEvent");
    auto recJets = flow.GetBranchReadOnly<vector<RecJet>>("recJets");
 
    ControlPlots raw("raw");
    vector<ControlPlots> calib { ControlPlots("nominal") };
    if (steering & DT::syst) {
        metainfo.Set<string>("variations", RecEvent::WeightVar, "Pref" + SysDown);
        metainfo.Set<string>("variations", RecEvent::WeightVar, "Pref" + SysUp);
 
        calib.push_back(ControlPlots("Pref" + SysDown));
        calib.push_back(ControlPlots("Pref" + SysUp));
    }
 
    auto file = config.get<fs::path>("corrections.prefiring.file");
    auto name = config.get<string>  ("corrections.prefiring.name");
    auto type = config.get<string>  ("corrections.prefiring.type");
    Prefiring::Functor apply(file, name, type, steering & DT::syst);
    metainfo.Set<fs::path>("corrections", "prefiring", "file", file);
    metainfo.Set<string>  ("corrections", "prefiring", "name", name);
    metainfo.Set<string>  ("corrections", "prefiring", "type", type);
 
    auto totRecWgt = [&](size_t i) {
        return (isMC ? genEvt->weights.front().v : 1.) * recEvt->weights.at(i).v;
    };
 
    for (DT::Looper looper(tIn); looper(); ++looper) {
        [[ maybe_unused ]]
        static auto& cout = steering & DT::verbose ? ::cout : DT::dev_null;
 
        raw(*recJets, totRecWgt(0));
        apply(*recEvt, *recJets);
        for (size_t i = 0; i < calib.size(); ++i)
            calib.at(i)(*recJets, totRecWgt(i));
 
        if (steering & DT::fill) tOut->Fill();
    }
 
    raw.Write(fOut);
    for (size_t i = 0; i < calib.size(); ++i)
        calib.at(i).Write(fOut);
 
    TDirectory * controlplots = fOut->mkdir("Prefiring");
    apply.Write(controlplots);
 
    metainfo.Set<bool>("git", "complete", true);
 
    cout << __func__ << ' ' << slice << " stop" << endl;
}

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

input	input ROOT files (histograms)
output	output ROOT file (histograms)
steering	parameters obtained from explicit options

{
    cout << __func__ << " start" << endl;
 
    auto fIn  = make_unique<TFile>(input.c_str(), "READ");
    auto fOut = make_unique<TFile>(output.c_str(), "RECREATE"); // TODO
 
    const vector<TString> eras {"16BCD", "16EF", "16FGH",
                                "17B", "17C", "17D", "17E", "17F"};
    auto N = eras.size();
 
    auto hA      = make_unique<TH1D>("a"     , ";;a"     , N, 0.5, 0.5 + N);
    auto hMu     = make_unique<TH1D>("mu"    , ";;#mu"   , N, 0.5, 0.5 + N);
    auto hSigma  = make_unique<TH1D>("sigma" , ";;#sigma", N, 0.5, 0.5 + N);
    auto hA2     = make_unique<TH1D>("a2"    , ";;a"     , N, 0.5, 0.5 + N);
    auto hMu2    = make_unique<TH1D>("mu2"   , ";;#mu"   , N, 0.5, 0.5 + N);
    auto hSigma2 = make_unique<TH1D>("sigma2", ";;#sigma", N, 0.5, 0.5 + N);
    for (TH1 * h: { hA .get(), hMu .get(), hSigma .get(),
                    hA2.get(), hMu2.get(), hSigma2.get() })
        h->SetDirectory(fOut.get());
 
    auto maxEta = 2.0;
    auto maxPt = 6500/cosh(maxEta);
    auto maxPt2 = 6500/cosh(maxEta+0.5);
 
    auto aMin = 0.0, aMax = 1.3;
    auto muMin = 100., muMax = 200.;
    auto sigmaMin = 0.1, sigmaMax = 2.;
 
    hA->SetMinimum(aMin);
    hA->SetMaximum(aMax);
    hMu->SetMinimum(muMin);
    hMu->SetMaximum(muMax);
    hSigma->SetMinimum(sigmaMin);
    hSigma->SetMaximum(sigmaMax);
    for (size_t i = 0; i < N; ++i) {
        TString era = eras.at(i);
 
        hA    ->GetXaxis()->SetBinLabel(i+1,era);
        hMu   ->GetXaxis()->SetBinLabel(i+1,era);
        hSigma->GetXaxis()->SetBinLabel(i+1,era);
 
        TString hname = "L1prefiring_jetpt_20" + era;
        auto h2 = unique_ptr<TH2>(fIn->Get<TH2>(hname));
        if (h2.get() == nullptr) {
            TString what = "Map " + hname + " can't be found";
            BOOST_THROW_EXCEPTION( DE::BadInput(what.Data(), h2) );
        }
        h2->SetDirectory(0);
        int yn = h2->GetXaxis()->FindBin(-maxEta-0.1),
            yp = h2->GetXaxis()->FindBin( maxEta+0.1),
            yn2 = h2->GetXaxis()->FindBin(-maxEta-0.6),
            yp2 = h2->GetXaxis()->FindBin( maxEta+0.6);
        auto h = unique_ptr<TH1>(h2->ProjectionY(Form("%d", rand()), yn, yn));
        h->Add(                  h2->ProjectionY(Form("%d", rand()), yp, yp));
        h->Scale(0.5);
        auto k = unique_ptr<TH1>(h2->ProjectionY(Form("%d", rand()), yn2, yn2));
        k->Add(                  h2->ProjectionY(Form("%d", rand()), yp2, yp2));
        k->Scale(0.5);
 
        auto g = Prefiring::Copy(h, maxPt);
        auto g2 = Prefiring::Copy(k, maxPt2);
        g->SetLineColor(kBlue);
        g->SetMarkerColor(kBlue);
        g2->SetLineColor(kRed);
        g2->SetMarkerColor(kRed);
 
        auto f  = new TF1("L1prefiring_jetpt_20" + era + "_ybin5", Prefiring::ansatz, 30, maxPt, 3);
        auto f2 = new TF1("L1prefiring_jetpt_20" + era + "_ybin6", Prefiring::ansatz, 30, maxPt2, 3);
        f->SetTitle("2.0 < |y| < 2.5");
        f->SetLineColor(g->GetLineColor());
        f->SetLineWidth(1);
        f2->SetTitle("2.5 < |y| < 3.0");
        f2->SetLineColor(g2->GetLineColor());
        f2->SetLineWidth(1);
 
        f->FixParameter(0,0.5*(aMin+aMax));           f2->FixParameter(0,0.5*(aMin+aMax));
        f->SetParLimits(0,aMin,aMax);                 f2->SetParLimits(0,aMin,aMax);
 
        f->SetParameter(1,0.5*(muMin+muMax));         f2->SetParameter(1,0.5*(muMin+muMax));
        f->SetParLimits(1,muMin,muMax);               f2->SetParLimits(1,muMin,muMax);
 
        f->SetParameter(2,0.5*(sigmaMin+sigmaMax));   f2->SetParameter(2,0.5*(sigmaMin+sigmaMax));
        f->SetParLimits(2,sigmaMin,sigmaMax);         f2->SetParLimits(2,sigmaMin,sigmaMax);
 
        auto r = g->Fit(f, "RS");
        auto r2 = g2->Fit(f2, "RS");
 
        fOut->cd();
        f ->Write();
        f2->Write();
 
        auto a     = r->Parameter(0), Ea     = r->Error(0),
             mu    = r->Parameter(1), Emu    = r->Error(1),
             sigma = r->Parameter(2), Esigma = r->Error(2);
 
        auto fUp = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt, 3);
        auto fDn = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt, 3);
        fUp->SetLineStyle(kDotted);
        fUp->SetLineColor(kBlue);
        fUp->SetLineWidth(1);
        fUp->SetParameter(0,a*1.2);
        fUp->SetParameter(1,mu);
        fUp->SetParameter(2,sigma);
        fDn->SetLineStyle(kDotted);
        fDn->SetLineColor(kBlue);
        fDn->SetLineWidth(1);
        fDn->SetParameter(0,a*0.8);
        fDn->SetParameter(1,mu);
        fDn->SetParameter(2,sigma);
 
        hA    ->SetBinContent(i+1,a);         hA    ->SetBinError(i+1,Ea);
        hMu   ->SetBinContent(i+1,mu);        hMu   ->SetBinError(i+1,Emu);
        hSigma->SetBinContent(i+1,sigma);     hSigma->SetBinError(i+1,Esigma);
 
        a     = r2->Parameter(0); Ea     = r2->Error(0);
        mu    = r2->Parameter(1); Emu    = r2->Error(1);
        sigma = r2->Parameter(2); Esigma = r2->Error(2);
 
        auto f2Up = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt2, 3);
        auto f2Dn = make_unique<TF1>(Form("%d", rand()), Prefiring::ansatz, 30, maxPt2, 3);
        f2Up->SetLineStyle(kDotted);
        f2Up->SetLineColor(kRed);
        f2Up->SetLineWidth(1);
        f2Up->SetParameter(0,a*1.2);
        f2Up->SetParameter(1,mu);
        f2Up->SetParameter(2,sigma);
        f2Dn->SetLineStyle(kDotted);
        f2Dn->SetLineColor(kRed);
        f2Dn->SetLineWidth(1);
        f2Dn->SetParameter(0,a*0.8);
        f2Dn->SetParameter(1,mu);
        f2Dn->SetParameter(2,sigma);
 
        hA2    ->SetBinContent(i+1,a);         hA2    ->SetBinError(i+1,Ea);
        hMu2   ->SetBinContent(i+1,mu);        hMu2   ->SetBinError(i+1,Emu);
        hSigma2->SetBinContent(i+1,sigma);     hSigma2->SetBinError(i+1,Esigma);
    }
 
    cout << __func__ << " stop" << endl;
}