CMS 3D CMS Logo

AnalysisRootpleProducerOnlyMC Class Reference

#include <QCDAnalysis/UEAnalysis/interface/AnalysisRootpleProducerOnlyMC.h>

Inheritance diagram for AnalysisRootpleProducerOnlyMC:

edm::EDAnalyzer

List of all members.

Public Member Functions

 AnalysisRootpleProducerOnlyMC (const edm::ParameterSet &)
virtual void analyze (const edm::Event &, const edm::EventSetup &)
virtual void beginJob (const edm::EventSetup &)
virtual void endJob ()
void fillChargedJet (float, float, float, float)
void fillEventInfo (int)
void fillInclusiveJet (float, float, float, float)
void fillMCParticles (float, float, float, float)
void store ()
virtual ~AnalysisRootpleProducerOnlyMC ()

Private Attributes

TTree * AnalysisTree
Handle< vector< GenParticle > > CandHandleMC
TClonesArray * ChargedJet
Handle< GenJetCollection > ChgGenJetsHandle
InputTag chgGenPartCollName
InputTag chgJetCollName
float EtaCJ [NCJMAX]
float EtaEHJ [NEHJMAX]
float EtaIJ [NIJMAX]
float EtaMC [NMCPMAX]
float EtaTJ [NTJMAX]
float EtaTK [NTKMAX]
int EventKind
Handle< HepMCProductEvtHandle
edm::Service< TFileServicefs
InputTag genJetCollName
Handle< GenJetCollection > GenJetsHandle
TClonesArray * InclusiveJet
InputTag mcEvent
float MomentumCJ [NCJMAX]
float MomentumEHJ [NEHJMAX]
float MomentumIJ [NIJMAX]
float MomentumMC [NMCPMAX]
float MomentumTJ [NTJMAX]
float MomentumTK [NTKMAX]
TClonesArray * MonteCarlo
int NumberCaloJet
int NumberChargedJet
int NumberInclusiveJet
int NumberMCParticles
int NumberTracks
int NumberTracksJet
float PhiCJ [NCJMAX]
float PhiEHJ [NEHJMAX]
float PhiIJ [NIJMAX]
float PhiMC [NMCPMAX]
float PhiTJ [NTJMAX]
float PhiTK [NTKMAX]
float piG
float TransverseMomentumCJ [NCJMAX]
float TransverseMomentumEHJ [NEHJMAX]
float TransverseMomentumIJ [NIJMAX]
float TransverseMomentumMC [NMCPMAX]
float TransverseMomentumTJ [NTJMAX]
float TransverseMomentumTK [NTKMAX]

Static Private Attributes

static const int NCJMAX = 10000
static const int NEHJMAX = 10000
static const int NIJMAX = 10000
static const int NMCPMAX = 10000
static const int NTJMAX = 10000
static const int NTKMAX = 10000


Detailed Description

Definition at line 35 of file AnalysisRootpleProducerOnlyMC.h.


Constructor & Destructor Documentation

AnalysisRootpleProducerOnlyMC::AnalysisRootpleProducerOnlyMC ( const edm::ParameterSet pset  )  [explicit]

Definition at line 68 of file AnalysisRootpleProducerOnlyMC.cc.

References chgGenPartCollName, chgJetCollName, genJetCollName, edm::ParameterSet::getUntrackedParameter(), mcEvent, NumberChargedJet, NumberInclusiveJet, NumberMCParticles, and piG.

00069 {
00070   mcEvent = pset.getUntrackedParameter<InputTag>("MCEvent",std::string(""));
00071   genJetCollName = pset.getUntrackedParameter<InputTag>("GenJetCollectionName",std::string(""));
00072   chgJetCollName = pset.getUntrackedParameter<InputTag>("ChgGenJetCollectionName",std::string(""));
00073   chgGenPartCollName = pset.getUntrackedParameter<InputTag>("ChgGenPartCollectionName",std::string(""));
00074 
00075   piG = acos(-1.);
00076   NumberMCParticles=0;
00077   NumberInclusiveJet=0;
00078   NumberChargedJet=0;
00079 }

virtual AnalysisRootpleProducerOnlyMC::~AnalysisRootpleProducerOnlyMC (  )  [inline, virtual]

Definition at line 41 of file AnalysisRootpleProducerOnlyMC.h.

00041 {} 


Member Function Documentation

void AnalysisRootpleProducerOnlyMC::analyze ( const edm::Event e,
const edm::EventSetup  
) [virtual]

Implements edm::EDAnalyzer.

Definition at line 129 of file AnalysisRootpleProducerOnlyMC.cc.

References CandHandleMC, ChargedJet, ChgGenJetsHandle, chgGenPartCollName, chgJetCollName, EventKind, EvtHandle, fillChargedJet(), fillInclusiveJet(), fillMCParticles(), genJetCollName, GenJetsHandle, edm::Event::getByLabel(), InclusiveJet, it, mcEvent, MonteCarlo, and store().

00130 {
00131 
00132   e.getByLabel( mcEvent           , EvtHandle        ) ;
00133   e.getByLabel( chgGenPartCollName, CandHandleMC     );
00134   e.getByLabel( chgJetCollName    , ChgGenJetsHandle );
00135   e.getByLabel( genJetCollName    , GenJetsHandle    );
00136   
00137   const HepMC::GenEvent* Evt = EvtHandle->GetEvent() ;
00138   
00139   EventKind = Evt->signal_process_id();
00140   
00141   std::vector<math::XYZTLorentzVector> GenPart;
00142   std::vector<GenJet> ChgGenJetContainer;
00143   std::vector<GenJet> GenJetContainer;
00144   
00145   GenPart.clear();
00146   ChgGenJetContainer.clear();
00147   GenJetContainer.clear();
00148   
00149   ChargedJet->Clear();
00150   InclusiveJet->Clear();
00151   MonteCarlo->Clear();
00152 
00153   if (ChgGenJetsHandle->size()){
00154 
00155     for ( GenJetCollection::const_iterator it(ChgGenJetsHandle->begin()), itEnd(ChgGenJetsHandle->end());
00156           it!=itEnd; ++it)
00157       {
00158         ChgGenJetContainer.push_back(*it);
00159       }
00160     
00161     std::stable_sort(ChgGenJetContainer.begin(),ChgGenJetContainer.end(),GenJetSort());
00162     
00163     std::vector<GenJet>::const_iterator it(ChgGenJetContainer.begin()), itEnd(ChgGenJetContainer.end());
00164     for ( int iChargedJet(0); it != itEnd; ++it, ++iChargedJet)
00165       {
00166         fillChargedJet(it->p(),it->pt(),it->eta(),it->phi());
00167         new((*ChargedJet)[iChargedJet]) TLorentzVector(it->px(), it->py(), it->pz(), it->energy());
00168       }
00169   }
00170   
00171   if (GenJetsHandle->size()){
00172     
00173     for ( GenJetCollection::const_iterator it(GenJetsHandle->begin()), itEnd(GenJetsHandle->end());
00174           it!=itEnd; ++it )
00175       {
00176         GenJetContainer.push_back(*it);
00177       }
00178 
00179     std::stable_sort(GenJetContainer.begin(),GenJetContainer.end(),GenJetSort());
00180 
00181     std::vector<GenJet>::const_iterator it(GenJetContainer.begin()), itEnd(GenJetContainer.end()); 
00182     for ( int iInclusiveJet(0); it != itEnd; ++it, ++iInclusiveJet)
00183       {
00184         fillInclusiveJet(it->p(),it->pt(),it->eta(),it->phi());
00185         new((*InclusiveJet)[iInclusiveJet]) TLorentzVector(it->px(), it->py(), it->pz(), it->energy());
00186       }
00187   }
00188   
00189   if (CandHandleMC->size()){
00190 
00191     for (vector<GenParticle>::const_iterator it(CandHandleMC->begin()), itEnd(CandHandleMC->end());
00192          it != itEnd;it++)
00193       {
00194         GenPart.push_back(it->p4());
00195       }
00196 
00197     std::stable_sort(GenPart.begin(),GenPart.end(),GreaterPt());
00198 
00199     std::vector<math::XYZTLorentzVector>::const_iterator it(GenPart.begin()), itEnd(GenPart.end());
00200     for( int iMonteCarlo(0); it != itEnd; ++it, ++iMonteCarlo )
00201       {
00202         fillMCParticles(it->P(),it->Pt(),it->Eta(),it->Phi());
00203         new((*MonteCarlo)[iMonteCarlo]) TLorentzVector(it->Px(), it->Py(), it->Pz(), it->E()); 
00204       }
00205   }
00206 
00207   store();
00208 }

void AnalysisRootpleProducerOnlyMC::beginJob ( const edm::EventSetup  )  [virtual]

Reimplemented from edm::EDAnalyzer.

Definition at line 82 of file AnalysisRootpleProducerOnlyMC.cc.

References AnalysisTree, ChargedJet, EtaCJ, EtaIJ, EtaMC, EventKind, fs, InclusiveJet, MomentumCJ, MomentumIJ, MomentumMC, MonteCarlo, NumberChargedJet, NumberInclusiveJet, NumberMCParticles, PhiCJ, PhiIJ, PhiMC, TransverseMomentumCJ, TransverseMomentumIJ, and TransverseMomentumMC.

00083 {
00084   // use TFileService for output to root file 
00085   AnalysisTree = fs->make<TTree>("AnalysisTree","MBUE Analysis Tree ");
00086  
00087   // process type
00088   AnalysisTree->Branch("EventKind",&EventKind,"EventKind/I");
00089   
00090   // store p, pt, eta, phi for particles and jets
00091 
00092   // GenParticles at hadron level
00093   AnalysisTree->Branch("NumberMCParticles",&NumberMCParticles,"NumberMCParticles/I");
00094   AnalysisTree->Branch("MomentumMC",MomentumMC,"MomentumMC[NumberMCParticles]/F");
00095   AnalysisTree->Branch("TransverseMomentumMC",TransverseMomentumMC,"TransverseMomentumMC[NumberMCParticles]/F");
00096   AnalysisTree->Branch("EtaMC",EtaMC,"EtaMC[NumberMCParticles]/F");
00097   AnalysisTree->Branch("PhiMC",PhiMC,"PhiMC[NumberMCParticles]/F");
00098 
00099   // GenJets
00100   AnalysisTree->Branch("NumberInclusiveJet",&NumberInclusiveJet,"NumberInclusiveJet/I");
00101   AnalysisTree->Branch("MomentumIJ",MomentumIJ,"MomentumIJ[NumberInclusiveJet]/F");
00102   AnalysisTree->Branch("TrasverseMomentumIJ",TransverseMomentumIJ,"TransverseMomentumIJ[NumberInclusiveJet]/F");
00103   AnalysisTree->Branch("EtaIJ",EtaIJ,"EtaIJ[NumberInclusiveJet]/F");
00104   AnalysisTree->Branch("PhiIJ",PhiIJ,"PhiIJ[NumberInclusiveJet]/F");
00105   
00106   // jets from charged GenParticles
00107   AnalysisTree->Branch("NumberChargedJet",&NumberChargedJet,"NumberChargedJet/I");
00108   AnalysisTree->Branch("MomentumCJ",MomentumCJ,"MomentumCJ[NumberChargedJet]/F");
00109   AnalysisTree->Branch("TrasverseMomentumCJ",TransverseMomentumCJ,"TransverseMomentumCJ[NumberChargedJet]/F");
00110   AnalysisTree->Branch("EtaCJ",EtaCJ,"EtaCJ[NumberChargedJet]/F");
00111   AnalysisTree->Branch("PhiCJ",PhiCJ,"PhiCJ[NumberChargedJet]/F");
00112   
00113 
00114   // alternative storage method:
00115   // save TClonesArrays of TLorentzVectors
00116   // i.e. store 4-vectors of particles and jets
00117 
00118   MonteCarlo = new TClonesArray("TLorentzVector", 10000);
00119   AnalysisTree->Branch("MonteCarlo", "TClonesArray", &MonteCarlo, 128000, 0);
00120   
00121   InclusiveJet = new TClonesArray("TLorentzVector", 10000);
00122   AnalysisTree->Branch("InclusiveJet", "TClonesArray", &InclusiveJet, 128000, 0);
00123 
00124   ChargedJet = new TClonesArray("TLorentzVector", 10000);
00125   AnalysisTree->Branch("ChargedJet", "TClonesArray", &ChargedJet, 128000, 0);
00126 }

void AnalysisRootpleProducerOnlyMC::endJob ( void   )  [virtual]

Reimplemented from edm::EDAnalyzer.

Definition at line 210 of file AnalysisRootpleProducerOnlyMC.cc.

00211 {
00212 }

void AnalysisRootpleProducerOnlyMC::fillChargedJet ( float  p,
float  pt,
float  eta,
float  phi 
)

Definition at line 59 of file AnalysisRootpleProducerOnlyMC.cc.

References EtaCJ, MomentumCJ, NumberChargedJet, PhiCJ, and TransverseMomentumCJ.

Referenced by analyze().

void AnalysisRootpleProducerOnlyMC::fillEventInfo ( int  e  ) 

Definition at line 36 of file AnalysisRootpleProducerOnlyMC.cc.

References EventKind.

00037 {
00038   EventKind = e;
00039 }

void AnalysisRootpleProducerOnlyMC::fillInclusiveJet ( float  p,
float  pt,
float  eta,
float  phi 
)

Definition at line 50 of file AnalysisRootpleProducerOnlyMC.cc.

References EtaIJ, MomentumIJ, NumberInclusiveJet, PhiIJ, and TransverseMomentumIJ.

Referenced by analyze().

void AnalysisRootpleProducerOnlyMC::fillMCParticles ( float  p,
float  pt,
float  eta,
float  phi 
)

Definition at line 41 of file AnalysisRootpleProducerOnlyMC.cc.

References EtaMC, MomentumMC, NumberMCParticles, PhiMC, and TransverseMomentumMC.

Referenced by analyze().

void AnalysisRootpleProducerOnlyMC::store (  ) 

Definition at line 27 of file AnalysisRootpleProducerOnlyMC.cc.

References AnalysisTree, NumberChargedJet, NumberInclusiveJet, and NumberMCParticles.

Referenced by analyze().

00028 {
00029   AnalysisTree->Fill();
00030 
00031   NumberMCParticles=0;
00032   NumberInclusiveJet=0;
00033   NumberChargedJet=0;
00034 }


Member Data Documentation

TTree* AnalysisRootpleProducerOnlyMC::AnalysisTree [private]

Definition at line 70 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and store().

Handle< vector<GenParticle> > AnalysisRootpleProducerOnlyMC::CandHandleMC [private]

Definition at line 61 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze().

TClonesArray* AnalysisRootpleProducerOnlyMC::ChargedJet [private]

Definition at line 90 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze(), and beginJob().

Handle< GenJetCollection > AnalysisRootpleProducerOnlyMC::ChgGenJetsHandle [private]

Definition at line 63 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze().

InputTag AnalysisRootpleProducerOnlyMC::chgGenPartCollName [private]

Definition at line 58 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), and analyze().

InputTag AnalysisRootpleProducerOnlyMC::chgJetCollName [private]

Definition at line 57 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), and analyze().

float AnalysisRootpleProducerOnlyMC::EtaCJ[NCJMAX] [private]

Definition at line 84 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillChargedJet().

float AnalysisRootpleProducerOnlyMC::EtaEHJ[NEHJMAX] [private]

Definition at line 86 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::EtaIJ[NIJMAX] [private]

Definition at line 83 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillInclusiveJet().

float AnalysisRootpleProducerOnlyMC::EtaMC[NMCPMAX] [private]

Definition at line 81 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillMCParticles().

float AnalysisRootpleProducerOnlyMC::EtaTJ[NTJMAX] [private]

Definition at line 85 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::EtaTK[NTKMAX] [private]

Definition at line 82 of file AnalysisRootpleProducerOnlyMC.h.

int AnalysisRootpleProducerOnlyMC::EventKind [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze(), beginJob(), and fillEventInfo().

Handle< HepMCProduct > AnalysisRootpleProducerOnlyMC::EvtHandle [private]

Definition at line 60 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze().

edm::Service<TFileService> AnalysisRootpleProducerOnlyMC::fs [private]

Definition at line 68 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob().

InputTag AnalysisRootpleProducerOnlyMC::genJetCollName [private]

Definition at line 56 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), and analyze().

Handle< GenJetCollection > AnalysisRootpleProducerOnlyMC::GenJetsHandle [private]

Definition at line 62 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze().

TClonesArray* AnalysisRootpleProducerOnlyMC::InclusiveJet [private]

Definition at line 89 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze(), and beginJob().

InputTag AnalysisRootpleProducerOnlyMC::mcEvent [private]

Definition at line 55 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), and analyze().

float AnalysisRootpleProducerOnlyMC::MomentumCJ[NCJMAX] [private]

Definition at line 84 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillChargedJet().

float AnalysisRootpleProducerOnlyMC::MomentumEHJ[NEHJMAX] [private]

Definition at line 86 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::MomentumIJ[NIJMAX] [private]

Definition at line 83 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillInclusiveJet().

float AnalysisRootpleProducerOnlyMC::MomentumMC[NMCPMAX] [private]

Definition at line 81 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillMCParticles().

float AnalysisRootpleProducerOnlyMC::MomentumTJ[NTJMAX] [private]

Definition at line 85 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::MomentumTK[NTKMAX] [private]

Definition at line 82 of file AnalysisRootpleProducerOnlyMC.h.

TClonesArray* AnalysisRootpleProducerOnlyMC::MonteCarlo [private]

Definition at line 88 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by analyze(), and beginJob().

const int AnalysisRootpleProducerOnlyMC::NCJMAX = 10000 [static, private]

Definition at line 75 of file AnalysisRootpleProducerOnlyMC.h.

const int AnalysisRootpleProducerOnlyMC::NEHJMAX = 10000 [static, private]

Definition at line 77 of file AnalysisRootpleProducerOnlyMC.h.

const int AnalysisRootpleProducerOnlyMC::NIJMAX = 10000 [static, private]

Definition at line 74 of file AnalysisRootpleProducerOnlyMC.h.

const int AnalysisRootpleProducerOnlyMC::NMCPMAX = 10000 [static, private]

Definition at line 72 of file AnalysisRootpleProducerOnlyMC.h.

const int AnalysisRootpleProducerOnlyMC::NTJMAX = 10000 [static, private]

Definition at line 76 of file AnalysisRootpleProducerOnlyMC.h.

const int AnalysisRootpleProducerOnlyMC::NTKMAX = 10000 [static, private]

Definition at line 73 of file AnalysisRootpleProducerOnlyMC.h.

int AnalysisRootpleProducerOnlyMC::NumberCaloJet [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

int AnalysisRootpleProducerOnlyMC::NumberChargedJet [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), beginJob(), fillChargedJet(), and store().

int AnalysisRootpleProducerOnlyMC::NumberInclusiveJet [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), beginJob(), fillInclusiveJet(), and store().

int AnalysisRootpleProducerOnlyMC::NumberMCParticles [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC(), beginJob(), fillMCParticles(), and store().

int AnalysisRootpleProducerOnlyMC::NumberTracks [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

int AnalysisRootpleProducerOnlyMC::NumberTracksJet [private]

Definition at line 79 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::PhiCJ[NCJMAX] [private]

Definition at line 84 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillChargedJet().

float AnalysisRootpleProducerOnlyMC::PhiEHJ[NEHJMAX] [private]

Definition at line 86 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::PhiIJ[NIJMAX] [private]

Definition at line 83 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillInclusiveJet().

float AnalysisRootpleProducerOnlyMC::PhiMC[NMCPMAX] [private]

Definition at line 81 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillMCParticles().

float AnalysisRootpleProducerOnlyMC::PhiTJ[NTJMAX] [private]

Definition at line 85 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::PhiTK[NTKMAX] [private]

Definition at line 82 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::piG [private]

Definition at line 66 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by AnalysisRootpleProducerOnlyMC().

float AnalysisRootpleProducerOnlyMC::TransverseMomentumCJ[NCJMAX] [private]

Definition at line 84 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillChargedJet().

float AnalysisRootpleProducerOnlyMC::TransverseMomentumEHJ[NEHJMAX] [private]

Definition at line 86 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::TransverseMomentumIJ[NIJMAX] [private]

Definition at line 83 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillInclusiveJet().

float AnalysisRootpleProducerOnlyMC::TransverseMomentumMC[NMCPMAX] [private]

Definition at line 81 of file AnalysisRootpleProducerOnlyMC.h.

Referenced by beginJob(), and fillMCParticles().

float AnalysisRootpleProducerOnlyMC::TransverseMomentumTJ[NTJMAX] [private]

Definition at line 85 of file AnalysisRootpleProducerOnlyMC.h.

float AnalysisRootpleProducerOnlyMC::TransverseMomentumTK[NTKMAX] [private]

Definition at line 82 of file AnalysisRootpleProducerOnlyMC.h.


The documentation for this class was generated from the following files:
Generated on Tue Jun 9 18:14:40 2009 for CMSSW by  doxygen 1.5.4