78 StatusCode writeROOT(
const string&
name,
int type,
float realeff_mean,
float fakeeff_mean,
float eff_spread,
float eff_delta_with_pt);
87 std::map<CP::SystematicVariation, float> &syst_map,
88 std::map<CP::SystematicVariation, float> &syst_err_map,
96 std::map<CP::SystematicVariation, float> &syst_map,
97 std::map<CP::SystematicVariation, float> &syst_poserr_map,
98 std::map<CP::SystematicVariation, float> &syst_negerr_map,
103 StatusCode doMerge
ATLAS_NOT_THREAD_SAFE(
const std::vector<std::string> & input,
const std::string &
name,
fbtTestToyMC_config &
config,
TH1F* h_lep_pt,
float &lep_pt,
TH1F* h_lep_eta,
float &lep_eta,
TH2F* h_lep_pt_eta,
float &fakes,
float &poserr,
float &negerr,
int icase);
107 double comboProb(
const vector<FakeBkgTools::ParticleData> & leptons_data,
const std::bitset<64> & tights,
const std::bitset<64> &reals) ;
141 config.realeff_mean = 0.90;
142 config.fakeeff_mean = 0.20;
144 config.eff_delta_with_pt = 0.;
148 config.saveFileNameBase =
"saveProgress";
149 config.mergeFileNameBase =
"saveProgress";
150 config.test_merge =
false;
151 config.test_histo =
false;
152 config.test_systematics =
false;
153 config.poisson_fluctuations =
false;
155 config.selection =
">=1T";
156 config.process =
">=1F[T]";
158 if (
config.verbose) cout <<
"maxnbaseline = " <<
config.maxnbaseline << endl;
168 if (
config.verbose) cout <<
"maxnbaseline = " <<
config.maxnbaseline << endl;
172 int nevents_thiscase;
175 TTree *ntuple =
new TTree(
"FakeBkg",
"Variables from toy MC");
176 ntuple->Branch(
"nevents", &nevents_thiscase,
"nevents/I");
178 ntuple->Branch(
"nevents_sel", &nevents_sel,
"nevents_sel/I");
179 Float_t fake_lep_frac;
180 ntuple->Branch(
"fake_lep_frac", &fake_lep_frac,
"fake_lep_frac/F");
182 ntuple->Branch(
"true_fakes", &true_fakes,
"true_fakes/F");
184 ntuple->Branch(
"asm_fakes", &asm_fakes,
"asm_fakes/F");
186 ntuple->Branch(
"asm_err", &asm_err,
"asm_err/F");
188 ntuple->Branch(
"fkf_fakes", &fkf_fakes,
"fkf_fakes/F");
190 ntuple->Branch(
"fkf_err", &fkf_err,
"fkf_err/F");
191 Float_t lhoodMM_fakes;
192 ntuple->Branch(
"lhoodMM_fakes", &lhoodMM_fakes,
"lhoodMM_fakes/F");
193 Float_t lhoodMM_poserr;
194 ntuple->Branch(
"lhoodMM_poserr", &lhoodMM_poserr,
"lhoodMM_poserr/F");
195 Float_t lhoodMM_negerr;
196 ntuple->Branch(
"lhoodMM_negerr", &lhoodMM_negerr,
"lhoodMM_negerr/F");
197 Float_t lhoodFF_fakes;
198 ntuple->Branch(
"lhoodFF_fakes", &lhoodFF_fakes,
"lhoodFF_fakes/F");
199 Float_t lhoodFF_poserr;
200 ntuple->Branch(
"lhoodFF_poserr", &lhoodFF_poserr,
"lhoodFF_poserr/F");
201 Float_t lhoodFF_negerr;
202 ntuple->Branch(
"lhoodFF_negerr", &lhoodFF_negerr,
"lhoodFF_negerr/F");
205 map<CP::SystematicVariation,float> lhoodMM_weight_map, lhoodMM_poserr_map, lhoodMM_negerr_map;
206 map<CP::SystematicVariation,float> lhoodFF_weight_map, lhoodFF_poserr_map, lhoodFF_negerr_map;
207 map<CP::SystematicVariation,float> asm_weight_map, asm_err_map;
208 map<CP::SystematicVariation,float> fkf_weight_map, fkf_err_map;
210 TRandom3
rand(242868252);
212 Double_t realeff_spread =
config.eff_spread;
213 Double_t fakeeff_spread =
config.eff_spread;
215 std::vector<bool> isTight;
216 std::vector<std::string> input = {
config.outputdirname+
"efficiencies_full" };
226 input.back().append(
".root");
233 for (
unsigned icase(0); icase <
config.ncases; icase++) {
236 TH1F* h_lep_pt_lhoodMM = 0;
237 TH1F* h_lep_eta_lhoodMM = 0;
238 TH2F* h_lep_pt_eta_lhoodMM = 0;
239 TH1F* h_lep_pt_lhoodFF = 0;
240 TH1F* h_lep_eta_lhoodFF = 0;
241 TH2F* h_lep_pt_eta_lhoodFF = 0;
242 TH1F* h_lep_pt_asm = 0;
243 TH1F* h_lep_eta_asm = 0;
244 TH2F* h_lep_pt_eta_asm = 0;
245 TH1F* h_lep_pt_fkf = 0;
246 TH1F* h_lep_eta_fkf = 0;
247 TH2F* h_lep_pt_eta_fkf = 0;
251 h_lep_pt_lhoodMM =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
253 h_lep_eta_lhoodMM =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
255 h_lep_pt_eta_lhoodMM =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
257 h_lep_pt_lhoodFF =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
259 h_lep_eta_lhoodFF =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
261 h_lep_pt_eta_lhoodFF =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
263 h_lep_pt_asm =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
265 h_lep_eta_asm =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
267 h_lep_pt_eta_asm =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
269 h_lep_pt_fkf =
new TH1F(
hname.c_str(),
hname.c_str(), 10, 0, 100);
271 h_lep_eta_fkf =
new TH1F(
hname.c_str(),
hname.c_str(), 10, -5, 5);
273 h_lep_pt_eta_fkf =
new TH2F(
hname.c_str(),
hname.c_str(), 10, 0, 100, 10, -5, 5);
276 float lep_pt = 0, lep_eta = 0;
280 cout <<
"Starting case " << icase << endl;
283 std::vector<CP::LhoodMM_tools*> lhmTool_sav;
285 std::vector<CP::LhoodMM_tools*> lhmTool_FF_sav;
287 std::vector<CP::AsymptMatrixTool*> asmTool_sav;
289 std::vector<CP::ApplyFakeFactor*> fkfTool_sav;
313 ANA_CHECK( lhmTool_FF.setProperty(
"DoFakeFactorFit",
true) );
318 for (
int iSave = 0; iSave <=
nSave; iSave++) {
321 lhmTool_sav.push_back(lhmTool_is);
330 lhmTool_FF_sav.push_back(lhmTool_FF_is);
332 ANA_CHECK( lhmTool_FF_sav[iSave]->setProperty(
"DoFakeFactorFit",
true) );
345 asmTool_sav.push_back(asmTool_is);
350 fkfTool_sav.push_back(fkfTool_is);
361 Double_t realeff_mean_thiscase =
rand.Gaus(
config.realeff_mean, realeff_spread);
362 if (realeff_mean_thiscase > 0.99) realeff_mean_thiscase = 0.99;
363 Double_t fakeeff_mean_thiscase =
rand.Gaus(
config.fakeeff_mean, fakeeff_spread);
364 if (fakeeff_mean_thiscase < 0.01) fakeeff_mean_thiscase = 0.01;
365 if (fakeeff_mean_thiscase > realeff_mean_thiscase) fakeeff_mean_thiscase = realeff_mean_thiscase-0.01;
367 fake_lep_frac =
rand.Uniform();
375 std::vector<LhoodMMEvent> mmevts;
379 std::vector<double> realeff, fakeeff;
381 vector<FinalState*>
fs;
383 for (
unsigned ilep(0); ilep <=
config.maxnbaseline; ilep++) {
386 fs.push_back(this_fs);
389 if (
config.poisson_fluctuations) {
390 nevents_thiscase =
rand.Poisson(
config.nevents);
392 nevents_thiscase =
config.nevents;
399 for (
int pass(0); pass<2; pass++) {
403 int nEventsMultFactor = 1;
404 if (pass > 0) nEventsMultFactor = 10;
405 for (Long64_t ievt(0); ievt<nevents_thiscase*nEventsMultFactor; ievt++) {
407 float nlep_frac = 1./(
config.maxnbaseline-
config.minnbaseline+1);
408 Int_t nlep_select =
config.minnbaseline+
rand.Uniform()/nlep_frac;
410 float extraweight(1.0);
416 vector<FakeBkgTools::ParticleData> leptons_data;
418 vector<bool> lep_real;
421 vector<float> lep_pts;
422 for (
int ilep(0); ilep < nlep_select; ilep++) {
423 lep_pts.push_back(100*
rand.Uniform());
426 sort(lep_pts.begin(), lep_pts.end(), std::greater<float>());
428 for (
int ilep(0); ilep < nlep_select; ilep++) {
432 if (
rand.Uniform() > fake_lep_frac) {
438 lep_real.push_back(isReal);
442 lep_pt = lep_pts[ilep];
443 lep_eta = -5. + 10*
rand.Uniform();
447 ANA_CHECK( lookupEfficiencies(*lepton, lepton_data) );
451 TightAcc(*lepton) =
true;
453 TightAcc(*lepton) =
false;
457 TightAcc(*lepton) =
true;
459 TightAcc(*lepton) =
false;
464 leptons_data.push_back(lepton_data);
473 ANA_CHECK( lhmTool_sav[savIndex]->addEvent(leptons, extraweight) );
474 ANA_CHECK( lhmTool_FF_sav[savIndex]->addEvent(leptons, extraweight) );
475 ANA_CHECK( asmTool_sav[savIndex]->addEvent(leptons, extraweight) );
476 ANA_CHECK( fkfTool_sav[savIndex]->addEvent(leptons, extraweight) );
483 for (
int ilep(0); ilep < nlep_select; ilep++) {
484 if (!lep_real[ilep]) {
491 if (pass == 1 && all_real) {
498 ANA_CHECK( fkfTool_sav[savIndex]->addEvent(leptons, -extraweight/nEventsMultFactor) );
499 ANA_CHECK( lhmTool_FF_sav[savIndex]->addEvent(leptons, -extraweight/nEventsMultFactor) );
505 vector<int> lep_indices;
506 std::bitset<64> tights, reals, charges;
507 for (
int ilep(0); ilep < nlep_select; ilep++) {
508 lep_indices.push_back(ilep);
509 if (lep_real[ilep]) reals.set(ilep);
513 for (
long comb(0); comb < (1<<nlep_select); ++comb) {
514 tights = std::bitset<64>(comb);
515 if (
fs[nlep_select]->accept_selection(tights, charges)) {
516 if (
fs[nlep_select]->accept_process(nlep_select, reals, tights)) {
517 true_fakes += extraweight*
comboProb(leptons_data, tights, reals);
523 for (
int ilep = 0; ilep < nlep_select; ilep++) {
524 if (TightAcc(*leptons[ilep])) tights.set(ilep);
526 if (
fs[nlep_select]->accept_selection(tights,charges) ) {
527 nevents_sel += extraweight;
543 fkfYield -= wgt/nEventsMultFactor;
544 fkfErr += wgt*wgt/(nEventsMultFactor);
550 if (ievt > 0 && (((
nSave*ievt)%nevents_thiscase ==0) || ievt == nevents_thiscase -1) ) {
553 std::unique_ptr<TFile> saveFile(TFile::Open(
saveFileName.c_str(),
"RECREATE"));
554 ANA_CHECK( lhmTool_sav[savIndex]->saveProgress(saveFile->mkdir(
"fakes")) );
558 saveFile = std::make_unique<TFile>(
saveFileName.c_str(),
"RECREATE");
559 ANA_CHECK( asmTool_sav[savIndex]->saveProgress(saveFile.get()) );
565 if (ievt > 0 && (((
nSave*ievt)%(nevents_thiscase*nEventsMultFactor) ==0) || ievt == nevents_thiscase*nEventsMultFactor -1) ) {
569 std::unique_ptr<TFile> saveFile(TFile::Open(
saveFileName.c_str(),
"RECREATE"));
570 ANA_CHECK( lhmTool_FF_sav[savIndex]->saveProgress(saveFile->mkdir(
"fakes")) );
574 saveFile = std::make_unique<TFile>(
saveFileName.c_str(),
"RECREATE");
575 ANA_CHECK( fkfTool_sav[savIndex]->saveProgress(saveFile.get()) );
579 float nfakes_tmp(0), err_tmp(0);
586 if (*
it !=
nullptr)
delete *
it;
591 asm_err = sqrt(asmErr);
592 asm_fakes = asmYield;
602 fkf_fakes = fkfYield;
604 if (
config.test_systematics) {
606 std::string systBrName, systBrNameF;
608 for(
auto& sysvar : sysvars){
611 float asm_syst_weight, asm_syst_err;
612 float fkf_syst_weight, fkf_syst_err;
613 float lhoodMM_syst_weight, lhoodMM_syst_poserr, lhoodMM_syst_negerr;
614 float lhoodFF_syst_weight, lhoodFF_syst_poserr, lhoodFF_syst_negerr;
619 ANA_CHECK(
setupSystBranchesAsym(
"lhoodMM", sysvar, lhoodMM_syst_weight, lhoodMM_syst_poserr, lhoodMM_syst_negerr, lhoodMM_weight_map, lhoodMM_poserr_map, lhoodMM_negerr_map, ntuple) );
620 ANA_CHECK(
setupSystBranchesAsym(
"lhoodFF", sysvar, lhoodFF_syst_weight, lhoodFF_syst_poserr, lhoodFF_syst_negerr, lhoodFF_weight_map, lhoodFF_poserr_map, lhoodFF_negerr_map, ntuple) );
624 asm_err_map.find(sysvar)->second,
625 asm_err_map.find(sysvar)->second));
629 fkf_err_map.find(sysvar)->second,
630 fkf_err_map.find(sysvar)->second) );
633 lhoodMM_poserr_map.find(sysvar)->second,
634 lhoodMM_negerr_map.find(sysvar)->second) );
637 lhoodFF_poserr_map.find(sysvar)->second,
638 lhoodFF_negerr_map.find(sysvar)->second) );
643 cout <<
"OUTPUT true_fakes = " << true_fakes << endl;
644 cout <<
"OUTPUT nfakes for lhoodMM = " << lhoodMM_fakes <<
" +" << lhoodMM_poserr <<
" -" << lhoodMM_negerr << endl;
645 cout <<
"OUTPUT nfakes for lhoodFF = " << lhoodFF_fakes <<
" +" << lhoodFF_poserr <<
" -" << lhoodFF_negerr << endl;
646 cout <<
"OUTPUT nfakes for asm = " << asm_fakes <<
" +- " << asm_err << endl;
647 cout <<
"OUTPUT nfakes for fkf = " << fkf_fakes <<
" +- " << fkf_err << endl;
648 cout <<
"OUTPUT true_fakes = " << true_fakes << endl;
651 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) );
652 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) );
653 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) );
654 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) );
661 h_lep_pt_lhoodMM->Write();
662 h_lep_eta_lhoodMM->Write();
663 h_lep_pt_eta_lhoodMM->Write();
664 h_lep_pt_lhoodFF->Write();
665 h_lep_eta_lhoodFF->Write();
666 h_lep_pt_eta_lhoodFF->Write();
667 h_lep_pt_asm->Write();
668 h_lep_eta_asm->Write();
669 h_lep_pt_eta_asm->Write();
670 h_lep_pt_fkf->Write();
671 h_lep_eta_fkf->Write();
672 h_lep_pt_eta_fkf->Write();
674 delete h_lep_pt_lhoodMM;
675 delete h_lep_eta_lhoodMM;
676 delete h_lep_pt_eta_lhoodMM;
677 delete h_lep_pt_lhoodFF;
678 delete h_lep_eta_lhoodFF;
679 delete h_lep_pt_eta_lhoodFF;
681 delete h_lep_eta_asm;
682 delete h_lep_pt_eta_asm;
684 delete h_lep_eta_fkf;
685 delete h_lep_pt_eta_fkf;
694 return StatusCode::SUCCESS;
698 TRandom3 rnd(235789);
702 std::unique_ptr<TFile>
file(TFile::Open(
name.c_str(),
"RECREATE"));
706 TH1D hElFake(
"FakeEfficiency_el_pt",
"FakeEfficiency_el_pt", nbin, 0., 1.*nbin);
707 TH1D hMuFake(
"FakeEfficiency_mu_pt",
"FakeEfficiency_mu_pt", nbin, 0., 1.*nbin);
708 TH1D hElReal(
"RealEfficiency_el_pt",
"RealEfficiency_el_pt", nbin, 0., 1.*nbin);
709 TH1D hMuReal(
"RealEfficiency_mu_pt",
"RealEfficiency_mu_pt", nbin, 0., 1.*nbin);
711 TH1D hElFake_bigSyst(
"FakeEfficiency_el_pt__bigSyst",
"FakeEfficiency_el_pt__bigSyst", nbin, 0., 1.*nbin);
712 TH1D hElFake_smallSyst(
"FakeEfficiency_el_pt__smallSyst",
"FakeEfficiency_el_pt__smallSyst", nbin, 0., 1.*nbin);
713 TH1D hMuFake_bigSyst(
"FakeEfficiency_mu_pt__bigSyst",
"FakeEfficiency_mu_pt__bigSyst", nbin, 0., 1.*nbin);
714 TH1D hMuFake_smallSyst(
"FakeEfficiency_mu_pt__smallSyst",
"FakeEfficiency_mu_pt__smallSyst", nbin, 0., 1.*nbin);
715 TH1D hElReal_bigSyst(
"RealEfficiency_el_pt__bigSyst",
"RealEfficiency_el_pt__bigSyst", nbin, 0., 1.*nbin);
716 TH1D hMuReal_bigSyst(
"RealEfficiency_mu_pt__bigSyst",
"RealEfficiency_mu_pt__bigSyst", nbin, 0., 1.*nbin);
718 for (
int ibin = 1; ibin <= nbin; ibin++) {
719 Double_t realeff = TMath::Min(rnd.Gaus(realeff_mean, eff_spread)+eff_delta_with_pt*(ibin-nbin/2), 0.99);
720 Double_t fakeeff = TMath::Max(rnd.Gaus(fakeeff_mean, eff_spread)-eff_delta_with_pt*(ibin-nbin/2), 0.01);
722 float minrfdiff = 0.10;
723 if (realeff - fakeeff < minrfdiff) {
724 if (realeff > minrfdiff) {
725 fakeeff = realeff - minrfdiff;
727 realeff = realeff + 0.10;
732 if (realeff < 0 || realeff > 1.) cout <<
"ERROR: Bad real efficiency value: " << realeff << endl;
733 if (fakeeff < 0 || fakeeff > 1.) cout <<
"ERROR: Bad fake efficiency value: " << fakeeff << endl;
734 if ((realeff - fakeeff) < minrfdiff) cout <<
"ERROR: Too small difference between real and fake efficiencies: " << realeff <<
" " << fakeeff << endl;
736 hElFake.SetBinContent(ibin, fakeeff);
737 hElFake.SetBinError(ibin, eff_spread);
739 hMuFake.SetBinContent(ibin, fakeeff);
740 hMuFake.SetBinError(ibin, eff_spread);
742 hElReal.SetBinContent(ibin, realeff);
743 hElReal.SetBinError(ibin, eff_spread);
745 hMuReal.SetBinContent(ibin, realeff);
746 hMuReal.SetBinError(ibin, eff_spread);
749 hElFake_bigSyst.SetBinContent(ibin, 0.20*hElFake.GetBinContent(ibin));
750 hElFake_smallSyst.SetBinContent(ibin, 0.02*hElFake.GetBinContent(ibin));
751 hMuFake_bigSyst.SetBinContent(ibin, 0.20*hMuFake.GetBinContent(ibin));
752 hMuFake_smallSyst.SetBinContent(ibin, 0.02*hMuFake.GetBinContent(ibin));
753 hElReal_bigSyst.SetBinContent(ibin, 0.20*hElReal.GetBinContent(ibin));
754 hMuReal_bigSyst.SetBinContent(ibin, 0.20*hMuReal.GetBinContent(ibin));
763 hElFake_bigSyst.Write();
764 hElReal_bigSyst.Write();
765 hElFake_smallSyst.Write();
766 hMuFake_bigSyst.Write();
767 hMuFake_smallSyst.Write();
768 hMuReal_bigSyst.Write();
771 return StatusCode::SUCCESS;
777 {
"ncases", required_argument, 0,
'c'},
778 {
"nevents", required_argument, 0,
'e'},
779 {
"minnb", required_argument, 0,
'm'},
780 {
"maxnb", required_argument, 0,
'n'},
781 {
"rmean", required_argument, 0,
'r'},
782 {
"fmean", required_argument, 0,
'f'},
783 {
"effsigma", required_argument, 0,
's'},
784 {
"effdeltawithpt", required_argument, 0,
'd'},
785 {
"sel", required_argument, 0,
'l'},
786 {
"proc", required_argument, 0,
'p'},
787 {
"test_save", no_argument, 0,
'S'},
788 {
"test_merge", required_argument, 0,
'M'},
789 {
"test_histo", required_argument, 0,
'H'},
790 {
"test_systematics", no_argument, 0,
'E'},
791 {
"verbose", no_argument, 0,
'v'},
792 {
"poisson", no_argument, 0,
'P'},
793 {
"help", no_argument, 0,
'h'},
801 while ((
c = getopt_long(
argc,
argv,
"c:e:m:n:r:f:s:d:l:p:S:M:HEPvh",
long_options, &longindex)) != -1) {
822 config.selection = optarg;
835 config.saveFileNameBase = optarg;
839 config.mergeFileNameBase = optarg;
845 config.test_systematics =
true;
851 config.poisson_fluctuations =
true;
868 while (
pos != string::npos) {
873 while (
pos != string::npos) {
877 return StatusCode::SUCCESS;
881 string outputdirname;
883 outputdirname =
"FakeBkgTools_toy_MC_nevents_";
885 outputdirname +=
"_ncases_";
887 outputdirname +=
"_minnbaseline_";
889 outputdirname +=
"_maxnbaseline_";
891 outputdirname +=
"_realeff_";
893 outputdirname +=
"_fakeeff_";
895 outputdirname +=
"_effspread_";
897 outputdirname +=
"_selection_";
898 outputdirname +=
config.selection;
899 outputdirname +=
"_process_";
900 outputdirname +=
config.process;
901 if (
config.poisson_fluctuations) {
902 outputdirname +=
"_poisson_";
908 pos = outputdirname.find(
">=");
909 while (
pos != string::npos) {
910 outputdirname.replace(
pos, 2,
"ge");
911 pos = outputdirname.find(
">=");
913 pos = outputdirname.find(
"<=");
914 while (
pos != string::npos) {
915 outputdirname.replace(
pos, 2,
"le");
916 pos = outputdirname.find(
"<=");
918 pos = outputdirname.find(
">");
919 while (
pos != string::npos) {
920 outputdirname.replace(
pos, 1,
"gt");
921 pos = outputdirname.find(
">");
923 pos = outputdirname.find(
"<");
924 while (
pos != string::npos) {
925 outputdirname.replace(
pos, 1,
"lt");
926 pos = outputdirname.find(
"<");
928 pos = outputdirname.find(
"=");
929 while (
pos != string::npos) {
930 outputdirname.replace(
pos, 1,
"eq");
931 pos = outputdirname.find(
"=");
933 pos = outputdirname.find(
",");
934 while (
pos != string::npos) {
935 outputdirname.replace(
pos, 1,
"");
936 pos = outputdirname.find(
",");
938 pos = outputdirname.find(
"[");
939 while (
pos != string::npos) {
940 outputdirname.replace(
pos, 1,
"");
941 pos = outputdirname.find(
"[");
943 pos = outputdirname.find(
"]");
944 while (
pos != string::npos) {
945 outputdirname.replace(
pos, 1,
"");
946 pos = outputdirname.find(
"]");
948 pos = outputdirname.find(
"\"");
949 while (
pos != string::npos) {
950 outputdirname.replace(
pos, 1,
"");
951 pos = outputdirname.find(
"\"");
953 pos = outputdirname.find(
"\\");
954 while (
pos != string::npos) {
955 outputdirname.replace(
pos, 1,
"");
956 pos = outputdirname.find(
"\\");
959 gSystem->mkdir(outputdirname.c_str());
960 rootfilename = outputdirname+
"/output.root";
962 config.outputdirname = outputdirname;
964 std::unique_ptr<TFile> f_out(TFile::Open(rootfilename.c_str(),
"RECREATE"));
983 return StatusCode::SUCCESS;
989 cout <<
"Um, no ROOT file!" << endl;
990 return StatusCode::FAILURE;
1003 if (
h_realeff_e == 0) cout <<
"No real e" << endl;
1004 if (
h_fakeeff_e == 0) cout <<
"No fake e" << endl;
1010 return StatusCode::SUCCESS;
1015 if (
h_realeff_e == 0) cout <<
"No real e" << endl;
1016 if (
h_fakeeff_e == 0) cout <<
"No fake e" << endl;
1026 return StatusCode::SUCCESS;
1029 double comboProb(
const vector<FakeBkgTools::ParticleData> & leptons_data,
const std::bitset<64> & tights,
const std::bitset<64> &reals) {
1032 for (
unsigned ilep = 0; ilep < leptons_data.size(); ilep++) {
1035 prob *= leptons_data[ilep].real_efficiency.nominal;
1037 prob *= leptons_data[ilep].fake_efficiency.nominal;
1041 prob *= (1.-leptons_data[ilep].real_efficiency.nominal);
1043 prob *= (1.-leptons_data[ilep].fake_efficiency.nominal);
1051 double comboProb_FF(
const vector<FakeBkgTools::ParticleData> & leptons_data,
const std::bitset<64> &tights,
const std::bitset<64> & reals) {
1055 for (
unsigned ilep = 0; ilep < leptons_data.size(); ilep++) {
1059 prob *= leptons_data[ilep].fake_efficiency.nominal;
1065 prob *= (1.-leptons_data[ilep].fake_efficiency.nominal);
1069 for (
unsigned ilep = 0; ilep < leptons_data.size(); ilep++) {
1071 float F = leptons_data[ilep].fake_efficiency.nominal/(1-leptons_data[ilep].fake_efficiency.nominal);
1082 std::map<CP::SystematicVariation, float> &syst_map,
1083 std::map<CP::SystematicVariation, float> &syst_err_map,
1085 syst_map.emplace(std::make_pair(sysvar,
weight));
1086 std::string systBrName = baseName;
1087 systBrName = systBrName+
"_fakes_"+sysvar.
name();
1088 std::string systBrNameF = systBrName+
"/F";
1089 ntuple->Branch(systBrName.c_str(), &(syst_map.find(sysvar)->second), systBrNameF.c_str());
1090 syst_err_map.emplace(std::make_pair(sysvar, weight_err));
1091 systBrName = baseName;
1092 systBrName = systBrName+
"_err_"+sysvar.
name();
1093 systBrNameF = systBrName+
"/F";
1094 ntuple->Branch(systBrName.c_str(), &(syst_err_map.find(sysvar)->second), systBrNameF.c_str());
1096 return StatusCode::SUCCESS;
1102 float &weight_poserr,
1103 float &weight_negerr,
1104 std::map<CP::SystematicVariation, float> &syst_map,
1105 std::map<CP::SystematicVariation, float> &syst_poserr_map,
1106 std::map<CP::SystematicVariation, float> &syst_negerr_map,
1108 syst_map.emplace(std::make_pair(sysvar,
weight));
1109 std::string systBrName = baseName;
1110 systBrName = systBrName+
"_fakes_"+sysvar.
name();
1111 std::string systBrNameF = systBrName+
"/F";
1112 ntuple->Branch(systBrName.c_str(), &(syst_map.find(sysvar)->second), systBrNameF.c_str());
1113 syst_poserr_map.emplace(std::make_pair(sysvar, weight_poserr));
1114 systBrName = baseName;
1115 systBrName = systBrName+
"_poserr_"+sysvar.
name();
1116 systBrNameF = systBrName+
"/F";
1117 ntuple->Branch(systBrName.c_str(), &(syst_poserr_map.find(sysvar)->second), systBrNameF.c_str());
1118 syst_negerr_map.emplace(std::make_pair(sysvar, weight_negerr));
1119 systBrName = baseName;
1120 systBrName = systBrName+
"_negerr_"+sysvar.
name();
1121 systBrNameF = systBrName+
"/F";
1122 ntuple->Branch(systBrName.c_str(), &(syst_negerr_map.find(sysvar)->second), systBrNameF.c_str());
1124 return StatusCode::SUCCESS;
1127 StatusCode doMerge
ATLAS_NOT_THREAD_SAFE(
const std::vector<std::string> &
input,
const std::string &
name,
fbtTestToyMC_config &
config,
TH1F* h_lep_pt,
float &lep_pt,
TH1F* h_lep_eta,
float &lep_eta,
TH2F* h_lep_pt_eta,
float &fakes,
float &poserr,
float &negerr,
int icase) {
1130 system(haddcmd.c_str());
1132 std::unique_ptr<CP::BaseFakeBkgTool>
tool;
1134 if (
name ==
"lhm" ||
name ==
"lhf") {
1135 std::string toolName =
"Lhood";
1136 if (
name ==
"lhm") {
1141 toolName +=
to_string(icase)+
"_tools_merge";
1142 tool = std::make_unique<CP::LhoodMM_tools>(toolName);
1143 if (
name ==
"lhf") {
1146 ANA_CHECK(
tool->setProperty(
"ProgressFileDirectory",
"fakes") );
1147 }
else if (
name ==
"asm") {
1148 tool = std::make_unique<CP::AsymptMatrixTool>(
"asm_tool_merge");
1149 }
else if (
name ==
"fkf") {
1150 tool = std::make_unique<CP::ApplyFakeFactor>(
"fkf_tool_merge");
1153 std::string mergeFileName =
config.mergeFileNameBase+
"_"+
name+
"_"+
to_string(icase)+
".root";
1154 std::cout << mergeFileName << std::endl;
1155 ANA_CHECK(
tool->setProperty(
"ProgressFileName", mergeFileName) );
1159 ANA_CHECK(
tool->register1DHistogram(h_lep_eta, &lep_eta) );
1160 ANA_CHECK(
tool->register2DHistogram(h_lep_pt_eta, &lep_pt, &lep_eta) );
1165 return StatusCode::SUCCESS;
1169 std::cout << std::endl <<
"fbtTestToyMC provides a toy MC model of a set of pseudoexperiments. It is intended to explore the statistical properties of the fake lepton background estimate for a given set of experiemental conditions (number of events, number of leptons required, typical values of the real and fake efficiencies, etc.)" << std::endl << std::endl;
1170 std::cout <<
"Options: " << std::endl;
1171 std::cout <<
" --ncases, -c: number of pseudoexperiments to run (default: 100)" << std::endl;
1172 std::cout <<
" --nevents, -e: number of events in each pseudoexperiment (default: 100)" << std::endl;
1173 std::cout <<
" --poisson, -P: use Poisson-distributed number of events in each pseudoexperiment, rather than a fixed number (default: false)" << std::endl;
1174 std::cout <<
" --minnb, -m: minimum number of baseline leptons per event (default: 1)" << std::endl;
1175 std::cout <<
" --maxnb, -m: maximum number of baseline leptons per event (default: 4)" << std::endl;
1176 std::cout <<
" --rmean, -r: average real lepton efficiency (default: 0.9)" << std::endl;
1177 std::cout <<
" --fmean, -f: average fake lepton efficiency (default: 0.2)" << std::endl;
1178 std::cout <<
" --effsigma, -s: standard deviation for lepton-to-lepton variation in real and fake efficiencies (default: 0.10)" << std::endl;
1179 std::cout <<
" --effdeltawithpt, -d: rate at which the real (fake) efficiency increases (decreases) with simulated pt. Note that the simulated pt range is 0 to 100" << std::endl;
1180 std::cout <<
" --sel, -l: selection string to be used by the tools (default \" >= 1T\") " << std::endl;
1181 std::cout <<
" --proc, -p: process string to be used by the tools (default \" >= 1F[T]\") " << std::endl;
1182 std::cout <<
" --test_save, -S: save output from subjobs in each pseudoexperiement. If set, requires an argument to specify the base name of the root files where the output will be saved. (default: false) " << std::endl;
1183 std::cout <<
" --test_merge, -M: merge output from subjobs in each pseudoexperiement. If set, requires an argument to specify the base name of the root files to be merged. This should match the base name used in a previous run with the --test_save option (default: false) " << std::endl;
1184 std::cout <<
" --test_histo, -H: test the filling of 1- and 2-dimensional histograms (default: false) " << std::endl;
1185 std::cout <<
" --test_systematics, -E: test the handling of systematic unceratinties (default: false) " << std::endl;
1186 std::cout <<
" --verbose, -v: enable verbose output (default: false) " << std::endl;
1187 std::cout <<
" --help, -h: print this help message" << std::endl;
1188 return StatusCode::SUCCESS;;