DAS  3.0
Das Analysis System
HTnFillerfinal

#include <HTn.h>

+ Inheritance diagram for HTnFiller:
+ Collaboration diagram for HTnFiller:

Public Member Functions

 HTnFiller (const HTn &obs, TTreeReader &reader)
 
std::list< int > fillRec (DistVariation &) override
 
void match () override
 
void fillMC (DistVariation &) override
 
- Public Member Functions inherited from Filler
virtual ~Filler ()=default
 

Public Attributes

HTn obs
 
std::optional< TTreeReaderArray< GenJet > > genJets
 
TTreeReaderArray< RecJetrecJets
 
std::optional< TTreeReaderValue< GenEvent > > gEv
 
TTreeReaderValue< RecEventrEv
 
std::vector< std::pair< GenJet, RecJet > > matches
 
std::vector< GenJetmisses
 
std::vector< RecJetfakes
 
std::optional< bool > matched
 

Constructor & Destructor Documentation

◆ HTnFiller()

HTnFiller ( const HTn obs,
TTreeReader &  reader 
)

Constructor.

26  : obs(obs)
27  , genJets(initOptionalBranch<decltype(genJets)>(reader, "genJets"))
28  , recJets(reader, "recJets")
29  , gEv(initOptionalBranch<decltype(gEv)>(reader, "genEvent"))
30  , rEv(reader, "recEvent")
31 {
32 }

Member Function Documentation

◆ fillMC()

void fillMC ( DistVariation v)
overridevirtual

See Filler::fillMC

Reimplemented from Filler.

128 {
129  if (!obs.isMC)
130  BOOST_THROW_EXCEPTION( runtime_error(__func__ + " should only be called for MC"s) );
131 
132  auto recdijet = selection<RecJet>(recJets, v);
133  auto [recn,recW] = getMultiplicity<RecJet>(recJets, v);
134  auto irec = getBinNumber<RecJet>(recdijet, recn, v, obs.recBinning);
135  recW *= rEv->Weight(v);
136  if (irec > 0) recW *= recdijet.Weight(v);
137 
138  auto gendijet = selection<GenJet>(*genJets, v);
139  auto [genn,genW] = getMultiplicity<GenJet>(*genJets, v);
140  auto igen = getBinNumber<GenJet>(gendijet, genn, v, obs.genBinning);
141  genW = (*gEv)->Weight(v);
142  if (igen > 0) genW *= gendijet.Weight(v);
143 
144  if (igen > 0) v.gen->Fill(igen, genW);
145 
146  if (*matched) {
147  if (irec > 0 && igen > 0) { v.RM->Fill(igen, irec, genW * recW );
148  v.missNoMatch->Fill(igen, genW * (1-recW)); }
149  else if (irec == 0 && igen > 0) v.missOut->Fill(igen, genW );
150  else if (irec > 0 && igen == 0) v.fakeOut->Fill( irec, genW * recW );
151  }
152  else {
153  if (igen > 0) v.missNoMatch->Fill(igen, genW );
154  if (irec > 0) v.fakeNoMatch->Fill(irec, genW * recW);
155  }
156 }

◆ fillRec()

list< int > fillRec ( DistVariation v)
overridevirtual

See Filler::fillRec

Reimplemented from Filler.

109 {
110  auto dijet = selection<RecJet>(recJets, v);
111  if (!dijet) return {};
112 
113  auto [n,w] = getMultiplicity(recJets, v);
114  int i = getBinNumber(dijet, n, v, obs.recBinning);
115  if (i == 0) return {};
116 
117  w *= rEv->Weight(v);
118  if (obs.isMC) w *= (*gEv)->Weight(v);
119  w *= dijet.Weight(v);
120 
121  v.tmp->Fill(i, w);
122  v.rec->Fill(i, w);
123 
124  return list<int>{i};
125 }

◆ match()

void match ( )
overridevirtual

Match the two pairs of leading jets (if any) and set matched member.

Reimplemented from Filler.

35 {
36  matched.reset();
37 
38  auto match = [this](size_t i, size_t j) {
39  const FourVector& g = genJets->At(i).p4,
40  r = recJets.At(j).p4;
41  using ROOT::Math::VectorUtil::DeltaR;
42  auto DR = DeltaR(g, r);
43  //cout << g << '\t' << r << '\t' << DR << '\t' << result << '\n';
44  return DR < obs.maxDR;
45  };
46 
47  // matching (swapping leading and subleading is allowed)
48  matched = genJets->GetSize() > 1 && recJets.GetSize() > 1
49  && ( (match(0,0) && match(1,1))
50  || (match(0,1) && match(1,0)) );
51 }

Member Data Documentation

◆ fakes

std::vector<RecJet> fakes

◆ genJets

std::optional<TTreeReaderArray<GenJet> > genJets

◆ gEv

std::optional<TTreeReaderValue<GenEvent> > gEv

◆ matched

std::optional<bool> matched

◆ matches

std::vector<std::pair<GenJet,RecJet> > matches

◆ misses

std::vector<GenJet> misses

◆ obs

HTn obs

Backreference to the observable.

◆ recJets

TTreeReaderArray<RecJet> recJets

◆ rEv

TTreeReaderValue<RecEvent> rEv

The documentation for this struct was generated from the following files:
getBinNumber
double getBinNumber(const Di< Jet, Jet > &MNjets, const Uncertainties::Variation &v, TUnfoldBinning *bng)
Definition: MNjets.cc:47
DAS::Di::Weight
double Weight(const Uncertainties::Variation &v=Uncertainties::nominal) const override
Definition: Di.h:73
DAS::Unfolding::DistVariation::rec
std::unique_ptr< TH1 > rec
reconstructed-level distribution
Definition: DistVariation.h:32
DAS::Unfolding::Rij::HTnFiller::rEv
TTreeReaderValue< RecEvent > rEv
Definition: HTn.h:76
DAS::Unfolding::Rij::HTnFiller::obs
HTn obs
Backreference to the observable.
Definition: HTn.h:71
DAS::JetEnergy::w
static const float w
Definition: common.h:47
DAS::Unfolding::Observable::genBinning
TUnfoldBinning * genBinning
particle-level binning
Definition: Observable.h:150
DAS::Unfolding::Rij::HTnFiller::recJets
TTreeReaderArray< RecJet > recJets
Definition: HTn.h:74
DAS::Unfolding::Rij::HTnFiller::gEv
std::optional< TTreeReaderValue< GenEvent > > gEv
Definition: HTn.h:75
DAS::Unfolding::Rij::HTnFiller::genJets
std::optional< TTreeReaderArray< GenJet > > genJets
Definition: HTn.h:73
DAS::Unfolding::Rij::HTnFiller::matched
std::optional< bool > matched
Definition: HTn.h:90
DAS::Unfolding::Observable::isMC
static bool isMC
Definition: Observable.h:144
DAS::Unfolding::initOptionalBranch
auto initOptionalBranch(TTreeReader &reader, const char *name)
Definition: Observable.h:39
recdijet
DAS::RecDijet recdijet
Definition: classes.h:34
DAS::Unfolding::DistVariation::gen
std::unique_ptr< TH1 > gen
generated-level distribution
Definition: DistVariation.h:34
DAS::Unfolding::DistVariation::tmp
std::unique_ptr< TH1 > tmp
temporary histogram help fill the covariance matrix
Definition: DistVariation.h:33
DAS::Unfolding::DistVariation::missNoMatch
std::unique_ptr< TH1 > missNoMatch
losses (unmatched entries)
Definition: DistVariation.h:35
DAS::Unfolding::DistVariation::RM
std::unique_ptr< TH2 > RM
response matrix
Definition: DistVariation.h:42
DAS::Unfolding::DistVariation::missOut
std::unique_ptr< TH1 > missOut
losses (migration out of phase space)
Definition: DistVariation.h:36
DAS::Unfolding::DistVariation::fakeNoMatch
std::unique_ptr< TH1 > fakeNoMatch
background (unmatched entries)
Definition: DistVariation.h:37
DAS::FourVector
ROOT::Math::LorentzVector< ROOT::Math::PtEtaPhiM4D< float > > FourVector
Definition: PhysicsObject.h:15
DAS::Unfolding::Observable::recBinning
TUnfoldBinning * recBinning
detector-level binning
Definition: Observable.h:149
DAS::Unfolding::DistVariation::fakeOut
std::unique_ptr< TH1 > fakeOut
background (migration out of phase space)
Definition: DistVariation.h:38
gendijet
DAS::GenDijet gendijet
Definition: classes.h:31
DAS::Unfolding::Rij::HTnFiller::match
void match() override
Match the two pairs of leading jets (if any) and set matched member.
Definition: HTn.cc:34
DAS::Unfolding::Observable::maxDR
static double maxDR
max Delta R
Definition: Observable.h:146