CMS 3D CMS Logo

List of all members | Public Member Functions | Private Attributes
calcTopMass Class Reference
Inheritance diagram for calcTopMass:
edm::one::EDAnalyzer< edm::one::SharedResources > edm::one::EDAnalyzerBase edm::EDConsumerBase

Public Member Functions

void analyze (const edm::Event &iEvent, const edm::EventSetup &iSetup) override
 
 calcTopMass (const edm::ParameterSet &)
 
 ~calcTopMass () override
 
- Public Member Functions inherited from edm::one::EDAnalyzer< edm::one::SharedResources >
 EDAnalyzer ()=default
 
 EDAnalyzer (const EDAnalyzer &)=delete
 
SerialTaskQueueglobalLuminosityBlocksQueue () final
 
SerialTaskQueueglobalRunsQueue () final
 
const EDAnalyzeroperator= (const EDAnalyzer &)=delete
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsInputProcessBlocks () const final
 
bool wantsProcessBlocks () const final
 
- Public Member Functions inherited from edm::one::EDAnalyzerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzerBase ()
 
ModuleDescription const & moduleDescription () const
 
bool wantsStreamLuminosityBlocks () const
 
bool wantsStreamRuns () const
 
 ~EDAnalyzerBase () override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
std::vector< ESProxyIndex > const & esGetTokenIndicesVector (edm::Transition iTrans) const
 
std::vector< ESRecordIndex > const & esGetTokenRecordIndicesVector (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::array< std::vector< ModuleDescription const *> *, NumBranchTypes > &modulesAll, std::vector< ModuleProcessName > &modulesInPreviousProcesses, ProductRegistry const &preg, std::map< std::string, ModuleDescription const *> const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void selectInputProcessBlocks (ProductRegistry const &productRegistry, ProcessBlockHelperBase const &processBlockHelperBase)
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Attributes

float bMass
 
float cMass
 
TH1F * hMassCorFl0
 
TH1F * hMassCorFlM
 
TH1F * hMassCorMix
 
TH1F * hMassNoCorr
 
std::string m_bT_CorrectorName
 
std::string m_cT_CorrectorName
 
std::string m_qT_CorrectorName
 
std::string m_tT_CorrectorName
 
float qMass
 
edm::InputTag sourceByRefer_
 
edm::InputTag sourcePartons_
 
edm::Handle< reco::JetMatchedPartonsCollectiontheJetPartonMatch
 

Additional Inherited Members

- Public Types inherited from edm::one::EDAnalyzerBase
typedef EDAnalyzerBase ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::one::EDAnalyzerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
template<BranchType B = InEvent>
EDConsumerBaseAdaptor< Bconsumes (edm::InputTag tag) noexcept
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<Transition Tr = Transition::Event>
constexpr auto esConsumes ()
 
template<Transition Tr = Transition::Event>
auto esConsumes (ESInputTag tag)
 
template<Transition Tr = Transition::Event>
ESGetTokenGeneric esConsumes (eventsetup::EventSetupRecordKey const &iRecord, eventsetup::DataKey const &iKey)
 Used with EventSetupRecord::doGet. More...
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
void resetItemsToGetFrom (BranchType iType)
 

Detailed Description

Definition at line 55 of file CalcTopMassExample.cc.

Constructor & Destructor Documentation

◆ calcTopMass()

calcTopMass::calcTopMass ( const edm::ParameterSet iConfig)
explicit

Definition at line 86 of file CalcTopMassExample.cc.

References compareTotals::fs, edm::ParameterSet::getParameter(), TFileService::kSharedResource, and AlCaHLTBitMon_QueryRunRegistry::string.

86  {
87  sourceByRefer_ = iConfig.getParameter<InputTag>("srcByReference");
88  m_qT_CorrectorName = iConfig.getParameter<std::string>("qTopCorrector");
89  m_cT_CorrectorName = iConfig.getParameter<std::string>("cTopCorrector");
90  m_bT_CorrectorName = iConfig.getParameter<std::string>("bTopCorrector");
91  m_tT_CorrectorName = iConfig.getParameter<std::string>("tTopCorrector");
92 
93  bMass = 4.5;
94  cMass = 1.5;
95  qMass = 0.3;
96 
97  usesResource(TFileService::kSharedResource);
99  hMassNoCorr = fs->make<TH1F>("hMassNoCorr", "", 100, 100, 300);
100  hMassCorFl0 = fs->make<TH1F>("hMassCorFl0", "", 100, 100, 300);
101  hMassCorFlM = fs->make<TH1F>("hMassCorFlM", "", 100, 100, 300);
102  hMassCorMix = fs->make<TH1F>("hMassCorMix", "", 100, 100, 300);
103 }
static const std::string kSharedResource
Definition: TFileService.h:76
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
std::string m_cT_CorrectorName
edm::InputTag sourceByRefer_
std::string m_qT_CorrectorName
std::string m_tT_CorrectorName
std::string m_bT_CorrectorName

◆ ~calcTopMass()

calcTopMass::~calcTopMass ( )
inlineoverride

Definition at line 58 of file CalcTopMassExample.cc.

58 {};

Member Function Documentation

◆ analyze()

void calcTopMass::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
overridevirtual

Implements edm::one::EDAnalyzerBase.

Definition at line 105 of file CalcTopMassExample.cc.

References EcnaPython_AdcPeg12_S1_10_R170298_1_0_150_Dee0::cerr, reco::JetCorrector::correction(), gather_cfg::cout, cppFunctionSkipper::exception, edm::Ref< C, T, F >::get(), JetCorrector::getJetCorrector(), iEvent, edm::Ref< C, T, F >::isNonnull(), dqmiolumiharvest::j, and reco::MatchedPartons::physicsDefinitionParton().

105  {
106  cout << "[calcTopMass] analysing event " << iEvent.id() << endl;
107  try {
109  } catch (std::exception &ce) {
110  cerr << "[calcTopMass] caught std::exception " << ce.what() << endl;
111  return;
112  }
113 
114  // get all correctors from top events
115  const JetCorrector *qTopCorrector = JetCorrector::getJetCorrector(m_qT_CorrectorName, iSetup);
116  const JetCorrector *cTopCorrector = JetCorrector::getJetCorrector(m_cT_CorrectorName, iSetup);
117  const JetCorrector *bTopCorrector = JetCorrector::getJetCorrector(m_bT_CorrectorName, iSetup);
118  const JetCorrector *tTopCorrector = JetCorrector::getJetCorrector(m_tT_CorrectorName, iSetup);
119 
120  TLorentzVector jetsNoCorr[6];
121  TLorentzVector jetsCorFl0[6];
122  TLorentzVector jetsCorFlM[6];
123  TLorentzVector jetsCorMix[6];
124 
125  bool isQuarkFound[6] = {false};
126 
127  for (JetMatchedPartonsCollection::const_iterator j = theJetPartonMatch->begin(); j != theJetPartonMatch->end(); j++) {
128  const math::XYZTLorentzVector theJet = (*j).first.get()->p4();
129  const MatchedPartons aMatch = (*j).second;
130  const GenParticleRef &thePhyDef = aMatch.physicsDefinitionParton();
131 
132  if (thePhyDef.isNonnull()) {
133  int particPDG = thePhyDef.get()->pdgId();
134 
135  if (particPDG == 5) { //b from t
136  double bTcorr = bTopCorrector->correction(theJet);
137  double tTcorr = tTopCorrector->correction(theJet);
138  jetsNoCorr[0].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
139  jetsCorFl0[0].SetPtEtaPhiM(theJet.pt() * bTcorr, theJet.eta(), theJet.phi(), 0);
140  jetsCorFlM[0].SetPtEtaPhiM(theJet.pt() * bTcorr, theJet.eta(), theJet.phi(), bMass);
141  jetsCorMix[0].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
142  isQuarkFound[0] = true;
143  } else if (particPDG == -5) { //bbar from tbar
144  double bTcorr = bTopCorrector->correction(theJet);
145  double tTcorr = tTopCorrector->correction(theJet);
146  jetsNoCorr[3].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
147  jetsCorFl0[3].SetPtEtaPhiM(theJet.pt() * bTcorr, theJet.eta(), theJet.phi(), 0);
148  jetsCorFlM[3].SetPtEtaPhiM(theJet.pt() * bTcorr, theJet.eta(), theJet.phi(), bMass);
149  jetsCorMix[3].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
150  isQuarkFound[3] = true;
151  } else if (particPDG == 2) { //W+ from t
152  double qTcorr = qTopCorrector->correction(theJet);
153  double tTcorr = tTopCorrector->correction(theJet);
154  jetsNoCorr[1].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
155  jetsCorFl0[1].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), 0);
156  jetsCorFlM[1].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), qMass);
157  jetsCorMix[1].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
158  isQuarkFound[1] = true;
159  } else if (particPDG == 4) { //W+ from t
160  double cTcorr = cTopCorrector->correction(theJet);
161  double tTcorr = tTopCorrector->correction(theJet);
162  jetsNoCorr[1].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
163  jetsCorFl0[1].SetPtEtaPhiM(theJet.pt() * cTcorr, theJet.eta(), theJet.phi(), 0);
164  jetsCorFlM[1].SetPtEtaPhiM(theJet.pt() * cTcorr, theJet.eta(), theJet.phi(), cMass);
165  jetsCorMix[1].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
166  isQuarkFound[1] = true;
167  } else if (particPDG == -1) { //W+ from t
168  double qTcorr = qTopCorrector->correction(theJet);
169  double tTcorr = tTopCorrector->correction(theJet);
170  jetsNoCorr[2].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
171  jetsCorFl0[2].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), 0);
172  jetsCorFlM[2].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), qMass);
173  jetsCorMix[2].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
174  isQuarkFound[2] = true;
175  } else if (particPDG == -3) { //W+ from t
176  double qTcorr = qTopCorrector->correction(theJet);
177  double tTcorr = tTopCorrector->correction(theJet);
178  jetsNoCorr[2].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
179  jetsCorFl0[2].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), 0);
180  jetsCorFlM[2].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), qMass);
181  jetsCorMix[2].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
182  isQuarkFound[2] = true;
183  } else if (particPDG == -2) { //W- from tbar
184  double qTcorr = qTopCorrector->correction(theJet);
185  double tTcorr = tTopCorrector->correction(theJet);
186  jetsNoCorr[4].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
187  jetsCorFl0[4].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), 0);
188  jetsCorFlM[4].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), qMass);
189  jetsCorMix[4].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
190  isQuarkFound[4] = true;
191  } else if (particPDG == -4) { //W- from tbar
192  double qTcorr = qTopCorrector->correction(theJet);
193  double tTcorr = tTopCorrector->correction(theJet);
194  jetsNoCorr[4].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
195  jetsCorFl0[4].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), 0);
196  jetsCorFlM[4].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), qMass);
197  jetsCorMix[4].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
198  isQuarkFound[4] = true;
199  } else if (particPDG == 1) { //W- from tbar
200  double qTcorr = qTopCorrector->correction(theJet);
201  double tTcorr = tTopCorrector->correction(theJet);
202  jetsNoCorr[5].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
203  jetsCorFl0[5].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), 0);
204  jetsCorFlM[5].SetPtEtaPhiM(theJet.pt() * qTcorr, theJet.eta(), theJet.phi(), qMass);
205  jetsCorMix[5].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
206  isQuarkFound[5] = true;
207  } else if (particPDG == 3) { //W- from tbar
208  double cTcorr = cTopCorrector->correction(theJet);
209  double tTcorr = tTopCorrector->correction(theJet);
210  jetsNoCorr[5].SetPtEtaPhiM(theJet.pt(), theJet.eta(), theJet.phi(), theJet.mass());
211  jetsCorFl0[5].SetPtEtaPhiM(theJet.pt() * cTcorr, theJet.eta(), theJet.phi(), 0);
212  jetsCorFlM[5].SetPtEtaPhiM(theJet.pt() * cTcorr, theJet.eta(), theJet.phi(), cMass);
213  jetsCorMix[5].SetPtEtaPhiM(theJet.pt() * tTcorr, theJet.eta(), theJet.phi(), 0);
214  isQuarkFound[5] = true;
215  }
216  }
217  }
218 
219  if (isQuarkFound[0] && isQuarkFound[1] && isQuarkFound[2]) {
220  TLorentzVector *theTopPNoCorr = new TLorentzVector(jetsNoCorr[0] + jetsNoCorr[1] + jetsNoCorr[2]);
221  TLorentzVector *theTopPCorFl0 = new TLorentzVector(jetsCorFl0[0] + jetsCorFl0[1] + jetsCorFl0[2]);
222  TLorentzVector *theTopPCorFlM = new TLorentzVector(jetsCorFlM[0] + jetsCorFlM[1] + jetsCorFlM[2]);
223  TLorentzVector *theTopPCorMix = new TLorentzVector(jetsCorMix[0] + jetsCorMix[1] + jetsCorMix[2]);
224  hMassNoCorr->Fill(theTopPNoCorr->M());
225  hMassCorFl0->Fill(theTopPCorFl0->M());
226  hMassCorFlM->Fill(theTopPCorFlM->M());
227  hMassCorMix->Fill(theTopPCorMix->M());
228  }
229 
230  if (isQuarkFound[3] && isQuarkFound[4] && isQuarkFound[5]) {
231  TLorentzVector *theTopMNoCorr = new TLorentzVector(jetsNoCorr[3] + jetsNoCorr[4] + jetsNoCorr[5]);
232  TLorentzVector *theTopMCorFl0 = new TLorentzVector(jetsCorFl0[3] + jetsCorFl0[4] + jetsCorFl0[5]);
233  TLorentzVector *theTopMCorFlM = new TLorentzVector(jetsCorFlM[3] + jetsCorFlM[4] + jetsCorFlM[5]);
234  TLorentzVector *theTopMCorMix = new TLorentzVector(jetsCorMix[3] + jetsCorMix[4] + jetsCorMix[5]);
235  hMassNoCorr->Fill(theTopMNoCorr->M());
236  hMassCorFl0->Fill(theTopMCorFl0->M());
237  hMassCorFlM->Fill(theTopMCorFlM->M());
238  hMassCorMix->Fill(theTopMCorMix->M());
239  }
240 }
edm::Handle< reco::JetMatchedPartonsCollection > theJetPartonMatch
std::string m_cT_CorrectorName
edm::InputTag sourceByRefer_
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
const_iterator begin() const
double correction(const LorentzVector &fJet) const
get correction using Jet information only
Definition: JetCorrector.h:46
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
int iEvent
Definition: GenABIO.cc:224
const GenParticleRef & physicsDefinitionParton() const
const_iterator end() const
std::string m_qT_CorrectorName
static const JetCorrector * getJetCorrector(const std::string &fName, const edm::EventSetup &fSetup)
retrieve corrector from the event setup. troughs exception if something is missing ...
Definition: JetCorrector.cc:48
std::string m_tT_CorrectorName
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:232
std::string m_bT_CorrectorName

Member Data Documentation

◆ bMass

float calcTopMass::bMass
private

Definition at line 71 of file CalcTopMassExample.cc.

◆ cMass

float calcTopMass::cMass
private

Definition at line 72 of file CalcTopMassExample.cc.

◆ hMassCorFl0

TH1F* calcTopMass::hMassCorFl0
private

Definition at line 76 of file CalcTopMassExample.cc.

◆ hMassCorFlM

TH1F* calcTopMass::hMassCorFlM
private

Definition at line 77 of file CalcTopMassExample.cc.

◆ hMassCorMix

TH1F* calcTopMass::hMassCorMix
private

Definition at line 78 of file CalcTopMassExample.cc.

◆ hMassNoCorr

TH1F* calcTopMass::hMassNoCorr
private

Definition at line 75 of file CalcTopMassExample.cc.

◆ m_bT_CorrectorName

std::string calcTopMass::m_bT_CorrectorName
private

Definition at line 68 of file CalcTopMassExample.cc.

◆ m_cT_CorrectorName

std::string calcTopMass::m_cT_CorrectorName
private

Definition at line 67 of file CalcTopMassExample.cc.

◆ m_qT_CorrectorName

std::string calcTopMass::m_qT_CorrectorName
private

Definition at line 66 of file CalcTopMassExample.cc.

◆ m_tT_CorrectorName

std::string calcTopMass::m_tT_CorrectorName
private

Definition at line 69 of file CalcTopMassExample.cc.

◆ qMass

float calcTopMass::qMass
private

Definition at line 73 of file CalcTopMassExample.cc.

◆ sourceByRefer_

edm::InputTag calcTopMass::sourceByRefer_
private

Definition at line 63 of file CalcTopMassExample.cc.

◆ sourcePartons_

edm::InputTag calcTopMass::sourcePartons_
private

Definition at line 62 of file CalcTopMassExample.cc.

◆ theJetPartonMatch

edm::Handle<reco::JetMatchedPartonsCollection> calcTopMass::theJetPartonMatch
private

Definition at line 64 of file CalcTopMassExample.cc.