Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-08-23 03:24:51

0001 #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h"
0002 #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h"
0003 #include "Calibration/IsolatedParticles/interface/GenSimInfo.h"
0004 
0005 #include "DataFormats/DetId/interface/DetId.h"
0006 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0007 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0008 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0009 #include "DataFormats/GeometrySurface/interface/GloballyPositioned.h"
0010 
0011 #include "DataFormats/Math/interface/deltaPhi.h"
0012 #include "DataFormats/Math/interface/LorentzVector.h"
0013 
0014 #include "FWCore/Framework/interface/Frameworkfwd.h"
0015 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0016 #include "FWCore/Framework/interface/EventSetup.h"
0017 #include "FWCore/Framework/interface/Event.h"
0018 #include "FWCore/Framework/interface/MakerMacros.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0020 //TFile Service
0021 #include "FWCore/ServiceRegistry/interface/Service.h"
0022 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0023 
0024 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0025 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0026 #include "Geometry/Records/interface/CaloTopologyRecord.h"
0027 #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
0028 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0029 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0030 
0031 #include "MagneticField/Engine/interface/MagneticField.h"
0032 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0033 
0034 #include "RecoCaloTools/Navigation/interface/CaloNavigator.h"
0035 
0036 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0037 #include "SimGeneral/HepPDTRecord/interface/ParticleDataTable.h"
0038 // track associator
0039 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0040 
0041 // root objects
0042 #include "TROOT.h"
0043 #include "TSystem.h"
0044 #include "TFile.h"
0045 #include "TH1I.h"
0046 #include "TH2D.h"
0047 #include "TTree.h"
0048 
0049 #include <cmath>
0050 #include <iostream>
0051 #include <iomanip>
0052 #include <list>
0053 #include <vector>
0054 
0055 class StudyCaloGen : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0056 public:
0057   explicit StudyCaloGen(const edm::ParameterSet &);
0058 
0059   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0060 
0061 private:
0062   void beginJob() override;
0063   void analyze(const edm::Event &, const edm::EventSetup &) override;
0064   void endJob() override {}
0065 
0066   void fillTrack(
0067       GlobalPoint &posVec, math::XYZTLorentzVector &momVec, GlobalPoint &posECAL, int pdgId, bool okECAL, bool accpet);
0068   void fillIsolatedTrack(math::XYZTLorentzVector &momVec, GlobalPoint &posECAL, int pdgId);
0069   void bookHistograms();
0070   void clearTreeVectors();
0071   int particleCode(int);
0072 
0073   static constexpr int NPBins_ = 3;
0074   static constexpr int NEtaBins_ = 4;
0075   static constexpr int PBins_ = 32, EtaBins_ = 60, Particles = 12;
0076   int nEventProc;
0077   double genPartPBins_[NPBins_ + 1], genPartEtaBins_[NEtaBins_ + 1];
0078   double ptMin_, etaMax_, pCutIsolate_;
0079   bool a_Isolation_;
0080   std::string genSrc_;
0081 
0082   edm::EDGetTokenT<edm::HepMCProduct> tok_hepmc_;
0083   edm::EDGetTokenT<reco::GenParticleCollection> tok_genParticles_;
0084 
0085   edm::ESGetToken<CaloGeometry, CaloGeometryRecord> tok_geom_;
0086   edm::ESGetToken<CaloTopology, CaloTopologyRecord> tok_caloTopology_;
0087   edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> tok_topo_;
0088   edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> tok_magField_;
0089   edm::ESGetToken<HepPDT::ParticleDataTable, PDTRecord> tok_pdt_;
0090 
0091   bool useHepMC_;
0092   double a_coneR_, a_charIsoR_, a_neutIsoR_, a_mipR_;
0093   int verbosity_;
0094 
0095   TH1I *h_NEventProc;
0096   TH2D *h_pEta[Particles];
0097   TTree *tree_;
0098 
0099   std::vector<double> *t_isoTrkPAll;
0100   std::vector<double> *t_isoTrkPtAll;
0101   std::vector<double> *t_isoTrkPhiAll;
0102   std::vector<double> *t_isoTrkEtaAll;
0103   std::vector<double> *t_isoTrkPdgIdAll;
0104   std::vector<double> *t_isoTrkDEtaAll;
0105   std::vector<double> *t_isoTrkDPhiAll;
0106 
0107   std::vector<double> *t_isoTrkP;
0108   std::vector<double> *t_isoTrkPt;
0109   std::vector<double> *t_isoTrkEne;
0110   std::vector<double> *t_isoTrkEta;
0111   std::vector<double> *t_isoTrkPhi;
0112   std::vector<double> *t_isoTrkEtaEC;
0113   std::vector<double> *t_isoTrkPhiEC;
0114   std::vector<double> *t_isoTrkPdgId;
0115 
0116   std::vector<double> *t_maxNearP31x31;
0117   std::vector<double> *t_cHadronEne31x31, *t_cHadronEne31x31_1, *t_cHadronEne31x31_2, *t_cHadronEne31x31_3;
0118   std::vector<double> *t_nHadronEne31x31;
0119   std::vector<double> *t_photonEne31x31;
0120   std::vector<double> *t_eleEne31x31;
0121   std::vector<double> *t_muEne31x31;
0122 
0123   std::vector<double> *t_maxNearP25x25;
0124   std::vector<double> *t_cHadronEne25x25, *t_cHadronEne25x25_1, *t_cHadronEne25x25_2, *t_cHadronEne25x25_3;
0125   std::vector<double> *t_nHadronEne25x25;
0126   std::vector<double> *t_photonEne25x25;
0127   std::vector<double> *t_eleEne25x25;
0128   std::vector<double> *t_muEne25x25;
0129 
0130   std::vector<double> *t_maxNearP21x21;
0131   std::vector<double> *t_cHadronEne21x21, *t_cHadronEne21x21_1, *t_cHadronEne21x21_2, *t_cHadronEne21x21_3;
0132   std::vector<double> *t_nHadronEne21x21;
0133   std::vector<double> *t_photonEne21x21;
0134   std::vector<double> *t_eleEne21x21;
0135   std::vector<double> *t_muEne21x21;
0136 
0137   std::vector<double> *t_maxNearP15x15;
0138   std::vector<double> *t_cHadronEne15x15, *t_cHadronEne15x15_1, *t_cHadronEne15x15_2, *t_cHadronEne15x15_3;
0139   std::vector<double> *t_nHadronEne15x15;
0140   std::vector<double> *t_photonEne15x15;
0141   std::vector<double> *t_eleEne15x15;
0142   std::vector<double> *t_muEne15x15;
0143 
0144   std::vector<double> *t_maxNearP11x11;
0145   std::vector<double> *t_cHadronEne11x11, *t_cHadronEne11x11_1, *t_cHadronEne11x11_2, *t_cHadronEne11x11_3;
0146   std::vector<double> *t_nHadronEne11x11;
0147   std::vector<double> *t_photonEne11x11;
0148   std::vector<double> *t_eleEne11x11;
0149   std::vector<double> *t_muEne11x11;
0150 
0151   std::vector<double> *t_maxNearP9x9;
0152   std::vector<double> *t_cHadronEne9x9, *t_cHadronEne9x9_1, *t_cHadronEne9x9_2, *t_cHadronEne9x9_3;
0153   std::vector<double> *t_nHadronEne9x9;
0154   std::vector<double> *t_photonEne9x9;
0155   std::vector<double> *t_eleEne9x9;
0156   std::vector<double> *t_muEne9x9;
0157 
0158   std::vector<double> *t_maxNearP7x7;
0159   std::vector<double> *t_cHadronEne7x7, *t_cHadronEne7x7_1, *t_cHadronEne7x7_2, *t_cHadronEne7x7_3;
0160   std::vector<double> *t_nHadronEne7x7;
0161   std::vector<double> *t_photonEne7x7;
0162   std::vector<double> *t_eleEne7x7;
0163   std::vector<double> *t_muEne7x7;
0164 
0165   std::vector<double> *t_maxNearP3x3;
0166   std::vector<double> *t_cHadronEne3x3, *t_cHadronEne3x3_1, *t_cHadronEne3x3_2, *t_cHadronEne3x3_3;
0167   std::vector<double> *t_nHadronEne3x3;
0168   std::vector<double> *t_photonEne3x3;
0169   std::vector<double> *t_eleEne3x3;
0170   std::vector<double> *t_muEne3x3;
0171 
0172   std::vector<double> *t_maxNearP1x1;
0173   std::vector<double> *t_cHadronEne1x1, *t_cHadronEne1x1_1, *t_cHadronEne1x1_2, *t_cHadronEne1x1_3;
0174   std::vector<double> *t_nHadronEne1x1;
0175   std::vector<double> *t_photonEne1x1;
0176   std::vector<double> *t_eleEne1x1;
0177   std::vector<double> *t_muEne1x1;
0178 
0179   std::vector<double> *t_maxNearPHC1x1;
0180   std::vector<double> *t_cHadronEneHC1x1, *t_cHadronEneHC1x1_1, *t_cHadronEneHC1x1_2, *t_cHadronEneHC1x1_3;
0181   std::vector<double> *t_nHadronEneHC1x1;
0182   std::vector<double> *t_photonEneHC1x1;
0183   std::vector<double> *t_eleEneHC1x1;
0184   std::vector<double> *t_muEneHC1x1;
0185 
0186   std::vector<double> *t_maxNearPHC3x3;
0187   std::vector<double> *t_cHadronEneHC3x3, *t_cHadronEneHC3x3_1, *t_cHadronEneHC3x3_2, *t_cHadronEneHC3x3_3;
0188   std::vector<double> *t_nHadronEneHC3x3;
0189   std::vector<double> *t_photonEneHC3x3;
0190   std::vector<double> *t_eleEneHC3x3;
0191   std::vector<double> *t_muEneHC3x3;
0192 
0193   std::vector<double> *t_maxNearPHC5x5;
0194   std::vector<double> *t_cHadronEneHC5x5, *t_cHadronEneHC5x5_1, *t_cHadronEneHC5x5_2, *t_cHadronEneHC5x5_3;
0195   std::vector<double> *t_nHadronEneHC5x5;
0196   std::vector<double> *t_photonEneHC5x5;
0197   std::vector<double> *t_eleEneHC5x5;
0198   std::vector<double> *t_muEneHC5x5;
0199 
0200   std::vector<double> *t_maxNearPHC7x7;
0201   std::vector<double> *t_cHadronEneHC7x7, *t_cHadronEneHC7x7_1, *t_cHadronEneHC7x7_2, *t_cHadronEneHC7x7_3;
0202   std::vector<double> *t_nHadronEneHC7x7;
0203   std::vector<double> *t_photonEneHC7x7;
0204   std::vector<double> *t_eleEneHC7x7;
0205   std::vector<double> *t_muEneHC7x7;
0206 
0207   std::vector<double> *t_maxNearPR;
0208   std::vector<double> *t_cHadronEneR, *t_cHadronEneR_1, *t_cHadronEneR_2, *t_cHadronEneR_3;
0209   std::vector<double> *t_nHadronEneR;
0210   std::vector<double> *t_photonEneR;
0211   std::vector<double> *t_eleEneR;
0212   std::vector<double> *t_muEneR;
0213 
0214   std::vector<double> *t_maxNearPIsoR;
0215   std::vector<double> *t_cHadronEneIsoR, *t_cHadronEneIsoR_1, *t_cHadronEneIsoR_2, *t_cHadronEneIsoR_3;
0216   std::vector<double> *t_nHadronEneIsoR;
0217   std::vector<double> *t_photonEneIsoR;
0218   std::vector<double> *t_eleEneIsoR;
0219   std::vector<double> *t_muEneIsoR;
0220 
0221   std::vector<double> *t_maxNearPHCR;
0222   std::vector<double> *t_cHadronEneHCR, *t_cHadronEneHCR_1, *t_cHadronEneHCR_2, *t_cHadronEneHCR_3;
0223   std::vector<double> *t_nHadronEneHCR;
0224   std::vector<double> *t_photonEneHCR;
0225   std::vector<double> *t_eleEneHCR;
0226   std::vector<double> *t_muEneHCR;
0227 
0228   std::vector<double> *t_maxNearPIsoHCR;
0229   std::vector<double> *t_cHadronEneIsoHCR, *t_cHadronEneIsoHCR_1, *t_cHadronEneIsoHCR_2, *t_cHadronEneIsoHCR_3;
0230   std::vector<double> *t_nHadronEneIsoHCR;
0231   std::vector<double> *t_photonEneIsoHCR;
0232   std::vector<double> *t_eleEneIsoHCR;
0233   std::vector<double> *t_muEneIsoHCR;
0234 
0235   spr::genSimInfo isoinfo1x1, isoinfo3x3, isoinfo7x7, isoinfo9x9, isoinfo11x11;
0236   spr::genSimInfo isoinfo15x15, isoinfo21x21, isoinfo25x25, isoinfo31x31;
0237   spr::genSimInfo isoinfoHC1x1, isoinfoHC3x3, isoinfoHC5x5, isoinfoHC7x7;
0238   spr::genSimInfo isoinfoR, isoinfoIsoR, isoinfoHCR, isoinfoIsoHCR;
0239 };
0240 
0241 StudyCaloGen::StudyCaloGen(const edm::ParameterSet &iConfig)
0242     : ptMin_(iConfig.getUntrackedParameter<double>("PTMin", 1.0)),
0243       etaMax_(iConfig.getUntrackedParameter<double>("MaxChargedHadronEta", 2.5)),
0244       pCutIsolate_(iConfig.getUntrackedParameter<double>("PMaxIsolation", 20.0)),
0245       a_Isolation_(iConfig.getUntrackedParameter<bool>("UseConeIsolation", false)),
0246       genSrc_(iConfig.getUntrackedParameter("GenSrc", std::string("generatorSmeared"))),
0247       useHepMC_(iConfig.getUntrackedParameter<bool>("UseHepMC", false)),
0248       a_coneR_(iConfig.getUntrackedParameter<double>("ConeRadius", 34.98)),
0249       a_mipR_(iConfig.getUntrackedParameter<double>("ConeRadiusMIP", 14.0)),
0250       verbosity_(iConfig.getUntrackedParameter<int>("Verbosity", 0)) {
0251   usesResource(TFileService::kSharedResource);
0252 
0253   a_charIsoR_ = a_coneR_ + 28.9;
0254   a_neutIsoR_ = a_charIsoR_ * 0.726;
0255 
0256   tok_hepmc_ = consumes<edm::HepMCProduct>(edm::InputTag(genSrc_));
0257   tok_genParticles_ = consumes<reco::GenParticleCollection>(edm::InputTag(genSrc_));
0258 
0259   if (!strcmp("Dummy", genSrc_.c_str())) {
0260     if (useHepMC_)
0261       genSrc_ = "generatorSmeared";
0262     else
0263       genSrc_ = "genParticles";
0264   }
0265   edm::LogVerbatim("IsoTrack") << "Generator Source " << genSrc_ << " Use HepMC " << useHepMC_ << " ptMin " << ptMin_
0266                                << " etaMax " << etaMax_ << "\n a_coneR " << a_coneR_ << " a_charIsoR " << a_charIsoR_
0267                                << " a_neutIsoR " << a_neutIsoR_ << " a_mipR " << a_mipR_ << " debug " << verbosity_
0268                                << "\nIsolation Flag " << a_Isolation_ << " with cut " << pCutIsolate_ << " GeV";
0269 
0270   tok_geom_ = esConsumes<CaloGeometry, CaloGeometryRecord>();
0271   tok_caloTopology_ = esConsumes<CaloTopology, CaloTopologyRecord>();
0272   tok_topo_ = esConsumes<HcalTopology, HcalRecNumberingRecord>();
0273   tok_magField_ = esConsumes<MagneticField, IdealMagneticFieldRecord>();
0274   tok_pdt_ = esConsumes<HepPDT::ParticleDataTable, PDTRecord>();
0275 }
0276 
0277 void StudyCaloGen::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0278   edm::ParameterSetDescription desc;
0279   desc.addUntracked<std::string>("GenSrc", "genParticles");
0280   desc.addUntracked<bool>("UseHepMC", false);
0281   desc.addUntracked<double>("ChargedHadronSeedP", 1.0);
0282   desc.addUntracked<double>("PTMin", 1.0);
0283   desc.addUntracked<double>("MaxChargedHadronEta", 2.5);
0284   desc.addUntracked<double>("ConeRadius", 34.98);
0285   desc.addUntracked<double>("ConeRadiusMIP", 14.0);
0286   desc.addUntracked<bool>("UseConeIsolation", true);
0287   desc.addUntracked<double>("PMaxIsolation", 5.0);
0288   desc.addUntracked<int>("Verbosity", 0);
0289   descriptions.add("studyCaloGen", desc);
0290 }
0291 
0292 void StudyCaloGen::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0293   clearTreeVectors();
0294 
0295   nEventProc++;
0296   const MagneticField *bField = &iSetup.getData(tok_magField_);
0297 
0298   // get particle data table
0299   const HepPDT::ParticleDataTable *pdt = &iSetup.getData(tok_pdt_);
0300 
0301   // get handle to HEPMCProduct
0302   edm::Handle<edm::HepMCProduct> hepmc;
0303   edm::Handle<reco::GenParticleCollection> genParticles;
0304   if (useHepMC_)
0305     iEvent.getByToken(tok_hepmc_, hepmc);
0306   else
0307     iEvent.getByToken(tok_genParticles_, genParticles);
0308 
0309   const CaloGeometry *geo = &iSetup.getData(tok_geom_);
0310   const CaloTopology *caloTopology = &iSetup.getData(tok_caloTopology_);
0311   const HcalTopology *theHBHETopology = &iSetup.getData(tok_topo_);
0312 
0313   GlobalPoint posVec, posECAL;
0314   math::XYZTLorentzVector momVec;
0315   if (verbosity_ > 0)
0316     edm::LogVerbatim("IsoTrack") << "event number " << iEvent.id().event();
0317   if (useHepMC_) {
0318     const HepMC::GenEvent *myGenEvent = hepmc->GetEvent();
0319     std::vector<spr::propagatedGenTrackID> trackIDs = spr::propagateCALO(myGenEvent, pdt, geo, bField, etaMax_, false);
0320 
0321     for (unsigned int indx = 0; indx < trackIDs.size(); ++indx) {
0322       int charge = trackIDs[indx].charge;
0323       HepMC::GenEvent::particle_const_iterator p = trackIDs[indx].trkItr;
0324       momVec = math::XYZTLorentzVector(
0325           (*p)->momentum().px(), (*p)->momentum().py(), (*p)->momentum().pz(), (*p)->momentum().e());
0326       if (verbosity_ > 1)
0327         edm::LogVerbatim("IsoTrack") << "trkIndx " << indx << " pdgid " << trackIDs[indx].pdgId << " charge " << charge
0328                                      << " momVec " << momVec;
0329       // only stable particles avoiding electrons and muons
0330       if (trackIDs[indx].ok && (std::abs(trackIDs[indx].pdgId) < 11 || std::abs(trackIDs[indx].pdgId) >= 21)) {
0331         // consider particles within a phased space
0332         if (momVec.Pt() > ptMin_ && std::abs(momVec.eta()) < etaMax_) {
0333           posVec = GlobalPoint(0.1 * (*p)->production_vertex()->position().x(),
0334                                0.1 * (*p)->production_vertex()->position().y(),
0335                                0.1 * (*p)->production_vertex()->position().z());
0336           posECAL = trackIDs[indx].pointECAL;
0337           fillTrack(posVec, momVec, posECAL, trackIDs[indx].pdgId, trackIDs[indx].okECAL, true);
0338           if (verbosity_ > 1)
0339             edm::LogVerbatim("IsoTrack") << "posECAL " << posECAL << " okECAL " << trackIDs[indx].okECAL << "okHCAL "
0340                                          << trackIDs[indx].okHCAL;
0341           if (trackIDs[indx].okECAL) {
0342             if (std::abs(charge) > 0) {
0343               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 0, 0, isoinfo1x1, false);
0344               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 1, 1, isoinfo3x3, false);
0345               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 3, 3, isoinfo7x7, false);
0346               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 4, 4, isoinfo9x9, false);
0347               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 5, 5, isoinfo11x11, false);
0348               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 7, 7, isoinfo15x15, false);
0349               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 10, 10, isoinfo21x21, false);
0350               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 12, 12, isoinfo25x25, false);
0351               spr::eGenSimInfo(trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 15, 15, isoinfo31x31, false);
0352               spr::eGenSimInfo(trackIDs[indx].detIdECAL,
0353                                p,
0354                                trackIDs,
0355                                geo,
0356                                caloTopology,
0357                                a_mipR_,
0358                                trackIDs[indx].directionECAL,
0359                                isoinfoR,
0360                                false);
0361               spr::eGenSimInfo(trackIDs[indx].detIdECAL,
0362                                p,
0363                                trackIDs,
0364                                geo,
0365                                caloTopology,
0366                                a_neutIsoR_,
0367                                trackIDs[indx].directionECAL,
0368                                isoinfoIsoR,
0369                                false);
0370               if (trackIDs[indx].okHCAL) {
0371                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 0, 0, isoinfoHC1x1, false);
0372                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 1, 1, isoinfoHC3x3, false);
0373                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 2, 2, isoinfoHC5x5, false);
0374                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 3, 3, isoinfoHC7x7, false);
0375                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL,
0376                                  p,
0377                                  trackIDs,
0378                                  geo,
0379                                  theHBHETopology,
0380                                  a_coneR_,
0381                                  trackIDs[indx].directionHCAL,
0382                                  isoinfoHCR,
0383                                  false);
0384                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL,
0385                                  p,
0386                                  trackIDs,
0387                                  geo,
0388                                  theHBHETopology,
0389                                  a_charIsoR_,
0390                                  trackIDs[indx].directionHCAL,
0391                                  isoinfoIsoHCR,
0392                                  false);
0393               }
0394 
0395               bool saveTrack = true;
0396               if (a_Isolation_)
0397                 saveTrack = (isoinfoR.maxNearP < pCutIsolate_);
0398               else
0399                 saveTrack = (isoinfo7x7.maxNearP < pCutIsolate_);
0400               if (saveTrack)
0401                 fillIsolatedTrack(momVec, posECAL, trackIDs[indx].pdgId);
0402             }
0403           }
0404         } else {  // stabale particles within |eta|=2.5
0405           fillTrack(posVec, momVec, posECAL, 0, false, false);
0406         }
0407       }
0408     }
0409 
0410     unsigned int indx;
0411     HepMC::GenEvent::particle_const_iterator p;
0412     for (p = myGenEvent->particles_begin(), indx = 0; p != myGenEvent->particles_end(); ++p, ++indx) {
0413       int pdgId = ((*p)->pdg_id());
0414       int ix = particleCode(pdgId);
0415       if (ix >= 0) {
0416         double pp = (*p)->momentum().rho();
0417         double eta = (*p)->momentum().eta();
0418         h_pEta[ix]->Fill(pp, eta);
0419       }
0420     }
0421   } else {  // loop over gen particles
0422     std::vector<spr::propagatedGenParticleID> trackIDs =
0423         spr::propagateCALO(genParticles, pdt, geo, bField, etaMax_, (verbosity_ > 0));
0424 
0425     for (unsigned int indx = 0; indx < trackIDs.size(); ++indx) {
0426       int charge = trackIDs[indx].charge;
0427       reco::GenParticleCollection::const_iterator p = trackIDs[indx].trkItr;
0428 
0429       momVec = math::XYZTLorentzVector(p->momentum().x(), p->momentum().y(), p->momentum().z(), p->energy());
0430       if (verbosity_ > 1)
0431         edm::LogVerbatim("IsoTrack") << "trkIndx " << indx << " pdgid " << trackIDs[indx].pdgId << " charge " << charge
0432                                      << " momVec " << momVec;
0433       // only stable particles avoiding electrons and muons
0434       if (trackIDs[indx].ok && std::abs(trackIDs[indx].pdgId) > 21) {
0435         // consider particles within a phased space
0436         if (verbosity_ > 1)
0437           edm::LogVerbatim("IsoTrack") << " pt " << momVec.Pt() << " eta " << momVec.eta();
0438         if (momVec.Pt() > ptMin_ && std::abs(momVec.eta()) < etaMax_) {
0439           posVec = GlobalPoint(p->vertex().x(), p->vertex().y(), p->vertex().z());
0440           posECAL = trackIDs[indx].pointECAL;
0441           if (verbosity_ > 0)
0442             edm::LogVerbatim("IsoTrack") << "posECAL " << posECAL << " okECAL " << trackIDs[indx].okECAL << "okHCAL "
0443                                          << trackIDs[indx].okHCAL;
0444           fillTrack(posVec, momVec, posECAL, trackIDs[indx].pdgId, trackIDs[indx].okECAL, true);
0445           if (trackIDs[indx].okECAL) {
0446             if (std::abs(charge) > 0) {
0447               spr::eGenSimInfo(
0448                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 0, 0, isoinfo1x1, verbosity_ > 1);
0449               spr::eGenSimInfo(
0450                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 1, 1, isoinfo3x3, verbosity_ > 0);
0451               spr::eGenSimInfo(
0452                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 3, 3, isoinfo7x7, verbosity_ > 1);
0453               spr::eGenSimInfo(
0454                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 4, 4, isoinfo9x9, verbosity_ > 1);
0455               spr::eGenSimInfo(
0456                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 5, 5, isoinfo11x11, verbosity_ > 1);
0457               spr::eGenSimInfo(
0458                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 7, 7, isoinfo15x15, verbosity_ > 1);
0459               spr::eGenSimInfo(
0460                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 10, 10, isoinfo21x21, verbosity_ > 1);
0461               spr::eGenSimInfo(
0462                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 12, 12, isoinfo25x25, verbosity_ > 1);
0463               spr::eGenSimInfo(
0464                   trackIDs[indx].detIdECAL, p, trackIDs, geo, caloTopology, 15, 15, isoinfo31x31, verbosity_ > 1);
0465               spr::eGenSimInfo(trackIDs[indx].detIdECAL,
0466                                p,
0467                                trackIDs,
0468                                geo,
0469                                caloTopology,
0470                                a_mipR_,
0471                                trackIDs[indx].directionECAL,
0472                                isoinfoR,
0473                                verbosity_ > 1);
0474               spr::eGenSimInfo(trackIDs[indx].detIdECAL,
0475                                p,
0476                                trackIDs,
0477                                geo,
0478                                caloTopology,
0479                                a_neutIsoR_,
0480                                trackIDs[indx].directionECAL,
0481                                isoinfoIsoR,
0482                                verbosity_ > 1);
0483               if (trackIDs[indx].okHCAL) {
0484                 spr::hGenSimInfo(
0485                     trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 0, 0, isoinfoHC1x1, verbosity_ > 1);
0486                 spr::hGenSimInfo(
0487                     trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 1, 1, isoinfoHC3x3, verbosity_ > 1);
0488                 spr::hGenSimInfo(
0489                     trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 2, 2, isoinfoHC5x5, verbosity_ > 1);
0490                 spr::hGenSimInfo(
0491                     trackIDs[indx].detIdHCAL, p, trackIDs, theHBHETopology, 3, 3, isoinfoHC7x7, verbosity_ > 1);
0492                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL,
0493                                  p,
0494                                  trackIDs,
0495                                  geo,
0496                                  theHBHETopology,
0497                                  a_coneR_,
0498                                  trackIDs[indx].directionHCAL,
0499                                  isoinfoHCR,
0500                                  verbosity_ > 1);
0501                 spr::hGenSimInfo(trackIDs[indx].detIdHCAL,
0502                                  p,
0503                                  trackIDs,
0504                                  geo,
0505                                  theHBHETopology,
0506                                  a_charIsoR_,
0507                                  trackIDs[indx].directionHCAL,
0508                                  isoinfoIsoHCR,
0509                                  verbosity_ > 1);
0510               }
0511 
0512               bool saveTrack = true;
0513               if (a_Isolation_)
0514                 saveTrack = (isoinfoIsoR.maxNearP < pCutIsolate_);
0515               else
0516                 saveTrack = (isoinfo7x7.maxNearP < pCutIsolate_);
0517               if (saveTrack)
0518                 fillIsolatedTrack(momVec, posECAL, trackIDs[indx].pdgId);
0519             }
0520           }
0521         } else {  // stabale particles within |eta|=2.5
0522           fillTrack(posVec, momVec, posECAL, 0, false, false);
0523         }
0524       }
0525     }  // loop over gen particles
0526 
0527     unsigned int indx;
0528     reco::GenParticleCollection::const_iterator p;
0529     for (p = genParticles->begin(), indx = 0; p != genParticles->end(); ++p, ++indx) {
0530       int pdgId = (p->pdgId());
0531       int ix = particleCode(pdgId);
0532       if (ix >= 0) {
0533         double pp = (p->momentum()).R();
0534         double eta = (p->momentum()).Eta();
0535         h_pEta[ix]->Fill(pp, eta);
0536       }
0537     }
0538   }
0539 
0540   //t_nEvtProc->push_back(nEventProc);
0541   h_NEventProc->SetBinContent(1, nEventProc);
0542   tree_->Fill();
0543 }
0544 
0545 void StudyCaloGen::beginJob() {
0546   nEventProc = 0;
0547 
0548   double tempgen_TH[NPBins_ + 1] = {0.0, 5.0, 12.0, 300.0};
0549   for (int i = 0; i <= NPBins_; i++)
0550     genPartPBins_[i] = tempgen_TH[i];
0551 
0552   double tempgen_Eta[NEtaBins_ + 1] = {0.0, 0.5, 1.1, 1.7, 2.3};
0553   for (int i = 0; i <= NEtaBins_; i++)
0554     genPartEtaBins_[i] = tempgen_Eta[i];
0555 
0556   bookHistograms();
0557 }
0558 
0559 void StudyCaloGen::fillTrack(
0560     GlobalPoint &posVec, math::XYZTLorentzVector &momVec, GlobalPoint &posECAL, int pdgId, bool okECAL, bool accept) {
0561   if (accept) {
0562     t_isoTrkPAll->push_back(momVec.P());
0563     t_isoTrkPtAll->push_back(momVec.Pt());
0564     t_isoTrkPhiAll->push_back(momVec.phi());
0565     t_isoTrkEtaAll->push_back(momVec.eta());
0566     t_isoTrkPdgIdAll->push_back(pdgId);
0567     if (okECAL) {
0568       double phi1 = momVec.phi();
0569       double phi2 = (posECAL - posVec).phi();
0570       double dphi = reco::deltaPhi(phi1, phi2);
0571       double deta = momVec.eta() - (posECAL - posVec).eta();
0572       t_isoTrkDPhiAll->push_back(dphi);
0573       t_isoTrkDEtaAll->push_back(deta);
0574     } else {
0575       t_isoTrkDPhiAll->push_back(999.0);
0576       t_isoTrkDEtaAll->push_back(999.0);
0577     }
0578   } else {
0579     t_isoTrkDPhiAll->push_back(-999.0);
0580     t_isoTrkDEtaAll->push_back(-999.0);
0581   }
0582 }
0583 
0584 void StudyCaloGen::fillIsolatedTrack(math::XYZTLorentzVector &momVec, GlobalPoint &posECAL, int pdgId) {
0585   t_isoTrkP->push_back(momVec.P());
0586   t_isoTrkPt->push_back(momVec.Pt());
0587   t_isoTrkEne->push_back(momVec.E());
0588   t_isoTrkEta->push_back(momVec.eta());
0589   t_isoTrkPhi->push_back(momVec.phi());
0590   t_isoTrkEtaEC->push_back(posECAL.eta());
0591   t_isoTrkPhiEC->push_back(posECAL.phi());
0592   t_isoTrkPdgId->push_back(pdgId);
0593 
0594   t_maxNearP31x31->push_back(isoinfo31x31.maxNearP);
0595   t_cHadronEne31x31->push_back(isoinfo31x31.cHadronEne);
0596   t_cHadronEne31x31_1->push_back(isoinfo31x31.cHadronEne_[0]);
0597   t_cHadronEne31x31_2->push_back(isoinfo31x31.cHadronEne_[1]);
0598   t_cHadronEne31x31_3->push_back(isoinfo31x31.cHadronEne_[2]);
0599   t_nHadronEne31x31->push_back(isoinfo31x31.nHadronEne);
0600   t_photonEne31x31->push_back(isoinfo31x31.photonEne);
0601   t_eleEne31x31->push_back(isoinfo31x31.eleEne);
0602   t_muEne31x31->push_back(isoinfo31x31.muEne);
0603 
0604   t_maxNearP25x25->push_back(isoinfo25x25.maxNearP);
0605   t_cHadronEne25x25->push_back(isoinfo25x25.cHadronEne);
0606   t_cHadronEne25x25_1->push_back(isoinfo25x25.cHadronEne_[0]);
0607   t_cHadronEne25x25_2->push_back(isoinfo25x25.cHadronEne_[1]);
0608   t_cHadronEne25x25_3->push_back(isoinfo25x25.cHadronEne_[2]);
0609   t_nHadronEne25x25->push_back(isoinfo25x25.nHadronEne);
0610   t_photonEne25x25->push_back(isoinfo25x25.photonEne);
0611   t_eleEne25x25->push_back(isoinfo25x25.eleEne);
0612   t_muEne25x25->push_back(isoinfo25x25.muEne);
0613 
0614   t_maxNearP21x21->push_back(isoinfo21x21.maxNearP);
0615   t_cHadronEne21x21->push_back(isoinfo21x21.cHadronEne);
0616   t_cHadronEne21x21_1->push_back(isoinfo21x21.cHadronEne_[0]);
0617   t_cHadronEne21x21_2->push_back(isoinfo21x21.cHadronEne_[1]);
0618   t_cHadronEne21x21_3->push_back(isoinfo21x21.cHadronEne_[2]);
0619   t_nHadronEne21x21->push_back(isoinfo21x21.nHadronEne);
0620   t_photonEne21x21->push_back(isoinfo21x21.photonEne);
0621   t_eleEne21x21->push_back(isoinfo21x21.eleEne);
0622   t_muEne21x21->push_back(isoinfo21x21.muEne);
0623 
0624   t_maxNearP15x15->push_back(isoinfo15x15.maxNearP);
0625   t_cHadronEne15x15->push_back(isoinfo15x15.cHadronEne);
0626   t_cHadronEne15x15_1->push_back(isoinfo15x15.cHadronEne_[0]);
0627   t_cHadronEne15x15_2->push_back(isoinfo15x15.cHadronEne_[1]);
0628   t_cHadronEne15x15_3->push_back(isoinfo15x15.cHadronEne_[2]);
0629   t_nHadronEne15x15->push_back(isoinfo15x15.nHadronEne);
0630   t_photonEne15x15->push_back(isoinfo15x15.photonEne);
0631   t_eleEne15x15->push_back(isoinfo15x15.eleEne);
0632   t_muEne15x15->push_back(isoinfo15x15.muEne);
0633 
0634   t_maxNearP11x11->push_back(isoinfo11x11.maxNearP);
0635   t_cHadronEne11x11->push_back(isoinfo11x11.cHadronEne);
0636   t_cHadronEne11x11_1->push_back(isoinfo11x11.cHadronEne_[0]);
0637   t_cHadronEne11x11_2->push_back(isoinfo11x11.cHadronEne_[1]);
0638   t_cHadronEne11x11_3->push_back(isoinfo11x11.cHadronEne_[2]);
0639   t_nHadronEne11x11->push_back(isoinfo11x11.nHadronEne);
0640   t_photonEne11x11->push_back(isoinfo11x11.photonEne);
0641   t_eleEne11x11->push_back(isoinfo11x11.eleEne);
0642   t_muEne11x11->push_back(isoinfo11x11.muEne);
0643 
0644   t_maxNearP9x9->push_back(isoinfo9x9.maxNearP);
0645   t_cHadronEne9x9->push_back(isoinfo9x9.cHadronEne);
0646   t_cHadronEne9x9_1->push_back(isoinfo9x9.cHadronEne_[0]);
0647   t_cHadronEne9x9_2->push_back(isoinfo9x9.cHadronEne_[1]);
0648   t_cHadronEne9x9_3->push_back(isoinfo9x9.cHadronEne_[2]);
0649   t_nHadronEne9x9->push_back(isoinfo9x9.nHadronEne);
0650   t_photonEne9x9->push_back(isoinfo9x9.photonEne);
0651   t_eleEne9x9->push_back(isoinfo9x9.eleEne);
0652   t_muEne9x9->push_back(isoinfo9x9.muEne);
0653 
0654   t_maxNearP7x7->push_back(isoinfo7x7.maxNearP);
0655   t_cHadronEne7x7->push_back(isoinfo7x7.cHadronEne);
0656   t_cHadronEne7x7_1->push_back(isoinfo7x7.cHadronEne_[0]);
0657   t_cHadronEne7x7_2->push_back(isoinfo7x7.cHadronEne_[1]);
0658   t_cHadronEne7x7_3->push_back(isoinfo7x7.cHadronEne_[2]);
0659   t_nHadronEne7x7->push_back(isoinfo7x7.nHadronEne);
0660   t_photonEne7x7->push_back(isoinfo7x7.photonEne);
0661   t_eleEne7x7->push_back(isoinfo7x7.eleEne);
0662   t_muEne7x7->push_back(isoinfo7x7.muEne);
0663 
0664   t_maxNearP3x3->push_back(isoinfo3x3.maxNearP);
0665   t_cHadronEne3x3->push_back(isoinfo3x3.cHadronEne);
0666   t_cHadronEne3x3_1->push_back(isoinfo3x3.cHadronEne_[0]);
0667   t_cHadronEne3x3_2->push_back(isoinfo3x3.cHadronEne_[1]);
0668   t_cHadronEne3x3_3->push_back(isoinfo3x3.cHadronEne_[2]);
0669   t_nHadronEne3x3->push_back(isoinfo3x3.nHadronEne);
0670   t_photonEne3x3->push_back(isoinfo3x3.photonEne);
0671   t_eleEne3x3->push_back(isoinfo3x3.eleEne);
0672   t_muEne3x3->push_back(isoinfo3x3.muEne);
0673 
0674   t_maxNearP1x1->push_back(isoinfo1x1.maxNearP);
0675   t_cHadronEne1x1->push_back(isoinfo1x1.cHadronEne);
0676   t_cHadronEne1x1_1->push_back(isoinfo1x1.cHadronEne_[0]);
0677   t_cHadronEne1x1_2->push_back(isoinfo1x1.cHadronEne_[1]);
0678   t_cHadronEne1x1_3->push_back(isoinfo1x1.cHadronEne_[2]);
0679   t_nHadronEne1x1->push_back(isoinfo1x1.nHadronEne);
0680   t_photonEne1x1->push_back(isoinfo1x1.photonEne);
0681   t_eleEne1x1->push_back(isoinfo1x1.eleEne);
0682   t_muEne1x1->push_back(isoinfo1x1.muEne);
0683 
0684   t_maxNearPHC1x1->push_back(isoinfoHC1x1.maxNearP);
0685   t_cHadronEneHC1x1->push_back(isoinfoHC1x1.cHadronEne);
0686   t_cHadronEneHC1x1_1->push_back(isoinfoHC1x1.cHadronEne_[0]);
0687   t_cHadronEneHC1x1_2->push_back(isoinfoHC1x1.cHadronEne_[1]);
0688   t_cHadronEneHC1x1_3->push_back(isoinfoHC1x1.cHadronEne_[2]);
0689   t_nHadronEneHC1x1->push_back(isoinfoHC1x1.nHadronEne);
0690   t_photonEneHC1x1->push_back(isoinfoHC1x1.photonEne);
0691   t_eleEneHC1x1->push_back(isoinfoHC1x1.eleEne);
0692   t_muEneHC1x1->push_back(isoinfoHC1x1.muEne);
0693 
0694   t_maxNearPHC3x3->push_back(isoinfoHC3x3.maxNearP);
0695   t_cHadronEneHC3x3->push_back(isoinfoHC3x3.cHadronEne);
0696   t_cHadronEneHC3x3_1->push_back(isoinfoHC3x3.cHadronEne_[0]);
0697   t_cHadronEneHC3x3_2->push_back(isoinfoHC3x3.cHadronEne_[1]);
0698   t_cHadronEneHC3x3_3->push_back(isoinfoHC3x3.cHadronEne_[2]);
0699   t_nHadronEneHC3x3->push_back(isoinfoHC3x3.nHadronEne);
0700   t_photonEneHC3x3->push_back(isoinfoHC3x3.photonEne);
0701   t_eleEneHC3x3->push_back(isoinfoHC3x3.eleEne);
0702   t_muEneHC3x3->push_back(isoinfoHC3x3.muEne);
0703 
0704   t_maxNearPHC5x5->push_back(isoinfoHC5x5.maxNearP);
0705   t_cHadronEneHC5x5->push_back(isoinfoHC5x5.cHadronEne);
0706   t_cHadronEneHC5x5_1->push_back(isoinfoHC5x5.cHadronEne_[0]);
0707   t_cHadronEneHC5x5_2->push_back(isoinfoHC5x5.cHadronEne_[1]);
0708   t_cHadronEneHC5x5_3->push_back(isoinfoHC5x5.cHadronEne_[2]);
0709   t_nHadronEneHC5x5->push_back(isoinfoHC5x5.nHadronEne);
0710   t_photonEneHC5x5->push_back(isoinfoHC5x5.photonEne);
0711   t_eleEneHC5x5->push_back(isoinfoHC5x5.eleEne);
0712   t_muEneHC5x5->push_back(isoinfoHC5x5.muEne);
0713 
0714   t_maxNearPHC7x7->push_back(isoinfoHC7x7.maxNearP);
0715   t_cHadronEneHC7x7->push_back(isoinfoHC7x7.cHadronEne);
0716   t_cHadronEneHC7x7_1->push_back(isoinfoHC7x7.cHadronEne_[0]);
0717   t_cHadronEneHC7x7_2->push_back(isoinfoHC7x7.cHadronEne_[1]);
0718   t_cHadronEneHC7x7_3->push_back(isoinfoHC7x7.cHadronEne_[2]);
0719   t_nHadronEneHC7x7->push_back(isoinfoHC7x7.nHadronEne);
0720   t_photonEneHC7x7->push_back(isoinfoHC7x7.photonEne);
0721   t_eleEneHC7x7->push_back(isoinfoHC7x7.eleEne);
0722   t_muEneHC7x7->push_back(isoinfoHC7x7.muEne);
0723 
0724   t_maxNearPR->push_back(isoinfoR.maxNearP);
0725   t_cHadronEneR->push_back(isoinfoR.cHadronEne);
0726   t_cHadronEneR_1->push_back(isoinfoR.cHadronEne_[0]);
0727   t_cHadronEneR_2->push_back(isoinfoR.cHadronEne_[1]);
0728   t_cHadronEneR_3->push_back(isoinfoR.cHadronEne_[2]);
0729   t_nHadronEneR->push_back(isoinfoR.nHadronEne);
0730   t_photonEneR->push_back(isoinfoR.photonEne);
0731   t_eleEneR->push_back(isoinfoR.eleEne);
0732   t_muEneR->push_back(isoinfoR.muEne);
0733 
0734   t_maxNearPIsoR->push_back(isoinfoIsoR.maxNearP);
0735   t_cHadronEneIsoR->push_back(isoinfoIsoR.cHadronEne);
0736   t_cHadronEneIsoR_1->push_back(isoinfoIsoR.cHadronEne_[0]);
0737   t_cHadronEneIsoR_2->push_back(isoinfoIsoR.cHadronEne_[1]);
0738   t_cHadronEneIsoR_3->push_back(isoinfoIsoR.cHadronEne_[2]);
0739   t_nHadronEneIsoR->push_back(isoinfoIsoR.nHadronEne);
0740   t_photonEneIsoR->push_back(isoinfoIsoR.photonEne);
0741   t_eleEneIsoR->push_back(isoinfoIsoR.eleEne);
0742   t_muEneIsoR->push_back(isoinfoIsoR.muEne);
0743 
0744   t_maxNearPHCR->push_back(isoinfoHCR.maxNearP);
0745   t_cHadronEneHCR->push_back(isoinfoHCR.cHadronEne);
0746   t_cHadronEneHCR_1->push_back(isoinfoHCR.cHadronEne_[0]);
0747   t_cHadronEneHCR_2->push_back(isoinfoHCR.cHadronEne_[1]);
0748   t_cHadronEneHCR_3->push_back(isoinfoHCR.cHadronEne_[2]);
0749   t_nHadronEneHCR->push_back(isoinfoHCR.nHadronEne);
0750   t_photonEneHCR->push_back(isoinfoHCR.photonEne);
0751   t_eleEneHCR->push_back(isoinfoHCR.eleEne);
0752   t_muEneHCR->push_back(isoinfoHCR.muEne);
0753 
0754   t_maxNearPIsoHCR->push_back(isoinfoIsoHCR.maxNearP);
0755   t_cHadronEneIsoHCR->push_back(isoinfoIsoHCR.cHadronEne);
0756   t_cHadronEneIsoHCR_1->push_back(isoinfoIsoHCR.cHadronEne_[0]);
0757   t_cHadronEneIsoHCR_2->push_back(isoinfoIsoHCR.cHadronEne_[1]);
0758   t_cHadronEneIsoHCR_3->push_back(isoinfoIsoHCR.cHadronEne_[2]);
0759   t_nHadronEneIsoHCR->push_back(isoinfoIsoHCR.nHadronEne);
0760   t_photonEneIsoHCR->push_back(isoinfoIsoHCR.photonEne);
0761   t_eleEneIsoHCR->push_back(isoinfoIsoHCR.eleEne);
0762   t_muEneIsoHCR->push_back(isoinfoIsoHCR.muEne);
0763 }
0764 
0765 void StudyCaloGen::bookHistograms() {
0766   edm::Service<TFileService> fs;
0767   //char hname[100], htit[100];
0768 
0769   h_NEventProc = fs->make<TH1I>("h_NEventProc", "h_NEventProc", 2, -0.5, 0.5);
0770 
0771   double pBin[PBins_ + 1] = {0.0,   2.0,   4.0,   6.0,   8.0,   10.0,  20.0,  30.0,  40.0,  50.0,  60.0,
0772                              70.0,  80.0,  90.0,  100.0, 150.0, 200.0, 250.0, 300.0, 350.0, 400.0, 450.0,
0773                              500.0, 550.0, 600.0, 650.0, 700.0, 750.0, 800.0, 850.0, 900.0, 950.0, 1000.0};
0774   double etaBin[EtaBins_ + 1] = {-3.0, -2.9, -2.8, -2.7, -2.6, -2.5, -2.4, -2.3, -2.2, -2.1, -2.0, -1.9, -1.8,
0775                                  -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1, -1.0, -0.9, -0.8, -0.7, -0.6, -0.5,
0776                                  -0.4, -0.3, -0.2, -0.1, 0.0,  0.1,  0.2,  0.3,  0.4,  0.5,  0.6,  0.7,  0.8,
0777                                  0.9,  1.0,  1.1,  1.2,  1.3,  1.4,  1.5,  1.6,  1.7,  1.8,  1.9,  2.0,  2.1,
0778                                  2.2,  2.3,  2.4,  2.5,  2.6,  2.7,  2.8,  2.9,  3.0};
0779   std::string particle[Particles] = {
0780       "electron", "positron", "#gamma", "#pi^+", "#pi^-", "K^+", "K^-", "p", "n", "pbar", "nbar", "K^0_L"};
0781   TFileDirectory dir1 = fs->mkdir("pEta");
0782   char name[20], title[50];
0783   for (int i = 0; i < Particles; ++i) {
0784     sprintf(name, "pEta%d", i);
0785     sprintf(title, "#eta vs momentum for %s", particle[i].c_str());
0786     h_pEta[i] = dir1.make<TH2D>(name, title, PBins_, pBin, EtaBins_, etaBin);
0787   }
0788 
0789   // build the tree
0790   tree_ = fs->make<TTree>("StudyCaloGen", "StudyCaloGen");
0791 
0792   t_isoTrkPAll = new std::vector<double>();
0793   t_isoTrkPtAll = new std::vector<double>();
0794   t_isoTrkPhiAll = new std::vector<double>();
0795   t_isoTrkEtaAll = new std::vector<double>();
0796   t_isoTrkDPhiAll = new std::vector<double>();
0797   t_isoTrkDEtaAll = new std::vector<double>();
0798   t_isoTrkPdgIdAll = new std::vector<double>();
0799 
0800   t_isoTrkP = new std::vector<double>();
0801   t_isoTrkPt = new std::vector<double>();
0802   t_isoTrkEne = new std::vector<double>();
0803   t_isoTrkEta = new std::vector<double>();
0804   t_isoTrkPhi = new std::vector<double>();
0805   t_isoTrkEtaEC = new std::vector<double>();
0806   t_isoTrkPhiEC = new std::vector<double>();
0807   t_isoTrkPdgId = new std::vector<double>();
0808 
0809   t_maxNearP31x31 = new std::vector<double>();
0810   t_cHadronEne31x31 = new std::vector<double>();
0811   t_cHadronEne31x31_1 = new std::vector<double>();
0812   t_cHadronEne31x31_2 = new std::vector<double>();
0813   t_cHadronEne31x31_3 = new std::vector<double>();
0814   t_nHadronEne31x31 = new std::vector<double>();
0815   t_photonEne31x31 = new std::vector<double>();
0816   t_eleEne31x31 = new std::vector<double>();
0817   t_muEne31x31 = new std::vector<double>();
0818 
0819   t_maxNearP25x25 = new std::vector<double>();
0820   t_cHadronEne25x25 = new std::vector<double>();
0821   t_cHadronEne25x25_1 = new std::vector<double>();
0822   t_cHadronEne25x25_2 = new std::vector<double>();
0823   t_cHadronEne25x25_3 = new std::vector<double>();
0824   t_nHadronEne25x25 = new std::vector<double>();
0825   t_photonEne25x25 = new std::vector<double>();
0826   t_eleEne25x25 = new std::vector<double>();
0827   t_muEne25x25 = new std::vector<double>();
0828 
0829   t_maxNearP21x21 = new std::vector<double>();
0830   t_cHadronEne21x21 = new std::vector<double>();
0831   t_cHadronEne21x21_1 = new std::vector<double>();
0832   t_cHadronEne21x21_2 = new std::vector<double>();
0833   t_cHadronEne21x21_3 = new std::vector<double>();
0834   t_nHadronEne21x21 = new std::vector<double>();
0835   t_photonEne21x21 = new std::vector<double>();
0836   t_eleEne21x21 = new std::vector<double>();
0837   t_muEne21x21 = new std::vector<double>();
0838 
0839   t_maxNearP15x15 = new std::vector<double>();
0840   t_cHadronEne15x15 = new std::vector<double>();
0841   t_cHadronEne15x15_1 = new std::vector<double>();
0842   t_cHadronEne15x15_2 = new std::vector<double>();
0843   t_cHadronEne15x15_3 = new std::vector<double>();
0844   t_nHadronEne15x15 = new std::vector<double>();
0845   t_photonEne15x15 = new std::vector<double>();
0846   t_eleEne15x15 = new std::vector<double>();
0847   t_muEne15x15 = new std::vector<double>();
0848 
0849   t_maxNearP11x11 = new std::vector<double>();
0850   t_cHadronEne11x11 = new std::vector<double>();
0851   t_cHadronEne11x11_1 = new std::vector<double>();
0852   t_cHadronEne11x11_2 = new std::vector<double>();
0853   t_cHadronEne11x11_3 = new std::vector<double>();
0854   t_nHadronEne11x11 = new std::vector<double>();
0855   t_photonEne11x11 = new std::vector<double>();
0856   t_eleEne11x11 = new std::vector<double>();
0857   t_muEne11x11 = new std::vector<double>();
0858 
0859   t_maxNearP9x9 = new std::vector<double>();
0860   t_cHadronEne9x9 = new std::vector<double>();
0861   t_cHadronEne9x9_1 = new std::vector<double>();
0862   t_cHadronEne9x9_2 = new std::vector<double>();
0863   t_cHadronEne9x9_3 = new std::vector<double>();
0864   t_nHadronEne9x9 = new std::vector<double>();
0865   t_photonEne9x9 = new std::vector<double>();
0866   t_eleEne9x9 = new std::vector<double>();
0867   t_muEne9x9 = new std::vector<double>();
0868 
0869   t_maxNearP7x7 = new std::vector<double>();
0870   t_cHadronEne7x7 = new std::vector<double>();
0871   t_cHadronEne7x7_1 = new std::vector<double>();
0872   t_cHadronEne7x7_2 = new std::vector<double>();
0873   t_cHadronEne7x7_3 = new std::vector<double>();
0874   t_nHadronEne7x7 = new std::vector<double>();
0875   t_photonEne7x7 = new std::vector<double>();
0876   t_eleEne7x7 = new std::vector<double>();
0877   t_muEne7x7 = new std::vector<double>();
0878 
0879   t_maxNearP3x3 = new std::vector<double>();
0880   t_cHadronEne3x3 = new std::vector<double>();
0881   t_cHadronEne3x3_1 = new std::vector<double>();
0882   t_cHadronEne3x3_2 = new std::vector<double>();
0883   t_cHadronEne3x3_3 = new std::vector<double>();
0884   t_nHadronEne3x3 = new std::vector<double>();
0885   t_photonEne3x3 = new std::vector<double>();
0886   t_eleEne3x3 = new std::vector<double>();
0887   t_muEne3x3 = new std::vector<double>();
0888 
0889   t_maxNearP1x1 = new std::vector<double>();
0890   t_cHadronEne1x1 = new std::vector<double>();
0891   t_cHadronEne1x1_1 = new std::vector<double>();
0892   t_cHadronEne1x1_2 = new std::vector<double>();
0893   t_cHadronEne1x1_3 = new std::vector<double>();
0894   t_nHadronEne1x1 = new std::vector<double>();
0895   t_photonEne1x1 = new std::vector<double>();
0896   t_eleEne1x1 = new std::vector<double>();
0897   t_muEne1x1 = new std::vector<double>();
0898 
0899   t_maxNearPHC1x1 = new std::vector<double>();
0900   t_cHadronEneHC1x1 = new std::vector<double>();
0901   t_cHadronEneHC1x1_1 = new std::vector<double>();
0902   t_cHadronEneHC1x1_2 = new std::vector<double>();
0903   t_cHadronEneHC1x1_3 = new std::vector<double>();
0904   t_nHadronEneHC1x1 = new std::vector<double>();
0905   t_photonEneHC1x1 = new std::vector<double>();
0906   t_eleEneHC1x1 = new std::vector<double>();
0907   t_muEneHC1x1 = new std::vector<double>();
0908 
0909   t_maxNearPHC3x3 = new std::vector<double>();
0910   t_cHadronEneHC3x3 = new std::vector<double>();
0911   t_cHadronEneHC3x3_1 = new std::vector<double>();
0912   t_cHadronEneHC3x3_2 = new std::vector<double>();
0913   t_cHadronEneHC3x3_3 = new std::vector<double>();
0914   t_nHadronEneHC3x3 = new std::vector<double>();
0915   t_photonEneHC3x3 = new std::vector<double>();
0916   t_eleEneHC3x3 = new std::vector<double>();
0917   t_muEneHC3x3 = new std::vector<double>();
0918 
0919   t_maxNearPHC5x5 = new std::vector<double>();
0920   t_cHadronEneHC5x5 = new std::vector<double>();
0921   t_cHadronEneHC5x5_1 = new std::vector<double>();
0922   t_cHadronEneHC5x5_2 = new std::vector<double>();
0923   t_cHadronEneHC5x5_3 = new std::vector<double>();
0924   t_nHadronEneHC5x5 = new std::vector<double>();
0925   t_photonEneHC5x5 = new std::vector<double>();
0926   t_eleEneHC5x5 = new std::vector<double>();
0927   t_muEneHC5x5 = new std::vector<double>();
0928 
0929   t_maxNearPHC7x7 = new std::vector<double>();
0930   t_cHadronEneHC7x7 = new std::vector<double>();
0931   t_cHadronEneHC7x7_1 = new std::vector<double>();
0932   t_cHadronEneHC7x7_2 = new std::vector<double>();
0933   t_cHadronEneHC7x7_3 = new std::vector<double>();
0934   t_nHadronEneHC7x7 = new std::vector<double>();
0935   t_photonEneHC7x7 = new std::vector<double>();
0936   t_eleEneHC7x7 = new std::vector<double>();
0937   t_muEneHC7x7 = new std::vector<double>();
0938 
0939   t_maxNearPR = new std::vector<double>();
0940   t_cHadronEneR = new std::vector<double>();
0941   t_cHadronEneR_1 = new std::vector<double>();
0942   t_cHadronEneR_2 = new std::vector<double>();
0943   t_cHadronEneR_3 = new std::vector<double>();
0944   t_nHadronEneR = new std::vector<double>();
0945   t_photonEneR = new std::vector<double>();
0946   t_eleEneR = new std::vector<double>();
0947   t_muEneR = new std::vector<double>();
0948 
0949   t_maxNearPIsoR = new std::vector<double>();
0950   t_cHadronEneIsoR = new std::vector<double>();
0951   t_cHadronEneIsoR_1 = new std::vector<double>();
0952   t_cHadronEneIsoR_2 = new std::vector<double>();
0953   t_cHadronEneIsoR_3 = new std::vector<double>();
0954   t_nHadronEneIsoR = new std::vector<double>();
0955   t_photonEneIsoR = new std::vector<double>();
0956   t_eleEneIsoR = new std::vector<double>();
0957   t_muEneIsoR = new std::vector<double>();
0958 
0959   t_maxNearPHCR = new std::vector<double>();
0960   t_cHadronEneHCR = new std::vector<double>();
0961   t_cHadronEneHCR_1 = new std::vector<double>();
0962   t_cHadronEneHCR_2 = new std::vector<double>();
0963   t_cHadronEneHCR_3 = new std::vector<double>();
0964   t_nHadronEneHCR = new std::vector<double>();
0965   t_photonEneHCR = new std::vector<double>();
0966   t_eleEneHCR = new std::vector<double>();
0967   t_muEneHCR = new std::vector<double>();
0968 
0969   t_maxNearPIsoHCR = new std::vector<double>();
0970   t_cHadronEneIsoHCR = new std::vector<double>();
0971   t_cHadronEneIsoHCR_1 = new std::vector<double>();
0972   t_cHadronEneIsoHCR_2 = new std::vector<double>();
0973   t_cHadronEneIsoHCR_3 = new std::vector<double>();
0974   t_nHadronEneIsoHCR = new std::vector<double>();
0975   t_photonEneIsoHCR = new std::vector<double>();
0976   t_eleEneIsoHCR = new std::vector<double>();
0977   t_muEneIsoHCR = new std::vector<double>();
0978 
0979   tree_->Branch("t_isoTrkPAll", "std::vector<double>", &t_isoTrkPAll);
0980   tree_->Branch("t_isoTrkPtAll", "std::vector<double>", &t_isoTrkPtAll);
0981   tree_->Branch("t_isoTrkPhiAll", "std::vector<double>", &t_isoTrkPhiAll);
0982   tree_->Branch("t_isoTrkEtaAll", "std::vector<double>", &t_isoTrkEtaAll);
0983   tree_->Branch("t_isoTrkDPhiAll", "std::vector<double>", &t_isoTrkDPhiAll);
0984   tree_->Branch("t_isoTrkDEtaAll", "std::vector<double>", &t_isoTrkDEtaAll);
0985   tree_->Branch("t_isoTrkPdgIdAll", "std::vector<double>", &t_isoTrkPdgIdAll);
0986 
0987   tree_->Branch("t_isoTrkP", "std::vector<double>", &t_isoTrkP);
0988   tree_->Branch("t_isoTrkPt", "std::vector<double>", &t_isoTrkPt);
0989   tree_->Branch("t_isoTrkEne", "std::vector<double>", &t_isoTrkEne);
0990   tree_->Branch("t_isoTrkEta", "std::vector<double>", &t_isoTrkEta);
0991   tree_->Branch("t_isoTrkPhi", "std::vector<double>", &t_isoTrkPhi);
0992   tree_->Branch("t_isoTrkEtaEC", "std::vector<double>", &t_isoTrkEtaEC);
0993   tree_->Branch("t_isoTrkPhiEC", "std::vector<double>", &t_isoTrkPhiEC);
0994   tree_->Branch("t_isoTrkPdgId", "std::vector<double>", &t_isoTrkPdgId);
0995 
0996   tree_->Branch("t_maxNearP31x31", "std::vector<double>", &t_maxNearP31x31);
0997   tree_->Branch("t_cHadronEne31x31", "std::vector<double>", &t_cHadronEne31x31);
0998   tree_->Branch("t_cHadronEne31x31_1", "std::vector<double>", &t_cHadronEne31x31_1);
0999   tree_->Branch("t_cHadronEne31x31_2", "std::vector<double>", &t_cHadronEne31x31_2);
1000   tree_->Branch("t_cHadronEne31x31_3", "std::vector<double>", &t_cHadronEne31x31_3);
1001   tree_->Branch("t_nHadronEne31x31", "std::vector<double>", &t_nHadronEne31x31);
1002   tree_->Branch("t_photonEne31x31", "std::vector<double>", &t_photonEne31x31);
1003   tree_->Branch("t_eleEne31x31", "std::vector<double>", &t_eleEne31x31);
1004   tree_->Branch("t_muEne31x31", "std::vector<double>", &t_muEne31x31);
1005 
1006   tree_->Branch("t_maxNearP25x25", "std::vector<double>", &t_maxNearP25x25);
1007   tree_->Branch("t_cHadronEne25x25", "std::vector<double>", &t_cHadronEne25x25);
1008   tree_->Branch("t_cHadronEne25x25_1", "std::vector<double>", &t_cHadronEne25x25_1);
1009   tree_->Branch("t_cHadronEne25x25_2", "std::vector<double>", &t_cHadronEne25x25_2);
1010   tree_->Branch("t_cHadronEne25x25_3", "std::vector<double>", &t_cHadronEne25x25_3);
1011   tree_->Branch("t_nHadronEne25x25", "std::vector<double>", &t_nHadronEne25x25);
1012   tree_->Branch("t_photonEne25x25", "std::vector<double>", &t_photonEne25x25);
1013   tree_->Branch("t_eleEne25x25", "std::vector<double>", &t_eleEne25x25);
1014   tree_->Branch("t_muEne25x25", "std::vector<double>", &t_muEne25x25);
1015 
1016   tree_->Branch("t_maxNearP21x21", "std::vector<double>", &t_maxNearP21x21);
1017   tree_->Branch("t_cHadronEne21x21", "std::vector<double>", &t_cHadronEne21x21);
1018   tree_->Branch("t_cHadronEne21x21_1", "std::vector<double>", &t_cHadronEne21x21_1);
1019   tree_->Branch("t_cHadronEne21x21_2", "std::vector<double>", &t_cHadronEne21x21_2);
1020   tree_->Branch("t_cHadronEne21x21_3", "std::vector<double>", &t_cHadronEne21x21_3);
1021   tree_->Branch("t_nHadronEne21x21", "std::vector<double>", &t_nHadronEne21x21);
1022   tree_->Branch("t_photonEne21x21", "std::vector<double>", &t_photonEne21x21);
1023   tree_->Branch("t_eleEne21x21", "std::vector<double>", &t_eleEne21x21);
1024   tree_->Branch("t_muEne21x21", "std::vector<double>", &t_muEne21x21);
1025 
1026   tree_->Branch("t_maxNearP15x15", "std::vector<double>", &t_maxNearP15x15);
1027   tree_->Branch("t_cHadronEne15x15", "std::vector<double>", &t_cHadronEne15x15);
1028   tree_->Branch("t_cHadronEne15x15_1", "std::vector<double>", &t_cHadronEne15x15_1);
1029   tree_->Branch("t_cHadronEne15x15_2", "std::vector<double>", &t_cHadronEne15x15_2);
1030   tree_->Branch("t_cHadronEne15x15_3", "std::vector<double>", &t_cHadronEne15x15_3);
1031   tree_->Branch("t_nHadronEne15x15", "std::vector<double>", &t_nHadronEne15x15);
1032   tree_->Branch("t_photonEne15x15", "std::vector<double>", &t_photonEne15x15);
1033   tree_->Branch("t_eleEne15x15", "std::vector<double>", &t_eleEne15x15);
1034   tree_->Branch("t_muEne15x15", "std::vector<double>", &t_muEne15x15);
1035 
1036   tree_->Branch("t_maxNearP11x11", "std::vector<double>", &t_maxNearP11x11);
1037   tree_->Branch("t_cHadronEne11x11", "std::vector<double>", &t_cHadronEne11x11);
1038   tree_->Branch("t_cHadronEne11x11_1", "std::vector<double>", &t_cHadronEne11x11_1);
1039   tree_->Branch("t_cHadronEne11x11_2", "std::vector<double>", &t_cHadronEne11x11_2);
1040   tree_->Branch("t_cHadronEne11x11_3", "std::vector<double>", &t_cHadronEne11x11_3);
1041   tree_->Branch("t_nHadronEne11x11", "std::vector<double>", &t_nHadronEne11x11);
1042   tree_->Branch("t_photonEne11x11", "std::vector<double>", &t_photonEne11x11);
1043   tree_->Branch("t_eleEne11x11", "std::vector<double>", &t_eleEne11x11);
1044   tree_->Branch("t_muEne11x11", "std::vector<double>", &t_muEne11x11);
1045 
1046   tree_->Branch("t_maxNearP9x9", "std::vector<double>", &t_maxNearP9x9);
1047   tree_->Branch("t_cHadronEne9x9", "std::vector<double>", &t_cHadronEne9x9);
1048   tree_->Branch("t_cHadronEne9x9_1", "std::vector<double>", &t_cHadronEne9x9_1);
1049   tree_->Branch("t_cHadronEne9x9_2", "std::vector<double>", &t_cHadronEne9x9_2);
1050   tree_->Branch("t_cHadronEne9x9_3", "std::vector<double>", &t_cHadronEne9x9_3);
1051   tree_->Branch("t_nHadronEne9x9", "std::vector<double>", &t_nHadronEne9x9);
1052   tree_->Branch("t_photonEne9x9", "std::vector<double>", &t_photonEne9x9);
1053   tree_->Branch("t_eleEne9x9", "std::vector<double>", &t_eleEne9x9);
1054   tree_->Branch("t_muEne9x9", "std::vector<double>", &t_muEne9x9);
1055 
1056   tree_->Branch("t_maxNearP7x7", "std::vector<double>", &t_maxNearP7x7);
1057   tree_->Branch("t_cHadronEne7x7", "std::vector<double>", &t_cHadronEne7x7);
1058   tree_->Branch("t_cHadronEne7x7_1", "std::vector<double>", &t_cHadronEne7x7_1);
1059   tree_->Branch("t_cHadronEne7x7_2", "std::vector<double>", &t_cHadronEne7x7_2);
1060   tree_->Branch("t_cHadronEne7x7_3", "std::vector<double>", &t_cHadronEne7x7_3);
1061   tree_->Branch("t_nHadronEne7x7", "std::vector<double>", &t_nHadronEne7x7);
1062   tree_->Branch("t_photonEne7x7", "std::vector<double>", &t_photonEne7x7);
1063   tree_->Branch("t_eleEne7x7", "std::vector<double>", &t_eleEne7x7);
1064   tree_->Branch("t_muEne7x7", "std::vector<double>", &t_muEne7x7);
1065 
1066   tree_->Branch("t_maxNearP3x3", "std::vector<double>", &t_maxNearP3x3);
1067   tree_->Branch("t_cHadronEne3x3", "std::vector<double>", &t_cHadronEne3x3);
1068   tree_->Branch("t_cHadronEne3x3_1", "std::vector<double>", &t_cHadronEne3x3_1);
1069   tree_->Branch("t_cHadronEne3x3_2", "std::vector<double>", &t_cHadronEne3x3_2);
1070   tree_->Branch("t_cHadronEne3x3_3", "std::vector<double>", &t_cHadronEne3x3_3);
1071   tree_->Branch("t_nHadronEne3x3", "std::vector<double>", &t_nHadronEne3x3);
1072   tree_->Branch("t_photonEne3x3", "std::vector<double>", &t_photonEne3x3);
1073   tree_->Branch("t_eleEne3x3", "std::vector<double>", &t_eleEne3x3);
1074   tree_->Branch("t_muEne3x3", "std::vector<double>", &t_muEne3x3);
1075 
1076   tree_->Branch("t_maxNearP1x1", "std::vector<double>", &t_maxNearP1x1);
1077   tree_->Branch("t_cHadronEne1x1", "std::vector<double>", &t_cHadronEne1x1);
1078   tree_->Branch("t_cHadronEne1x1_1", "std::vector<double>", &t_cHadronEne1x1_1);
1079   tree_->Branch("t_cHadronEne1x1_2", "std::vector<double>", &t_cHadronEne1x1_2);
1080   tree_->Branch("t_cHadronEne1x1_3", "std::vector<double>", &t_cHadronEne1x1_3);
1081   tree_->Branch("t_nHadronEne1x1", "std::vector<double>", &t_nHadronEne1x1);
1082   tree_->Branch("t_photonEne1x1", "std::vector<double>", &t_photonEne1x1);
1083   tree_->Branch("t_eleEne1x1", "std::vector<double>", &t_eleEne1x1);
1084   tree_->Branch("t_muEne1x1", "std::vector<double>", &t_muEne1x1);
1085 
1086   tree_->Branch("t_maxNearPHC1x1", "std::vector<double>", &t_maxNearPHC1x1);
1087   tree_->Branch("t_cHadronEneHC1x1", "std::vector<double>", &t_cHadronEneHC1x1);
1088   tree_->Branch("t_cHadronEneHC1x1_1", "std::vector<double>", &t_cHadronEneHC1x1_1);
1089   tree_->Branch("t_cHadronEneHC1x1_2", "std::vector<double>", &t_cHadronEneHC1x1_2);
1090   tree_->Branch("t_cHadronEneHC1x1_3", "std::vector<double>", &t_cHadronEneHC1x1_3);
1091   tree_->Branch("t_nHadronEneHC1x1", "std::vector<double>", &t_nHadronEneHC1x1);
1092   tree_->Branch("t_photonEneHC1x1", "std::vector<double>", &t_photonEneHC1x1);
1093   tree_->Branch("t_eleEneHC1x1", "std::vector<double>", &t_eleEneHC1x1);
1094   tree_->Branch("t_muEneHC1x1", "std::vector<double>", &t_muEneHC1x1);
1095 
1096   tree_->Branch("t_maxNearPHC3x3", "std::vector<double>", &t_maxNearPHC3x3);
1097   tree_->Branch("t_cHadronEneHC3x3", "std::vector<double>", &t_cHadronEneHC3x3);
1098   tree_->Branch("t_cHadronEneHC3x3_1", "std::vector<double>", &t_cHadronEneHC3x3_1);
1099   tree_->Branch("t_cHadronEneHC3x3_2", "std::vector<double>", &t_cHadronEneHC3x3_2);
1100   tree_->Branch("t_cHadronEneHC3x3_3", "std::vector<double>", &t_cHadronEneHC3x3_3);
1101   tree_->Branch("t_nHadronEneHC3x3", "std::vector<double>", &t_nHadronEneHC3x3);
1102   tree_->Branch("t_photonEneHC3x3", "std::vector<double>", &t_photonEneHC3x3);
1103   tree_->Branch("t_eleEneHC3x3", "std::vector<double>", &t_eleEneHC3x3);
1104   tree_->Branch("t_muEneHC3x3", "std::vector<double>", &t_muEneHC3x3);
1105 
1106   tree_->Branch("t_maxNearPHC5x5", "std::vector<double>", &t_maxNearPHC5x5);
1107   tree_->Branch("t_cHadronEneHC5x5", "std::vector<double>", &t_cHadronEneHC5x5);
1108   tree_->Branch("t_cHadronEneHC5x5_1", "std::vector<double>", &t_cHadronEneHC5x5_1);
1109   tree_->Branch("t_cHadronEneHC5x5_2", "std::vector<double>", &t_cHadronEneHC5x5_2);
1110   tree_->Branch("t_cHadronEneHC5x5_3", "std::vector<double>", &t_cHadronEneHC5x5_3);
1111   tree_->Branch("t_nHadronEneHC5x5", "std::vector<double>", &t_nHadronEneHC5x5);
1112   tree_->Branch("t_photonEneHC5x5", "std::vector<double>", &t_photonEneHC5x5);
1113   tree_->Branch("t_eleEneHC5x5", "std::vector<double>", &t_eleEneHC5x5);
1114   tree_->Branch("t_muEneHC5x5", "std::vector<double>", &t_muEneHC5x5);
1115 
1116   tree_->Branch("t_maxNearPHC7x7", "std::vector<double>", &t_maxNearPHC7x7);
1117   tree_->Branch("t_cHadronEneHC7x7", "std::vector<double>", &t_cHadronEneHC7x7);
1118   tree_->Branch("t_cHadronEneHC7x7_1", "std::vector<double>", &t_cHadronEneHC7x7_1);
1119   tree_->Branch("t_cHadronEneHC7x7_2", "std::vector<double>", &t_cHadronEneHC7x7_2);
1120   tree_->Branch("t_cHadronEneHC7x7_3", "std::vector<double>", &t_cHadronEneHC7x7_3);
1121   tree_->Branch("t_nHadronEneHC7x7", "std::vector<double>", &t_nHadronEneHC7x7);
1122   tree_->Branch("t_photonEneHC7x7", "std::vector<double>", &t_photonEneHC7x7);
1123   tree_->Branch("t_eleEneHC7x7", "std::vector<double>", &t_eleEneHC7x7);
1124   tree_->Branch("t_muEneHC7x7", "std::vector<double>", &t_muEneHC7x7);
1125 
1126   tree_->Branch("t_maxNearPR", "std::vector<double>", &t_maxNearPR);
1127   tree_->Branch("t_cHadronEneR", "std::vector<double>", &t_cHadronEneR);
1128   tree_->Branch("t_cHadronEneR_1", "std::vector<double>", &t_cHadronEneR_1);
1129   tree_->Branch("t_cHadronEneR_2", "std::vector<double>", &t_cHadronEneR_2);
1130   tree_->Branch("t_cHadronEneR_3", "std::vector<double>", &t_cHadronEneR_3);
1131   tree_->Branch("t_nHadronEneR", "std::vector<double>", &t_nHadronEneR);
1132   tree_->Branch("t_photonEneR", "std::vector<double>", &t_photonEneR);
1133   tree_->Branch("t_eleEneR", "std::vector<double>", &t_eleEneR);
1134   tree_->Branch("t_muEneR", "std::vector<double>", &t_muEneR);
1135 
1136   tree_->Branch("t_maxNearPIsoR", "std::vector<double>", &t_maxNearPIsoR);
1137   tree_->Branch("t_cHadronEneIsoR", "std::vector<double>", &t_cHadronEneIsoR);
1138   tree_->Branch("t_cHadronEneIsoR_1", "std::vector<double>", &t_cHadronEneIsoR_1);
1139   tree_->Branch("t_cHadronEneIsoR_2", "std::vector<double>", &t_cHadronEneIsoR_2);
1140   tree_->Branch("t_cHadronEneIsoR_3", "std::vector<double>", &t_cHadronEneIsoR_3);
1141   tree_->Branch("t_nHadronEneIsoR", "std::vector<double>", &t_nHadronEneIsoR);
1142   tree_->Branch("t_photonEneIsoR", "std::vector<double>", &t_photonEneIsoR);
1143   tree_->Branch("t_eleEneIsoR", "std::vector<double>", &t_eleEneIsoR);
1144   tree_->Branch("t_muEneIsoR", "std::vector<double>", &t_muEneIsoR);
1145 
1146   tree_->Branch("t_maxNearPHCR", "std::vector<double>", &t_maxNearPHCR);
1147   tree_->Branch("t_cHadronEneHCR", "std::vector<double>", &t_cHadronEneHCR);
1148   tree_->Branch("t_cHadronEneHCR_1", "std::vector<double>", &t_cHadronEneHCR_1);
1149   tree_->Branch("t_cHadronEneHCR_2", "std::vector<double>", &t_cHadronEneHCR_2);
1150   tree_->Branch("t_cHadronEneHCR_3", "std::vector<double>", &t_cHadronEneHCR_3);
1151   tree_->Branch("t_nHadronEneHCR", "std::vector<double>", &t_nHadronEneHCR);
1152   tree_->Branch("t_photonEneHCR", "std::vector<double>", &t_photonEneHCR);
1153   tree_->Branch("t_eleEneHCR", "std::vector<double>", &t_eleEneHCR);
1154   tree_->Branch("t_muEneHCR", "std::vector<double>", &t_muEneHCR);
1155 
1156   tree_->Branch("t_maxNearPIsoHCR", "std::vector<double>", &t_maxNearPIsoHCR);
1157   tree_->Branch("t_cHadronEneIsoHCR", "std::vector<double>", &t_cHadronEneIsoHCR);
1158   tree_->Branch("t_cHadronEneIsoHCR_1", "std::vector<double>", &t_cHadronEneIsoHCR_1);
1159   tree_->Branch("t_cHadronEneIsoHCR_2", "std::vector<double>", &t_cHadronEneIsoHCR_2);
1160   tree_->Branch("t_cHadronEneIsoHCR_3", "std::vector<double>", &t_cHadronEneIsoHCR_3);
1161   tree_->Branch("t_nHadronEneIsoHCR", "std::vector<double>", &t_nHadronEneIsoHCR);
1162   tree_->Branch("t_photonEneIsoHCR", "std::vector<double>", &t_photonEneIsoHCR);
1163   tree_->Branch("t_eleEneIsoHCR", "std::vector<double>", &t_eleEneIsoHCR);
1164   tree_->Branch("t_muEneIsoHCR", "std::vector<double>", &t_muEneIsoHCR);
1165 }
1166 
1167 void StudyCaloGen::clearTreeVectors() {
1168   // t_maxNearP31x31     ->clear();
1169   // t_nEvtProc          ->clear();
1170 
1171   t_isoTrkPAll->clear();
1172   t_isoTrkPtAll->clear();
1173   t_isoTrkPhiAll->clear();
1174   t_isoTrkEtaAll->clear();
1175   t_isoTrkDPhiAll->clear();
1176   t_isoTrkDEtaAll->clear();
1177   t_isoTrkPdgIdAll->clear();
1178 
1179   t_isoTrkP->clear();
1180   t_isoTrkPt->clear();
1181   t_isoTrkEne->clear();
1182   t_isoTrkEta->clear();
1183   t_isoTrkPhi->clear();
1184   t_isoTrkEtaEC->clear();
1185   t_isoTrkPhiEC->clear();
1186   t_isoTrkPdgId->clear();
1187 
1188   t_maxNearP31x31->clear();
1189   t_cHadronEne31x31->clear();
1190   t_cHadronEne31x31_1->clear();
1191   t_cHadronEne31x31_2->clear();
1192   t_cHadronEne31x31_3->clear();
1193   t_nHadronEne31x31->clear();
1194   t_photonEne31x31->clear();
1195   t_eleEne31x31->clear();
1196   t_muEne31x31->clear();
1197 
1198   t_maxNearP25x25->clear();
1199   t_cHadronEne25x25->clear();
1200   t_cHadronEne25x25_1->clear();
1201   t_cHadronEne25x25_2->clear();
1202   t_cHadronEne25x25_3->clear();
1203   t_nHadronEne25x25->clear();
1204   t_photonEne25x25->clear();
1205   t_eleEne25x25->clear();
1206   t_muEne25x25->clear();
1207 
1208   t_maxNearP21x21->clear();
1209   t_cHadronEne21x21->clear();
1210   t_cHadronEne21x21_1->clear();
1211   t_cHadronEne21x21_2->clear();
1212   t_cHadronEne21x21_3->clear();
1213   t_nHadronEne21x21->clear();
1214   t_photonEne21x21->clear();
1215   t_eleEne21x21->clear();
1216   t_muEne21x21->clear();
1217 
1218   t_maxNearP15x15->clear();
1219   t_cHadronEne15x15->clear();
1220   t_cHadronEne15x15_1->clear();
1221   t_cHadronEne15x15_2->clear();
1222   t_cHadronEne15x15_3->clear();
1223   t_nHadronEne15x15->clear();
1224   t_photonEne15x15->clear();
1225   t_eleEne15x15->clear();
1226   t_muEne15x15->clear();
1227 
1228   t_maxNearP11x11->clear();
1229   t_cHadronEne11x11->clear();
1230   t_cHadronEne11x11_1->clear();
1231   t_cHadronEne11x11_2->clear();
1232   t_cHadronEne11x11_3->clear();
1233   t_nHadronEne11x11->clear();
1234   t_photonEne11x11->clear();
1235   t_eleEne11x11->clear();
1236   t_muEne11x11->clear();
1237 
1238   t_maxNearP9x9->clear();
1239   t_cHadronEne9x9->clear();
1240   t_cHadronEne9x9_1->clear();
1241   t_cHadronEne9x9_2->clear();
1242   t_cHadronEne9x9_3->clear();
1243   t_nHadronEne9x9->clear();
1244   t_photonEne9x9->clear();
1245   t_eleEne9x9->clear();
1246   t_muEne9x9->clear();
1247 
1248   t_maxNearP7x7->clear();
1249   t_cHadronEne7x7->clear();
1250   t_cHadronEne7x7_1->clear();
1251   t_cHadronEne7x7_2->clear();
1252   t_cHadronEne7x7_3->clear();
1253   t_nHadronEne7x7->clear();
1254   t_photonEne7x7->clear();
1255   t_eleEne7x7->clear();
1256   t_muEne7x7->clear();
1257 
1258   t_maxNearP3x3->clear();
1259   t_cHadronEne3x3->clear();
1260   t_cHadronEne3x3_1->clear();
1261   t_cHadronEne3x3_2->clear();
1262   t_cHadronEne3x3_3->clear();
1263   t_nHadronEne3x3->clear();
1264   t_photonEne3x3->clear();
1265   t_eleEne3x3->clear();
1266   t_muEne3x3->clear();
1267 
1268   t_maxNearP1x1->clear();
1269   t_cHadronEne1x1->clear();
1270   t_cHadronEne1x1_1->clear();
1271   t_cHadronEne1x1_2->clear();
1272   t_cHadronEne1x1_3->clear();
1273   t_nHadronEne1x1->clear();
1274   t_photonEne1x1->clear();
1275   t_eleEne1x1->clear();
1276   t_muEne1x1->clear();
1277 
1278   t_maxNearPHC1x1->clear();
1279   t_cHadronEneHC1x1->clear();
1280   t_cHadronEneHC1x1_1->clear();
1281   t_cHadronEneHC1x1_2->clear();
1282   t_cHadronEneHC1x1_3->clear();
1283   t_nHadronEneHC1x1->clear();
1284   t_photonEneHC1x1->clear();
1285   t_eleEneHC1x1->clear();
1286   t_muEneHC1x1->clear();
1287 
1288   t_maxNearPHC3x3->clear();
1289   t_cHadronEneHC3x3->clear();
1290   t_cHadronEneHC3x3_1->clear();
1291   t_cHadronEneHC3x3_2->clear();
1292   t_cHadronEneHC3x3_3->clear();
1293   t_nHadronEneHC3x3->clear();
1294   t_photonEneHC3x3->clear();
1295   t_eleEneHC3x3->clear();
1296   t_muEneHC3x3->clear();
1297 
1298   t_maxNearPHC5x5->clear();
1299   t_cHadronEneHC5x5->clear();
1300   t_cHadronEneHC5x5_1->clear();
1301   t_cHadronEneHC5x5_2->clear();
1302   t_cHadronEneHC5x5_3->clear();
1303   t_nHadronEneHC5x5->clear();
1304   t_photonEneHC5x5->clear();
1305   t_eleEneHC5x5->clear();
1306   t_muEneHC5x5->clear();
1307 
1308   t_maxNearPHC7x7->clear();
1309   t_cHadronEneHC7x7->clear();
1310   t_cHadronEneHC7x7_1->clear();
1311   t_cHadronEneHC7x7_2->clear();
1312   t_cHadronEneHC7x7_3->clear();
1313   t_nHadronEneHC7x7->clear();
1314   t_photonEneHC7x7->clear();
1315   t_eleEneHC7x7->clear();
1316   t_muEneHC7x7->clear();
1317 
1318   t_maxNearPR->clear();
1319   t_cHadronEneR->clear();
1320   t_cHadronEneR_1->clear();
1321   t_cHadronEneR_2->clear();
1322   t_cHadronEneR_3->clear();
1323   t_nHadronEneR->clear();
1324   t_photonEneR->clear();
1325   t_eleEneR->clear();
1326   t_muEneR->clear();
1327 
1328   t_maxNearPIsoR->clear();
1329   t_cHadronEneIsoR->clear();
1330   t_cHadronEneIsoR_1->clear();
1331   t_cHadronEneIsoR_2->clear();
1332   t_cHadronEneIsoR_3->clear();
1333   t_nHadronEneIsoR->clear();
1334   t_photonEneIsoR->clear();
1335   t_eleEneIsoR->clear();
1336   t_muEneIsoR->clear();
1337 
1338   t_maxNearPHCR->clear();
1339   t_cHadronEneHCR->clear();
1340   t_cHadronEneHCR_1->clear();
1341   t_cHadronEneHCR_2->clear();
1342   t_cHadronEneHCR_3->clear();
1343   t_nHadronEneHCR->clear();
1344   t_photonEneHCR->clear();
1345   t_eleEneHCR->clear();
1346   t_muEneHCR->clear();
1347 
1348   t_maxNearPIsoHCR->clear();
1349   t_cHadronEneIsoHCR->clear();
1350   t_cHadronEneIsoHCR_1->clear();
1351   t_cHadronEneIsoHCR_2->clear();
1352   t_cHadronEneIsoHCR_3->clear();
1353   t_nHadronEneIsoHCR->clear();
1354   t_photonEneIsoHCR->clear();
1355   t_eleEneIsoHCR->clear();
1356   t_muEneIsoHCR->clear();
1357 }
1358 
1359 int StudyCaloGen::particleCode(int pdgId) {
1360   int partID[Particles] = {11, -11, 21, 211, -211, 321, -321, 2212, 2112, -2212, -2112, 130};
1361   int ix = -1;
1362   for (int ik = 0; ik < Particles; ++ik) {
1363     if (pdgId == partID[ik]) {
1364       ix = ik;
1365       break;
1366     }
1367   }
1368   return ix;
1369 }
1370 
1371 //define this as a plug-in
1372 DEFINE_FWK_MODULE(StudyCaloGen);