CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Types | Private Member Functions | Private Attributes
JetAnaPythia< Jet > Class Template Reference

#include <JetAnaPythia.h>

Inheritance diagram for JetAnaPythia< Jet >:
edm::EDAnalyzer edm::EDConsumerBase

Public Member Functions

 JetAnaPythia (edm::ParameterSet const &cfg)
 
- Public Member Functions inherited from edm::EDAnalyzer
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzer ()
 
ModuleDescription const & moduleDescription () const
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
 EDConsumerBase ()
 
ProductHolderIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
std::vector
< ProductHolderIndexAndSkipBit >
const & 
itemsToGetFromEvent () const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesDependentUpon (std::string const &iProcessName, std::string const &iModuleLabel, bool iPrint, std::vector< char const * > &oModuleLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
bool registeredToConsume (ProductHolderIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)
 
virtual ~EDConsumerBase ()
 

Private Types

typedef std::vector< Jet > JetCollection
 

Private Member Functions

void analyze (edm::Event const &e, edm::EventSetup const &iSetup)
 
void beginJob ()
 
void endJob ()
 
void FillHist1D (const TString &histName, const Double_t &x, const Double_t &wt)
 

Private Attributes

std::string anaLevel
 
bool debug
 
float diJetMass
 
float diPartMass
 
float etaJet1
 
float etaJet2
 
float etaPart1
 
float etaPart2
 
int eventsGen
 
std::string HistoFileName
 
std::string JetAlgorithm
 
TFile * m_file
 
std::map< TString, TH1 * > m_HistNames1D
 
TTree * mcTruthTree_
 
int nJets
 
int NJets
 
float pt_hat
 
std::vector< double > ptHatEdges
 
float ptJet1
 
float ptJet2
 
float ptPart1
 
float ptPart2
 
float weight
 
float xsec
 
std::vector< double > xsecGen
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
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 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)
 

Detailed Description

template<class Jet>
class JetAnaPythia< Jet >

Definition at line 16 of file JetAnaPythia.h.

Member Typedef Documentation

template<class Jet >
typedef std::vector<Jet> JetAnaPythia< Jet >::JetCollection
private

Definition at line 21 of file JetAnaPythia.h.

Constructor & Destructor Documentation

template<class Jet >
JetAnaPythia< Jet >::JetAnaPythia ( edm::ParameterSet const &  cfg)

Definition at line 27 of file JetAnaPythia.cc.

References debug, edm::ParameterSet::getParameter(), and AlCaHLTBitMon_QueryRunRegistry::string.

28 {
29  JetAlgorithm = cfg.getParameter<std::string> ("JetAlgorithm");
30  HistoFileName = cfg.getParameter<std::string> ("HistoFileName");
31  NJets = cfg.getParameter<int> ("NJets");
32  debug = cfg.getParameter<bool> ("debug");
33  eventsGen = cfg.getParameter<int> ("eventsGen");
34  anaLevel = cfg.getParameter<std::string> ("anaLevel");
35  xsecGen = cfg.getParameter< vector<double> > ("xsecGen");
36  ptHatEdges = cfg.getParameter< vector<double> > ("ptHatEdges");
37 
38 }
std::vector< double > xsecGen
Definition: JetAnaPythia.h:60
tuple cfg
Definition: looper.py:293
std::string JetAlgorithm
Definition: JetAnaPythia.h:42
std::string HistoFileName
Definition: JetAnaPythia.h:44
std::vector< double > ptHatEdges
Definition: JetAnaPythia.h:61
std::string anaLevel
Definition: JetAnaPythia.h:56

Member Function Documentation

template<class Jet >
void JetAnaPythia< Jet >::analyze ( edm::Event const &  e,
edm::EventSetup const &  iSetup 
)
privatevirtual

evt.getRun().getByLabel("generator", genInfoProduct );

Histograms for Dijet Mass Analysis ////

Histograms for Dijet Ratio Analysis: Inner region ///

Histograms for Dijet Ratio Analysis: Outer region ////

We are looking at dijet resonances. ////

We are looking at QCD ////

Diparton mass for dijet mass analysis ////

Diparton mass for dijet ratio analysis: inner region ///

Diparton mass for dijet ratio analysis: outer region ///

Implements edm::EDAnalyzer.

Definition at line 140 of file JetAnaPythia.cc.

References funct::abs(), gather_cfg::cout, debug, edm::Event::getByLabel(), i, cmsHarvester::index, fwrapper::jets, ResonanceBuilder::mass, AlCaHLTBitMon_ParallelJobs::p, reco::LeafCandidate::p4(), reco::LeafCandidate::pdgId(), lumiQueryAPI::q, alignCSCRings::r, reco::LeafCandidate::status(), and histoStyle::weight.

141 {
142  int notDone=1;
143  while(notDone)
144  { //while loop to allow us to tailor the analysis level for faster running.
145  TString hname;
146 
147  // Process Info
148 
149  //edm::Handle< double > genEventScale;
150  //evt.getByLabel("genEventScale", genEventScale );
151  //pt_hat = *genEventScale;
152 
154  evt.getByLabel("generatorSmeared", MCevt);
155  HepMC::GenEvent * myGenEvent = new HepMC::GenEvent(*(MCevt->GetEvent()));
156 
157  double pthat = myGenEvent->event_scale();
158  pt_hat = float(pthat);
159 
160  delete myGenEvent;
161 
162  if(anaLevel != "generating"){ //We are not generating events, so xsec is there
163  //edm::Handle< GenRunInfoProduct > genInfoProduct;
165  //xsec = (double)genInfoProduct->externalXSecLO();
166  xsec=0.0;
167  if( ptHatEdges.size()>xsecGen.size() ){
168  for(unsigned int i_pthat = 0; i_pthat < xsecGen.size(); ++i_pthat){
169  if( pthat >= ptHatEdges[i_pthat] && pthat < ptHatEdges[i_pthat+1])xsec=float(xsecGen[i_pthat]);
170  }
171  }
172  else
173  {
174  std::cout << "Number of PtHat bin edges too small. Xsec set to zero" << std::endl;
175  }
176  }
177  else
178  {
179  xsec = xsecGen[0]; //Generating events, no xsec in event, get xsec from user input
180  }
181  if(debug)std::cout << "cross section=" <<xsec << " pb" << std::endl;
183 
184  if(debug)std::cout << "pt_hat=" <<pt_hat << std::endl;
185  hname = "PtHat";
186  FillHist1D(hname, pt_hat, 1.0);
187  hname = "PtHatFine";
188  FillHist1D(hname, pt_hat, 1.0);
189  hname = "PtHatWt";
190  FillHist1D(hname, pt_hat, weight);
191  hname = "PtHatFineWt";
192  FillHist1D(hname, pt_hat, weight);
193  if(anaLevel=="PtHatOnly")break; //ptHatOnly should be very fast
194 
195  // Jet Info
196  math::XYZTLorentzVector p4jet[2];
197  float etajet[2];
200  evt.getByLabel(JetAlgorithm,jets);
201  typename JetCollection::const_iterator i_jet;
202  int index = 0;
203 
205  hname = "NumberOfJets";
206  nJets = jets->size();
207  FillHist1D(hname,nJets,1.0);
208 
209 
210  // Two Leading Jet Info
211  for(i_jet = jets->begin(); i_jet != jets->end() && index < 2; ++i_jet)
212  {
213  hname = "JetPt";
214  FillHist1D(hname,i_jet->pt(),1.0);
215  hname = "JetEta";
216  FillHist1D(hname,i_jet->eta(),1.0);
217  hname = "JetPhi";
218  FillHist1D(hname,i_jet->phi(),1.0);
219  p4jet[index] = i_jet->p4();
220  etajet[index] = i_jet->eta();
221  if(debug)std::cout << "jet " << index+1 <<": pt=" <<i_jet->pt() << ", eta=" <<etajet[index] << std::endl;
222  index++;
223  }
224 
225  // TTree variables //
226  etaJet1 = etajet[0];
227  etaJet2 = etajet[1];
228  ptJet1 = p4jet[0].pt();
229  ptJet2 = p4jet[1].pt();
230  diJetMass = (p4jet[0]+p4jet[1]).mass();
231 
233  if(index==2&&abs(etaJet1)<1.3&&abs(etaJet2)<1.3){
234  hname = "DijetMass";
235  FillHist1D(hname,diJetMass ,1.0);
236  hname = "DijetMassWt";
237  FillHist1D(hname,diJetMass ,weight);
238  }
239 
241  if(index==2&&abs(etaJet1)<0.7&&abs(etaJet2)<0.7){
242  hname = "DijetMassIn";
243  FillHist1D(hname,diJetMass ,1.0);
244  hname = "DijetMassInWt";
245  FillHist1D(hname,diJetMass ,weight);
246  }
248  if(index==2 && (abs(etaJet1)>0.7&&abs(etaJet1)<1.3)
249  && (abs(etaJet2)>0.7&&abs(etaJet2)<1.3) ){
250  hname = "DijetMassOut";
251  FillHist1D(hname, diJetMass ,1.0);
252  hname = "DijetMassOutWt";
253  FillHist1D(hname,diJetMass ,weight);
254  }
255  if(anaLevel=="Jets")break; //Jets level for samples without genParticles
256 
257 
258  // Parton Info
259  edm::Handle<std::vector<reco::GenParticle> > genParticlesHandle_;
260  evt.getByLabel("genParticles",genParticlesHandle_);
261  if(debug)for( size_t i = 0; i < genParticlesHandle_->size(); ++ i ) {
262  const reco::GenParticle & p = (*genParticlesHandle_)[i];
263  int id = p.pdgId();
264  int st = p.status();
265  math::XYZTLorentzVector genP4 = p.p4();
266  if(i>=2&&i<=8)std::cout << "particle " << i << ": id=" << id << ", status=" << st << ", mass=" << genP4.mass() << ", pt=" << genP4.pt() << ", eta=" << genP4.eta() << std::endl;
267  }
268  // Examine the 7th particle in pythia.
269  // It should be either a resonance (abs(id)>=32) or the first outgoing parton
270  // for the processes we will consider: dijet resonances, QCD, or QCD +contact interactions.
271  const reco::GenParticle & p = (*genParticlesHandle_)[6];
272  int id = p.pdgId();
273  math::XYZTLorentzVector resonance_p, parton1_p, parton2_p;
274  if(abs(id)>=32){
276  resonance_p = p.p4();
277  hname = "ResonanceMass";
278  FillHist1D(hname,resonance_p.mass() ,1.0);
279  const reco::GenParticle & q = (*genParticlesHandle_)[7];
280  parton1_p = q.p4();
281  const reco::GenParticle & r = (*genParticlesHandle_)[8];
282  parton2_p = r.p4();
283  if(debug)std::cout << "Resonance mass=" << resonance_p.mass() << ", parton 1 pt=" << parton1_p.pt() << ", parton 2 pt=" << parton2_p.pt() << ", diparton mass=" << (parton1_p+parton2_p).mass() << std::endl;
284  }
285  else
286  {
288  parton1_p = p.p4();
289  const reco::GenParticle & q = (*genParticlesHandle_)[7];
290  parton2_p = q.p4();
291  if(debug)std::cout << "parton 1 pt=" << parton1_p.pt() << ", parton 2 pt=" << parton2_p.pt() << ", diparton mass=" << (parton1_p+parton2_p).mass() << std::endl;
292  }
293 
294  etaPart1 = parton1_p.eta();
295  etaPart2 = parton2_p.eta();
296  ptPart1 = parton1_p.pt();
297  ptPart2 = parton2_p.pt();
298  diPartMass = (parton1_p+parton2_p).mass();
300  if(abs(etaPart1)<1.3&&abs(etaPart2)<1.3){
301  hname = "DipartonMass";
302  FillHist1D(hname,diPartMass ,1.0);
303  hname = "DipartonMassWt";
304  FillHist1D(hname,diPartMass ,weight);
305  }
307  if(abs(etaPart1)<0.7&&abs(etaPart2)<0.7){
308  hname = "DipartonMassIn";
309  FillHist1D(hname,diPartMass ,1.0);
310  hname = "DipartonMassInWt";
311  FillHist1D(hname,diPartMass ,weight);
312  }
314  if( (abs(etaPart1)>0.7&&abs(etaPart1)<1.3)
315  && (abs(etaPart2)>0.7&&abs(etaPart2)<1.3) ){
316  hname = "DipartonMassOut";
317  FillHist1D(hname,diPartMass ,1.0);
318  hname = "DipartonMassOutWt";
319  FillHist1D(hname,diPartMass ,weight);
320  }
321 
322  // Fill the TTree //
323  mcTruthTree_->Fill();
324 
325  notDone=0; //We are done, exit the while loop
326  }//end of while
327 
328 }
int i
Definition: DBlmapReader.cc:9
std::vector< double > xsecGen
Definition: JetAnaPythia.h:60
float diJetMass
Definition: JetAnaPythia.h:37
void FillHist1D(const TString &histName, const Double_t &x, const Double_t &wt)
float etaPart2
Definition: JetAnaPythia.h:35
std::string JetAlgorithm
Definition: JetAnaPythia.h:42
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
virtual int status() const final
status word
float etaPart1
Definition: JetAnaPythia.h:35
float diPartMass
Definition: JetAnaPythia.h:38
vector< PseudoJet > jets
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
TTree * mcTruthTree_
Definition: JetAnaPythia.h:28
std::vector< double > ptHatEdges
Definition: JetAnaPythia.h:61
std::string anaLevel
Definition: JetAnaPythia.h:56
virtual int pdgId() const final
PDG identifier.
tuple cout
Definition: gather_cfg.py:145
virtual const LorentzVector & p4() const final
four-momentum Lorentz vector
Definition: LeafCandidate.h:99
template<class Jet >
void JetAnaPythia< Jet >::beginJob ( void  )
privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 41 of file JetAnaPythia.cc.

References M_PI, and histoStyle::weight.

42 {
43  TString hname;
44  m_file = new TFile(HistoFileName.c_str(),"RECREATE");
46  const int nMassBins = 103;
47  double massBoundaries[nMassBins+1] = {1, 3, 6, 10, 16, 23, 31, 40, 50, 61, 74, 88, 103, 119, 137, 156, 176, 197, 220, 244, 270, 296, 325, 354, 386, 419, 453, 489, 526, 565, 606, 649, 693, 740, 788, 838, 890, 944, 1000, 1058, 1118, 1181, 1246, 1313, 1383, 1455, 1530, 1607, 1687, 1770, 1856, 1945, 2037, 2132, 2231, 2332, 2438, 2546, 2659, 2775, 2895, 3019, 3147, 3279, 3416, 3558, 3704, 3854, 4010, 4171, 4337, 4509, 4686, 4869, 5058, 5253, 5455, 5663, 5877, 6099, 6328, 6564, 6808, 7060, 7320, 7589, 7866, 8152, 8447, 8752, 9067, 9391, 9726, 10072, 10430, 10798, 11179, 11571, 11977, 12395, 12827, 13272, 13732, 14000};
48 
49  hname = "JetPt";
50  m_HistNames1D[hname] = new TH1F(hname,hname,500,0,5000);
51 
52  hname = "JetEta";
53  m_HistNames1D[hname] = new TH1F(hname,hname,120,-6,6);
54 
55  hname = "JetPhi";
56  m_HistNames1D[hname] = new TH1F(hname,hname,100,-M_PI,M_PI);
57 
58  hname = "NumberOfJets";
59  m_HistNames1D[hname] = new TH1F(hname,hname,100,0,100);
60 
61  hname = "DijetMass";
62  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
63 
64  hname = "DijetMassWt";
65  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
66  m_HistNames1D.find(hname)->second->Sumw2();
67 
68  hname = "DijetMassIn";
69  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
70 
71  hname = "DijetMassInWt";
72  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
73  m_HistNames1D.find(hname)->second->Sumw2();
74 
75  hname = "DijetMassOut";
76  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
77 
78  hname = "DijetMassOutWt";
79  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
80  m_HistNames1D.find(hname)->second->Sumw2();
81 
82  hname = "ResonanceMass";
83  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
84 
85  hname = "DipartonMass";
86  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
87 
88  hname = "DipartonMassWt";
89  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
90  m_HistNames1D.find(hname)->second->Sumw2();
91 
92  hname = "DipartonMassIn";
93  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
94 
95  hname = "DipartonMassInWt";
96  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
97  m_HistNames1D.find(hname)->second->Sumw2();
98 
99  hname = "DipartonMassOut";
100  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
101 
102  hname = "DipartonMassOutWt";
103  m_HistNames1D[hname] = new TH1F(hname,hname,nMassBins,massBoundaries);
104  m_HistNames1D.find(hname)->second->Sumw2();
105 
106  hname = "PtHat";
107  m_HistNames1D[hname] = new TH1F(hname,hname,1000,0,5000);
108 
109  hname = "PtHatWt";
110  m_HistNames1D[hname] = new TH1F(hname,hname,1000,0,5000);
111  m_HistNames1D.find(hname)->second->Sumw2();
112 
113  hname = "PtHatFine";
114  m_HistNames1D[hname] = new TH1F(hname,hname,5000,0,5000);
115 
116  hname = "PtHatFineWt";
117  m_HistNames1D[hname] = new TH1F(hname,hname,5000,0,5000);
118  m_HistNames1D.find(hname)->second->Sumw2();
119 
120  mcTruthTree_ = new TTree("mcTruthTree","mcTruthTree");
121  mcTruthTree_->Branch("xsec", &xsec, "xsec/F");
122  mcTruthTree_->Branch("weight", &weight, "weight/F");
123  mcTruthTree_->Branch("pt_hat", &pt_hat, "pt_hat/F");
124  mcTruthTree_->Branch("nJets", &nJets, "nJets/I");
125  mcTruthTree_->Branch("etaJet1", &etaJet1, "etaJet1/F");
126  mcTruthTree_->Branch("etaJet2", &etaJet2, "etaJet2/F");
127  mcTruthTree_->Branch("ptJet1", &ptJet1, "ptJet1/F");
128  mcTruthTree_->Branch("ptJet2", &ptJet2, "ptJet2/F");
129  mcTruthTree_->Branch("diJetMass", &diJetMass, "diJetMass/F");
130  mcTruthTree_->Branch("etaPart1", &etaPart1, "etaPart1/F");
131  mcTruthTree_->Branch("etaPart2", &etaPart2, "etaPart2/F");
132  mcTruthTree_->Branch("ptPart1", &ptPart1, "ptPart1/F");
133  mcTruthTree_->Branch("ptPart2", &ptPart2, "ptPart2/F");
134  mcTruthTree_->Branch("diPartMass", &diPartMass, "diPartMass/F");
135 
136 
137 }
float diJetMass
Definition: JetAnaPythia.h:37
float etaPart2
Definition: JetAnaPythia.h:35
std::string HistoFileName
Definition: JetAnaPythia.h:44
float etaPart1
Definition: JetAnaPythia.h:35
std::map< TString, TH1 * > m_HistNames1D
Definition: JetAnaPythia.h:26
TFile * m_file
Definition: JetAnaPythia.h:39
float diPartMass
Definition: JetAnaPythia.h:38
TTree * mcTruthTree_
Definition: JetAnaPythia.h:28
#define M_PI
template<class Jet >
void JetAnaPythia< Jet >::endJob ( void  )
privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 331 of file JetAnaPythia.cc.

332 {
334  if (m_file !=0)
335  {
336  m_file->cd();
337  mcTruthTree_->Write();
338  for (std::map<TString, TH1*>::iterator hid = m_HistNames1D.begin(); hid != m_HistNames1D.end(); hid++)
339  hid->second->Write();
340  delete m_file;
341  m_file = 0;
342  }
343 }
std::map< TString, TH1 * > m_HistNames1D
Definition: JetAnaPythia.h:26
TFile * m_file
Definition: JetAnaPythia.h:39
TTree * mcTruthTree_
Definition: JetAnaPythia.h:28
template<class Jet >
void JetAnaPythia< Jet >::FillHist1D ( const TString &  histName,
const Double_t &  x,
const Double_t &  wt 
)
private

Definition at line 346 of file JetAnaPythia.cc.

References gather_cfg::cout.

347 {
348  std::map<TString, TH1*>::iterator hid=m_HistNames1D.find(histName);
349  if (hid==m_HistNames1D.end())
350  std::cout << "%fillHist -- Could not find histogram with name: " << histName << std::endl;
351  else
352  hid->second->Fill(value,wt);
353 }
std::map< TString, TH1 * > m_HistNames1D
Definition: JetAnaPythia.h:26
tuple cout
Definition: gather_cfg.py:145

Member Data Documentation

template<class Jet >
std::string JetAnaPythia< Jet >::anaLevel
private

Analysis level string. Can speed up job by looking at less /// PtHatOnly: only get PtHat and make PtHat histos Jets: do histogram analysis of jets, but not partons all: do analysis of everything and make histos and root tree generating: analysis of everything, make histos and root tree

Definition at line 56 of file JetAnaPythia.h.

template<class Jet >
bool JetAnaPythia< Jet >::debug
private

Definition at line 48 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::diJetMass
private

Definition at line 37 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::diPartMass
private

Definition at line 38 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::etaJet1
private

Definition at line 33 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::etaJet2
private

Definition at line 33 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::etaPart1
private

Definition at line 35 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::etaPart2
private

Definition at line 35 of file JetAnaPythia.h.

template<class Jet >
int JetAnaPythia< Jet >::eventsGen
private

Definition at line 50 of file JetAnaPythia.h.

template<class Jet >
std::string JetAnaPythia< Jet >::HistoFileName
private

Definition at line 44 of file JetAnaPythia.h.

template<class Jet >
std::string JetAnaPythia< Jet >::JetAlgorithm
private

Definition at line 42 of file JetAnaPythia.h.

template<class Jet >
TFile* JetAnaPythia< Jet >::m_file
private

Definition at line 39 of file JetAnaPythia.h.

template<class Jet >
std::map<TString, TH1*> JetAnaPythia< Jet >::m_HistNames1D
private

Definition at line 26 of file JetAnaPythia.h.

template<class Jet >
TTree* JetAnaPythia< Jet >::mcTruthTree_
private

Definition at line 28 of file JetAnaPythia.h.

template<class Jet >
int JetAnaPythia< Jet >::nJets
private

Definition at line 32 of file JetAnaPythia.h.

template<class Jet >
int JetAnaPythia< Jet >::NJets
private

Definition at line 46 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::pt_hat
private

Definition at line 31 of file JetAnaPythia.h.

template<class Jet >
std::vector<double> JetAnaPythia< Jet >::ptHatEdges
private

Definition at line 61 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::ptJet1
private

Definition at line 34 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::ptJet2
private

Definition at line 34 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::ptPart1
private

Definition at line 36 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::ptPart2
private

Definition at line 36 of file JetAnaPythia.h.

template<class Jet >
float JetAnaPythia< Jet >::weight
private
template<class Jet >
float JetAnaPythia< Jet >::xsec
private

Definition at line 29 of file JetAnaPythia.h.

template<class Jet >
std::vector<double> JetAnaPythia< Jet >::xsecGen
private

Generator cross section Only 1 entry in case analysis level is "generating" //// Multiple entries when analyzing ///

Definition at line 60 of file JetAnaPythia.h.