164 if (
config.verbose) cout <<
"maxnbaseline = " <<
config.maxnbaseline << endl;
168 int nevents_thiscase;
171 TTree *ntuple =
new TTree(
"FakeBkg",
"Variables from toy MC");
172 ntuple->Branch(
"nevents", &nevents_thiscase,
"nevents/I");
174 ntuple->Branch(
"nevents_sel", &nevents_sel,
"nevents_sel/I");
175 Float_t fake_lep_frac;
176 ntuple->Branch(
"fake_lep_frac", &fake_lep_frac,
"fake_lep_frac/F");
178 ntuple->Branch(
"true_fakes", &true_fakes,
"true_fakes/F");
180 ntuple->Branch(
"asm_fakes", &asm_fakes,
"asm_fakes/F");
182 ntuple->Branch(
"asm_err", &asm_err,
"asm_err/F");
184 ntuple->Branch(
"fkf_fakes", &fkf_fakes,
"fkf_fakes/F");
186 ntuple->Branch(
"fkf_err", &fkf_err,
"fkf_err/F");
187 Float_t lhoodMM_fakes;
188 ntuple->Branch(
"lhoodMM_fakes", &lhoodMM_fakes,
"lhoodMM_fakes/F");
189 Float_t lhoodMM_poserr;
190 ntuple->Branch(
"lhoodMM_poserr", &lhoodMM_poserr,
"lhoodMM_poserr/F");
191 Float_t lhoodMM_negerr;
192 ntuple->Branch(
"lhoodMM_negerr", &lhoodMM_negerr,
"lhoodMM_negerr/F");
193 Float_t lhoodFF_fakes;
194 ntuple->Branch(
"lhoodFF_fakes", &lhoodFF_fakes,
"lhoodFF_fakes/F");
195 Float_t lhoodFF_poserr;
196 ntuple->Branch(
"lhoodFF_poserr", &lhoodFF_poserr,
"lhoodFF_poserr/F");
197 Float_t lhoodFF_negerr;
198 ntuple->Branch(
"lhoodFF_negerr", &lhoodFF_negerr,
"lhoodFF_negerr/F");
201 map<CP::SystematicVariation,float> lhoodMM_weight_map, lhoodMM_poserr_map, lhoodMM_negerr_map;
202 map<CP::SystematicVariation,float> lhoodFF_weight_map, lhoodFF_poserr_map, lhoodFF_negerr_map;
203 map<CP::SystematicVariation,float> asm_weight_map, asm_err_map;
204 map<CP::SystematicVariation,float> fkf_weight_map, fkf_err_map;
206 TRandom3
rand(242868252);
208 Double_t realeff_spread =
config.eff_spread;
209 Double_t fakeeff_spread =
config.eff_spread;
211 std::vector<bool> isTight;
212 std::vector<std::string>
input = {
config.outputdirname+
"efficiencies_full" };
222 input.back().append(
".root");
229 for (
unsigned icase(0); icase <
config.ncases; icase++) {
232 TH1F* h_lep_pt_lhoodMM = 0;
233 TH1F* h_lep_eta_lhoodMM = 0;
234 TH2F* h_lep_pt_eta_lhoodMM = 0;
235 TH1F* h_lep_pt_lhoodFF = 0;
236 TH1F* h_lep_eta_lhoodFF = 0;
237 TH2F* h_lep_pt_eta_lhoodFF = 0;
238 TH1F* h_lep_pt_asm = 0;
239 TH1F* h_lep_eta_asm = 0;
240 TH2F* h_lep_pt_eta_asm = 0;
241 TH1F* h_lep_pt_fkf = 0;
242 TH1F* h_lep_eta_fkf = 0;
243 TH2F* h_lep_pt_eta_fkf = 0;
247 h_lep_pt_lhoodMM =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
249 h_lep_eta_lhoodMM =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
251 h_lep_pt_eta_lhoodMM =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
253 h_lep_pt_lhoodFF =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
255 h_lep_eta_lhoodFF =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
257 h_lep_pt_eta_lhoodFF =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
259 h_lep_pt_asm =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
261 h_lep_eta_asm =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
263 h_lep_pt_eta_asm =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
265 h_lep_pt_fkf =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
267 h_lep_eta_fkf =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
269 h_lep_pt_eta_fkf =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
272 float lep_pt, lep_eta;
276 cout <<
"Starting case " << icase << endl;
279 std::vector<CP::LhoodMM_tools*> lhmTool_sav;
281 std::vector<CP::LhoodMM_tools*> lhmTool_FF_sav;
283 std::vector<CP::AsymptMatrixTool*> asmTool_sav;
285 std::vector<CP::ApplyFakeFactor*> fkfTool_sav;
292 ANA_CHECK( lhmTool.register1DHistogram(h_lep_pt_lhoodMM, &lep_pt) );
293 ANA_CHECK( lhmTool.register1DHistogram(h_lep_eta_lhoodMM, &lep_eta) );
294 ANA_CHECK( lhmTool.register2DHistogram(h_lep_pt_eta_lhoodMM, &lep_pt, &lep_eta) );
295 ANA_CHECK( lhmTool_FF.register1DHistogram(h_lep_pt_lhoodFF, &lep_pt) );
296 ANA_CHECK( lhmTool_FF.register1DHistogram(h_lep_eta_lhoodFF, &lep_eta) );
297 ANA_CHECK( lhmTool_FF.register2DHistogram(h_lep_pt_eta_lhoodFF, &lep_pt, &lep_eta) );
298 ANA_CHECK( fkfTool.register1DHistogram(h_lep_pt_fkf, &lep_pt) );
299 ANA_CHECK( fkfTool.register1DHistogram(h_lep_eta_fkf, &lep_eta) );
300 ANA_CHECK( fkfTool.register2DHistogram(h_lep_pt_eta_fkf, &lep_pt, &lep_eta) );
301 ANA_CHECK( asmTool.register1DHistogram(h_lep_pt_asm, &lep_pt) );
302 ANA_CHECK( asmTool.register1DHistogram(h_lep_eta_asm, &lep_eta) );
303 ANA_CHECK( asmTool.register2DHistogram(h_lep_pt_eta_asm, &lep_pt, &lep_eta) );
309 ANA_CHECK( lhmTool_FF.setProperty(
"DoFakeFactorFit",
true) );
314 for (
int iSave = 0; iSave <=
nSave; iSave++) {
317 lhmTool_sav.push_back(lhmTool_is);
326 lhmTool_FF_sav.push_back(lhmTool_FF_is);
328 ANA_CHECK( lhmTool_FF_sav[iSave]->setProperty(
"DoFakeFactorFit",
true) );
341 asmTool_sav.push_back(asmTool_is);
346 fkfTool_sav.push_back(fkfTool_is);
357 Double_t realeff_mean_thiscase =
rand.Gaus(
config.realeff_mean, realeff_spread);
358 if (realeff_mean_thiscase > 0.99) realeff_mean_thiscase = 0.99;
359 Double_t fakeeff_mean_thiscase =
rand.Gaus(
config.fakeeff_mean, fakeeff_spread);
360 if (fakeeff_mean_thiscase < 0.01) fakeeff_mean_thiscase = 0.01;
361 if (fakeeff_mean_thiscase > realeff_mean_thiscase) fakeeff_mean_thiscase = realeff_mean_thiscase-0.01;
363 fake_lep_frac =
rand.Uniform();
371 std::vector<LhoodMMEvent> mmevts;
375 std::vector<double> realeff, fakeeff;
377 vector<FinalState*>
fs;
379 for (
unsigned ilep(0); ilep <=
config.maxnbaseline; ilep++) {
382 fs.push_back(this_fs);
385 if (
config.poisson_fluctuations) {
386 nevents_thiscase =
rand.Poisson(
config.nevents);
388 nevents_thiscase =
config.nevents;
393 for (
int pass(0); pass<2; pass++) {
397 int nEventsMultFactor = 1;
398 if (pass > 0) nEventsMultFactor = 10;
399 for (Long64_t ievt(0); ievt<nevents_thiscase*nEventsMultFactor; ievt++) {
401 float nlep_frac = 1./(
config.maxnbaseline-
config.minnbaseline+1);
402 Int_t nlep_select =
config.minnbaseline+
rand.Uniform()/nlep_frac;
404 float extraweight(1.0);
410 vector<FakeBkgTools::ParticleData> leptons_data;
412 vector<bool> lep_real;
415 vector<float> lep_pts;
416 for (
int ilep(0); ilep < nlep_select; ilep++) {
417 lep_pts.push_back(100*
rand.Uniform());
420 sort(lep_pts.begin(), lep_pts.end(), std::greater<float>());
422 for (
int ilep(0); ilep < nlep_select; ilep++) {
426 if (
rand.Uniform() > fake_lep_frac) {
432 lep_real.push_back(isReal);
436 lep_pt = lep_pts[ilep];
437 lep_eta = -5. + 10*
rand.Uniform();
438 lepton->
setP4( 1000*lep_pt, lep_eta, 0, 0.511);
445 lepton->
auxdata<
char>(
"Tight") =
true;
447 lepton->
auxdata<
char>(
"Tight") =
false;
451 lepton->
auxdata<
char>(
"Tight") =
true;
453 lepton->
auxdata<
char>(
"Tight") =
false;
458 leptons_data.push_back(lepton_data);
462 ANA_CHECK( asmTool.addEvent(leptons, extraweight) );
463 ANA_CHECK( fkfTool.addEvent(leptons, extraweight) );
464 ANA_CHECK( lhmTool.addEvent(leptons, extraweight) );
465 ANA_CHECK( lhmTool_FF.addEvent(leptons, extraweight) );
467 ANA_CHECK( lhmTool_sav[savIndex]->addEvent(leptons, extraweight) );
468 ANA_CHECK( lhmTool_FF_sav[savIndex]->addEvent(leptons, extraweight) );
469 ANA_CHECK( asmTool_sav[savIndex]->addEvent(leptons, extraweight) );
470 ANA_CHECK( fkfTool_sav[savIndex]->addEvent(leptons, extraweight) );
477 for (
int ilep(0); ilep < nlep_select; ilep++) {
478 if (!lep_real[ilep]) {
485 if (pass == 1 && all_real) {
486 ANA_CHECK( fkfTool.addEvent(leptons, -extraweight/nEventsMultFactor) );
489 ANA_CHECK( lhmTool_FF.addEvent(leptons, -extraweight/nEventsMultFactor) );
492 ANA_CHECK( fkfTool_sav[savIndex]->addEvent(leptons, -extraweight/nEventsMultFactor) );
493 ANA_CHECK( lhmTool_FF_sav[savIndex]->addEvent(leptons, -extraweight/nEventsMultFactor) );
499 vector<int> lep_indices;
500 std::bitset<64> tights, reals, charges;
501 for (
int ilep(0); ilep < nlep_select; ilep++) {
502 lep_indices.push_back(ilep);
503 if (lep_real[ilep]) reals.set(ilep);
507 for (
long comb(0); comb < (1<<nlep_select); ++comb) {
508 tights = std::bitset<64>(comb);
509 if (
fs[nlep_select]->accept_selection(tights, charges)) {
510 if (
fs[nlep_select]->accept_process(nlep_select, reals, tights)) {
511 true_fakes += extraweight*
comboProb(leptons_data, tights, reals);
517 for (
int ilep = 0; ilep < nlep_select; ilep++) {
518 if (leptons[ilep]->auxdata<char>(
"Tight")) tights.set(ilep);
520 if (
fs[nlep_select]->accept_selection(tights,charges) ) {
521 nevents_sel += extraweight;
537 fkfYield -= wgt/nEventsMultFactor;
538 fkfErr += wgt*wgt/(nEventsMultFactor);
544 if (ievt > 0 && (((
nSave*ievt)%nevents_thiscase ==0) || ievt == nevents_thiscase -1) ) {
559 if (ievt > 0 && (((
nSave*ievt)%(nevents_thiscase*nEventsMultFactor) ==0) || ievt == nevents_thiscase*nEventsMultFactor -1) ) {
573 float nfakes_tmp(0), err_tmp(0);
574 ANA_CHECK( fkfTool.getTotalYield(nfakes_tmp, err_tmp, err_tmp) );;
580 if (*
it !=
nullptr)
delete *
it;
585 asm_err = sqrt(asmErr);
586 asm_fakes = asmYield;
587 ANA_CHECK( asmTool.getTotalYield(asmYield, asmErr, asmErr) );
590 ANA_CHECK( lhmTool.getTotalYield(lhoodMM_fakes, lhoodMM_poserr, lhoodMM_negerr) );
592 ANA_CHECK( lhmTool_FF.getTotalYield(lhoodFF_fakes, lhoodFF_poserr, lhoodFF_negerr) );
594 ANA_CHECK( fkfTool.getTotalYield(fkfYield, fkfErr, fkfErr) );
596 fkf_fakes = fkfYield;
598 if (
config.test_systematics) {
599 auto sysvars = asmTool.affectingSystematics();
600 std::string systBrName, systBrNameF;
602 for(
auto& sysvar : sysvars){
603 if (
config.verbose) asmTool.getSystDescriptor().printUncertaintyDescription(sysvar);
605 float asm_syst_weight, asm_syst_err;
606 float fkf_syst_weight, fkf_syst_err;
607 float lhoodMM_syst_weight, lhoodMM_syst_poserr, lhoodMM_syst_negerr;
608 float lhoodFF_syst_weight, lhoodFF_syst_poserr, lhoodFF_syst_negerr;
613 ANA_CHECK(
setupSystBranchesAsym(
"lhoodMM", sysvar, lhoodMM_syst_weight, lhoodMM_syst_poserr, lhoodMM_syst_negerr, lhoodMM_weight_map, lhoodMM_poserr_map, lhoodMM_negerr_map, ntuple) );
614 ANA_CHECK(
setupSystBranchesAsym(
"lhoodFF", sysvar, lhoodFF_syst_weight, lhoodFF_syst_poserr, lhoodFF_syst_negerr, lhoodFF_weight_map, lhoodFF_poserr_map, lhoodFF_negerr_map, ntuple) );
616 ANA_CHECK( asmTool.applySystematicVariation({sysvar}) );
617 ANA_CHECK( asmTool.getTotalYield(asm_weight_map.find(sysvar)->second,
618 asm_err_map.find(sysvar)->second,
619 asm_err_map.find(sysvar)->second));
621 ANA_CHECK( fkfTool.applySystematicVariation({sysvar}) );;
622 ANA_CHECK( fkfTool.getTotalYield(fkf_weight_map.find(sysvar)->second,
623 fkf_err_map.find(sysvar)->second,
624 fkf_err_map.find(sysvar)->second) );
625 ANA_CHECK( lhmTool.applySystematicVariation({sysvar}) );;
626 ANA_CHECK( lhmTool.getTotalYield(lhoodMM_weight_map.find(sysvar)->second,
627 lhoodMM_poserr_map.find(sysvar)->second,
628 lhoodMM_negerr_map.find(sysvar)->second) );
629 ANA_CHECK( lhmTool_FF.applySystematicVariation({sysvar}) );;
630 ANA_CHECK( lhmTool_FF.getTotalYield(lhoodFF_weight_map.find(sysvar)->second,
631 lhoodFF_poserr_map.find(sysvar)->second,
632 lhoodFF_negerr_map.find(sysvar)->second) );
637 cout <<
"OUTPUT true_fakes = " << true_fakes << endl;
638 cout <<
"OUTPUT nfakes for lhoodMM = " << lhoodMM_fakes <<
" +" << lhoodMM_poserr <<
" -" << lhoodMM_negerr << endl;
639 cout <<
"OUTPUT nfakes for lhoodFF = " << lhoodFF_fakes <<
" +" << lhoodFF_poserr <<
" -" << lhoodFF_negerr << endl;
640 cout <<
"OUTPUT nfakes for asm = " << asm_fakes <<
" +- " << asm_err << endl;
641 cout <<
"OUTPUT nfakes for fkf = " << fkf_fakes <<
" +- " << fkf_err << endl;
642 cout <<
"OUTPUT true_fakes = " << true_fakes << endl;
645 ANA_CHECK(
doMerge(
input,
"lhm",
config, h_lep_pt_lhoodMM, lep_pt, h_lep_eta_lhoodMM, lep_eta, h_lep_pt_eta_lhoodMM, lhoodMM_fakes, lhoodMM_poserr, lhoodMM_negerr, icase) );
646 ANA_CHECK(
doMerge(
input,
"lhf",
config, h_lep_pt_lhoodFF, lep_pt, h_lep_eta_lhoodFF, lep_eta, h_lep_pt_eta_lhoodFF, lhoodFF_fakes, lhoodFF_poserr, lhoodFF_negerr, icase) );
647 ANA_CHECK(
doMerge(
input,
"asm",
config, h_lep_pt_asm, lep_pt, h_lep_eta_asm, lep_eta, h_lep_pt_eta_asm, asm_fakes, asm_err, asm_err, icase) );
648 ANA_CHECK(
doMerge(
input,
"fkf",
config, h_lep_pt_fkf, lep_pt, h_lep_eta_fkf, lep_eta, h_lep_pt_eta_fkf, fkf_fakes, fkf_err, fkf_err, icase) );
655 h_lep_pt_lhoodMM->Write();
656 h_lep_eta_lhoodMM->Write();
657 h_lep_pt_eta_lhoodMM->Write();
658 h_lep_pt_lhoodFF->Write();
659 h_lep_eta_lhoodFF->Write();
660 h_lep_pt_eta_lhoodFF->Write();
661 h_lep_pt_asm->Write();
662 h_lep_eta_asm->Write();
663 h_lep_pt_eta_asm->Write();
664 h_lep_pt_fkf->Write();
665 h_lep_eta_fkf->Write();
666 h_lep_pt_eta_fkf->Write();
668 delete h_lep_pt_lhoodMM;
669 delete h_lep_eta_lhoodMM;
670 delete h_lep_pt_eta_lhoodMM;
671 delete h_lep_pt_lhoodFF;
672 delete h_lep_eta_lhoodFF;
673 delete h_lep_pt_eta_lhoodFF;
675 delete h_lep_eta_asm;
676 delete h_lep_pt_eta_asm;
678 delete h_lep_eta_fkf;
679 delete h_lep_pt_eta_fkf;
688 return StatusCode::SUCCESS;