|
DAS
3.0
Das Analysis System
|
|
template<class Muon > |
bool | keepEvent (const vector< Muon > &muons, float pt1, float pt2) |
|
void | applyDimuonSkim (const vector< fs::path > &inputs, const fs::path &output, const pt::ptree &config, const int steering, const DT::Slice slice={1, 0}) |
|
◆ applyDimuonSkim()
void DAS::Muons::applyDimuonSkim |
( |
const vector< fs::path > & |
inputs, |
|
|
const fs::path & |
output, |
|
|
const pt::ptree & |
config, |
|
|
const int |
steering, |
|
|
const DT::Slice |
slice = {1,0} |
|
) |
| |
Skims the input tree to keep only events with at least two muons, with configurable pT thresholds.
- 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 |
47 cout << __func__ <<
' ' << slice <<
" start" << endl;
50 auto tIn = flow.GetInputTree(
inputs, slice);
51 auto [fOut, tOut] = flow.GetOutput(output);
55 const bool isMC = metainfo.Get<
bool>(
"flags",
"isMC");
57 const auto pt1 =
config.get<
float>(
"skims.dimuon.pt1");
58 metainfo.Set<
float>(
"skims",
"dimuon",
"pt1", pt1);
60 const auto pt2 =
config.get<
float>(
"skims.dimuon.pt2");
61 metainfo.Set<
float>(
"skims",
"dimuon",
"pt2", pt2);
63 cout <<
"Skimming events: pt(mu1) > " << pt1 <<
"\t pt(mu2) > " << pt2 << endl;
66 auto genEvt = flow.GetBranchReadOnly<GenEvent>(
"genEvent",
DT::facultative);
68 auto recMuons = flow.GetBranchReadOnly<vector<RecMuon>>(
"recMuons",
DT::facultative);
70 if (recMuons !=
nullptr)
71 recDimuon = flow.GetBranchWriteOnly<
RecDimuon>(
"recDimuon");
73 auto genMuons = flow.GetBranchReadOnly<vector<GenMuon>>(
"genMuons",
DT::facultative);
75 if (genMuons !=
nullptr)
76 genDimuon = flow.GetBranchWriteOnly<
GenDimuon>(
"genDimuon");
78 if (!recMuons && !genMuons)
79 BOOST_THROW_EXCEPTION(
DE::BadInput(
"No muons in input tree", tIn) );
82 TH1 * hSumWgt =
nullptr;
83 if (
isMC) hSumWgt =
new TH1F(
"hSumWgt",
";;#sum_{i} w_{i}", 1, -0.5, 0.5);
85 for (
DT::Looper looper(tIn); looper(); ++looper) {
89 bool passesRec =
true, passesGen =
true;
92 passesRec =
keepEvent(*recMuons, pt1, pt2);
93 if (passesRec) *recDimuon = recMuons->at(0) + recMuons->at(1);
94 else recDimuon->clear();
97 passesGen =
keepEvent(*genMuons, pt1, pt2);
98 if (passesGen) *genDimuon = genMuons->at(0) + genMuons->at(1);
99 else genDimuon->clear();
102 if ((steering &
DT::fill) && (passesRec || passesGen)) tOut->Fill();
105 auto w = genEvt->weights.front();
106 hSumWgt->Fill(0.0,
w);
109 metainfo.Set<
bool>(
"git",
"complete",
true);
111 if (hSumWgt) hSumWgt->Write();
113 cout << __func__ <<
' ' << slice <<
" stop" << endl;
◆ keepEvent()
bool DAS::Muons::keepEvent |
( |
const vector< Muon > & |
muons, |
|
|
float |
pt1, |
|
|
float |
pt2 |
|
) |
| |
Checks whether an event should be kept.
33 return muons.size() >= 2 &&
muons[0].p4.Pt() > pt1 &&
muons[1].p4.Pt() > pt2;
static bool verbose
Definition: Step.h:40
string muons
Definition: Ntupliser_cfg.py:43
static const float w
Definition: common.h:51
bool keepEvent(const vector< Muon > &muons, float pt1, float pt2)
Checks whether an event should be kept.
Definition: applyDimuonSkim.cc:31
config
Definition: Ntupliser_cfg.py:260
string isMC
Definition: Ntupliser_cfg.py:59
def inputs
Definition: jercExample.py:118
Di< RecMuon, RecMuon > RecDimuon
Definition: Di.h:81
Di< GenMuon, GenMuon > GenDimuon
Definition: Di.h:78