CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_10_patch1/src/Calibration/IsolatedParticles/plugins/IsolatedTracksHcalScale.cc

Go to the documentation of this file.
00001 
00002 #include "DataFormats/TrackReco/interface/TrackBase.h"
00003 #include "Calibration/IsolatedParticles/plugins/IsolatedTracksHcalScale.h"
00004 
00005 #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h"
00006 #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h"
00007 #include "Calibration/IsolatedParticles/interface/eCone.h"
00008 #include "Calibration/IsolatedParticles/interface/eHCALMatrixExtra.h"
00009 #include "MagneticField/Engine/interface/MagneticField.h"
00010 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
00011 #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h"
00012 #include "Geometry/CaloTopology/interface/EcalTrigTowerConstituentsMap.h"
00013 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
00014 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
00015 
00016 IsolatedTracksHcalScale::IsolatedTracksHcalScale(const edm::ParameterSet& iConfig) {
00017 
00018   //now do what ever initialization is needed
00019   doMC                                = iConfig.getUntrackedParameter<bool>("DoMC", false); 
00020   myverbose                           = iConfig.getUntrackedParameter<int>("Verbosity", 5          );
00021   theTrackQuality                     = iConfig.getUntrackedParameter<std::string>("TrackQuality","highPurity");
00022   reco::TrackBase::TrackQuality trackQuality_=reco::TrackBase::qualityByName(theTrackQuality);
00023   selectionParameters.minPt           = iConfig.getUntrackedParameter<double>("MinTrackPt", 10.0);
00024   selectionParameters.minQuality      = trackQuality_;
00025   selectionParameters.maxDxyPV        = iConfig.getUntrackedParameter<double>("MaxDxyPV", 0.2);
00026   selectionParameters.maxDzPV         = iConfig.getUntrackedParameter<double>("MaxDzPV",  5.0);
00027   selectionParameters.maxChi2         = iConfig.getUntrackedParameter<double>("MaxChi2",  5.0);
00028   selectionParameters.maxDpOverP      = iConfig.getUntrackedParameter<double>("MaxDpOverP",  0.1);
00029   selectionParameters.minOuterHit     = iConfig.getUntrackedParameter<int>("MinOuterHit", 4);
00030   selectionParameters.minLayerCrossed = iConfig.getUntrackedParameter<int>("MinLayerCrossed", 8);
00031   selectionParameters.maxInMiss       = iConfig.getUntrackedParameter<int>("MaxInMiss", 0);
00032   selectionParameters.maxOutMiss      = iConfig.getUntrackedParameter<int>("MaxOutMiss", 0);
00033   a_coneR                             = iConfig.getUntrackedParameter<double>("ConeRadius",34.98);
00034   a_charIsoR                          = a_coneR + 28.9;
00035   a_neutIsoR                          = a_charIsoR*0.726;
00036   a_mipR                              = iConfig.getUntrackedParameter<double>("ConeRadiusMIP",14.0);
00037   tMinE_                              = iConfig.getUntrackedParameter<double>("TimeMinCutECAL", -500.);
00038   tMaxE_                              = iConfig.getUntrackedParameter<double>("TimeMaxCutECAL",  500.);
00039   
00040   if (myverbose>=0) {
00041     std::cout <<"Parameters read from config file \n" 
00042               <<" doMC "              << doMC
00043               <<"\t myverbose "       << myverbose        
00044               <<"\t minPt "           << selectionParameters.minPt   
00045               <<"\t theTrackQuality " << theTrackQuality
00046               <<"\t minQuality "      << selectionParameters.minQuality
00047               <<"\t maxDxyPV "        << selectionParameters.maxDxyPV          
00048               <<"\t maxDzPV "         << selectionParameters.maxDzPV          
00049               <<"\t maxChi2 "         << selectionParameters.maxChi2          
00050               <<"\t maxDpOverP "      << selectionParameters.maxDpOverP
00051               <<"\t minOuterHit "     << selectionParameters.minOuterHit
00052               <<"\t minLayerCrossed " << selectionParameters.minLayerCrossed
00053               <<"\t maxInMiss "       << selectionParameters.maxInMiss
00054               <<"\t maxOutMiss "      << selectionParameters.maxOutMiss
00055               <<"\t a_coneR "         << a_coneR          
00056               <<"\t a_charIsoR "      << a_charIsoR          
00057               <<"\t a_neutIsoR "      << a_neutIsoR          
00058               <<"\t a_mipR "          << a_mipR 
00059               <<"\t time Range ("     << tMinE_ << ":" << tMaxE_ << ")"
00060               << std::endl;
00061   }
00062   initL1 = false;
00063   
00064 }
00065 
00066 IsolatedTracksHcalScale::~IsolatedTracksHcalScale() {}
00067 
00068 void IsolatedTracksHcalScale::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
00069 
00070   edm::ESHandle<MagneticField> bFieldH;
00071   iSetup.get<IdealMagneticFieldRecord>().get(bFieldH);
00072   bField = bFieldH.product();
00073 
00074   // get handles to calogeometry and calotopology
00075   edm::ESHandle<CaloGeometry> pG;
00076   iSetup.get<CaloGeometryRecord>().get(pG);
00077   const CaloGeometry* geo = pG.product();
00078 
00079   edm::ESHandle<CaloTopology> theCaloTopology;
00080   iSetup.get<CaloTopologyRecord>().get(theCaloTopology); 
00081   const CaloTopology *caloTopology = theCaloTopology.product();
00082   
00083   /*  
00084   edm::ESHandle<HcalTopology> htopo;
00085   iSetup.get<IdealGeometryRecord>().get(htopo);
00086   const HcalTopology* theHBHETopology = htopo.product();
00087   */
00088 
00089   // Retrieve the good/bad ECAL channels from the DB
00090   edm::ESHandle<EcalChannelStatus> ecalChStatus;
00091   iSetup.get<EcalChannelStatusRcd>().get(ecalChStatus);
00092   const EcalChannelStatus* theEcalChStatus = ecalChStatus.product();
00093 
00094   /*  
00095   // Retrieve trigger tower map
00096   edm::ESHandle<EcalTrigTowerConstituentsMap> hTtmap;
00097   iSetup.get<IdealGeometryRecord>().get(hTtmap);
00098   const EcalTrigTowerConstituentsMap& ttMap = *hTtmap;
00099   */
00100 
00101   clearTreeVectors();
00102 
00103   nEventProc++;
00104 
00105   t_RunNo = iEvent.id().run();
00106   t_EvtNo = iEvent.id().event();
00107   t_Lumi  = iEvent.luminosityBlock();
00108   t_Bunch = iEvent.bunchCrossing();
00109   if (myverbose>0) std::cout << nEventProc << " Run " << t_RunNo << " Event " << t_EvtNo << " Lumi " << t_Lumi << " Bunch " << t_Bunch << std::endl;
00110 
00111   edm::Handle<reco::TrackCollection> trkCollection;
00112   iEvent.getByLabel("generalTracks", trkCollection);
00113 
00114   edm::Handle<reco::VertexCollection> recVtxs;
00115   iEvent.getByLabel("offlinePrimaryVertices",recVtxs);  
00116 
00117   // Get the beamspot
00118   edm::Handle<reco::BeamSpot> beamSpotH;
00119   iEvent.getByLabel("offlineBeamSpot", beamSpotH);
00120 
00121   math::XYZPoint leadPV(0,0,0);
00122   if (recVtxs->size()>0 && !((*recVtxs)[0].isFake())) {
00123     leadPV = math::XYZPoint( (*recVtxs)[0].x(),(*recVtxs)[0].y(), (*recVtxs)[0].z() );
00124   } else if (beamSpotH.isValid()) {
00125     leadPV = beamSpotH->position();
00126   }
00127 
00128   if (myverbose>0) {
00129     std::cout << "Primary Vertex " << leadPV;
00130     if (beamSpotH.isValid()) std::cout << " Beam Spot " << beamSpotH->position();
00131     std::cout << std::endl;
00132   }
00133 
00134   std::vector<spr::propagatedTrackDirection> trkCaloDirections;
00135   spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDirections, (myverbose>2));
00136   std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
00137   
00138   edm::Handle<EcalRecHitCollection> barrelRecHitsHandle;
00139   edm::Handle<EcalRecHitCollection> endcapRecHitsHandle;
00140   iEvent.getByLabel("ecalRecHit","EcalRecHitsEB",barrelRecHitsHandle);
00141   iEvent.getByLabel("ecalRecHit","EcalRecHitsEE",endcapRecHitsHandle);
00142 
00143   edm::Handle<HBHERecHitCollection> hbhe;
00144   iEvent.getByLabel("hbhereco",hbhe);
00145   const HBHERecHitCollection Hithbhe = *(hbhe.product());
00146 
00147   //get Handles to SimTracks and SimHits
00148   edm::Handle<edm::SimTrackContainer> SimTk;
00149   edm::SimTrackContainer::const_iterator simTrkItr;
00150   edm::Handle<edm::SimVertexContainer> SimVtx;
00151   edm::SimVertexContainer::const_iterator vtxItr = SimVtx->begin();
00152 
00153   //get Handles to PCaloHitContainers of eb/ee/hbhe
00154   edm::Handle<edm::PCaloHitContainer> pcaloeb;
00155   edm::Handle<edm::PCaloHitContainer> pcaloee;
00156   edm::Handle<edm::PCaloHitContainer> pcalohh;
00157 
00158   //associates tracker rechits/simhits to a track
00159   TrackerHitAssociator* associate=0;
00160  
00161   if (doMC) {
00162     iEvent.getByLabel("g4SimHits",SimTk);
00163     iEvent.getByLabel("g4SimHits",SimVtx);
00164     iEvent.getByLabel("g4SimHits", "EcalHitsEB", pcaloeb);
00165     iEvent.getByLabel("g4SimHits", "EcalHitsEE", pcaloee);
00166     iEvent.getByLabel("g4SimHits", "HcalHits", pcalohh);
00167     associate = new TrackerHitAssociator(iEvent);
00168   }
00169  
00170   unsigned int nTracks=0;
00171   for (trkDetItr = trkCaloDirections.begin(),nTracks=0; trkDetItr != trkCaloDirections.end(); trkDetItr++,nTracks++){
00172     const reco::Track* pTrack = &(*(trkDetItr->trkItr));
00173     if (spr::goodTrack(pTrack,leadPV,selectionParameters,(myverbose>2)) && trkDetItr->okECAL && trkDetItr->okHCAL) {
00174       int                nRH_eMipDR=0, nRH_eDR=0, nNearTRKs=0, nRecHitsCone=-99;
00175       double             distFromHotCell=-99.0, distFromHotCell2=-99.0;
00176       int                ietaHotCell=-99, iphiHotCell=-99;
00177       int                ietaHotCell2=-99, iphiHotCell2=-99;
00178       GlobalPoint        gposHotCell(0.,0.,0.), gposHotCell2(0.,0.,0.);
00179       std::vector<DetId> coneRecHitDetIds, coneRecHitDetIds2;
00180       std::pair<double, bool> e11x11_20SigP, e15x15_20SigP;
00181       double hCone = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL,
00182                                      a_coneR, trkDetItr->directionHCAL, nRecHitsCone,
00183                                      coneRecHitDetIds, distFromHotCell,
00184                                      ietaHotCell, iphiHotCell, gposHotCell);
00185       double eHCALDR = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL,
00186                                      a_charIsoR, trkDetItr->directionHCAL, nRecHitsCone,
00187                                      coneRecHitDetIds2, distFromHotCell2,
00188                                      ietaHotCell2, iphiHotCell2, gposHotCell2);
00189 
00190       double conehmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections, a_charIsoR, nNearTRKs, (myverbose>3));
00191 
00192       double eMipDR  = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
00193                                        trkDetItr->pointHCAL, trkDetItr->pointECAL,
00194                                        a_mipR, trkDetItr->directionECAL, nRH_eMipDR);
00195       double eECALDR = spr::eCone_ecal(geo, barrelRecHitsHandle,  endcapRecHitsHandle,
00196                                        trkDetItr->pointHCAL, trkDetItr->pointECAL,
00197                                        a_neutIsoR, trkDetItr->directionECAL, nRH_eDR);
00198       double eMipDR_1 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
00199                                         trkDetItr->pointHCAL, trkDetItr->pointECAL,
00200                                         a_mipR, trkDetItr->directionECAL, nRH_eMipDR,
00201                                         0.030, 0.150);
00202       double eECALDR_1= spr::eCone_ecal(geo, barrelRecHitsHandle,  endcapRecHitsHandle,
00203                                        trkDetItr->pointHCAL, trkDetItr->pointECAL,
00204                                         a_neutIsoR, trkDetItr->directionECAL, nRH_eDR,
00205                                         0.030, 0.150);
00206       double eMipDR_2 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
00207                                         trkDetItr->pointHCAL, trkDetItr->pointECAL,
00208                                         a_mipR, trkDetItr->directionECAL, nRH_eMipDR,
00209                                         0.060, 0.300);
00210       double eECALDR_2= spr::eCone_ecal(geo, barrelRecHitsHandle,  endcapRecHitsHandle,
00211                                        trkDetItr->pointHCAL, trkDetItr->pointECAL,
00212                                         a_neutIsoR, trkDetItr->directionECAL, nRH_eDR,
00213                                         0.060, 0.300);
00214 
00215       HcalDetId closestCell = (HcalDetId)(trkDetItr->detIdHCAL);
00216 
00217       edm::ESHandle<EcalSeverityLevelAlgo> sevlv;
00218       iSetup.get<EcalSeverityLevelAlgoRcd>().get(sevlv);
00219 
00220       e11x11_20SigP = spr::eECALmatrix(trkDetItr->detIdECAL,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),5,5,   0.060,  0.300, tMinE_,tMaxE_);
00221       e15x15_20SigP = spr::eECALmatrix(trkDetItr->detIdECAL,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),7,7,   0.060,  0.300, tMinE_,tMaxE_);
00222       
00223       // Fill the tree Branches here 
00224       t_trackP                ->push_back( pTrack->p() );
00225       t_trackPt               ->push_back( pTrack->pt() );
00226       t_trackEta              ->push_back( pTrack->momentum().eta() );
00227       t_trackPhi              ->push_back( pTrack->momentum().phi() );
00228       t_trackHcalEta          ->push_back( closestCell.ieta() );
00229       t_trackHcalPhi          ->push_back( closestCell.iphi() );
00230       t_hCone                 ->push_back( hCone);
00231       t_conehmaxNearP         ->push_back( conehmaxNearP);
00232       t_eMipDR                ->push_back( eMipDR);
00233       t_eECALDR               ->push_back( eECALDR);
00234       t_eHCALDR               ->push_back( eHCALDR);
00235       t_e11x11_20Sig          ->push_back( e11x11_20SigP.first );
00236       t_e15x15_20Sig          ->push_back( e15x15_20SigP.first );
00237       t_eMipDR_1              ->push_back( eMipDR_1);
00238       t_eECALDR_1             ->push_back( eECALDR_1);
00239       t_eMipDR_2              ->push_back( eMipDR_2);
00240       t_eECALDR_2             ->push_back( eECALDR_2);
00241 
00242       if (myverbose > 0) {
00243         std::cout << "Track p " << pTrack->p() << " pt " << pTrack->pt()
00244                   << " eta " << pTrack->momentum().eta() << " phi "
00245                   << pTrack->momentum().phi() << " ieta/iphi ("
00246                   << closestCell.ieta() << ", " << closestCell.iphi() 
00247                   << ") Energy in cone " << hCone << " Charge Isolation "
00248                   << conehmaxNearP << " eMIP (" << eMipDR << ", "
00249                   << eMipDR_1 << ", " << eMipDR_2 << ")"
00250                   << " Neutral isolation (ECAL) (" << eECALDR-eMipDR << ", "
00251                   << eECALDR_1-eMipDR_1 << ", " << eECALDR_2-eMipDR_2 << ")"
00252                   << " (ECAL NxN) " << e15x15_20SigP.first-e11x11_20SigP.first
00253                   << " (HCAL) " << eHCALDR-hCone << std::endl;
00254       }
00255 
00256       if (doMC) {
00257         int nSimHits = -999;
00258         double hsim;
00259         std::map<std::string, double> hsimInfo;
00260         std::vector<int> multiplicity;
00261         hsim = spr::eCone_hcal(geo, pcalohh, trkDetItr->pointHCAL, 
00262                                trkDetItr->pointECAL, a_coneR, 
00263                                trkDetItr->directionHCAL, nSimHits);
00264         hsimInfo = spr::eHCALSimInfoCone(iEvent, pcalohh, SimTk, SimVtx, 
00265                                          pTrack, *associate, geo, 
00266                                          trkDetItr->pointHCAL, 
00267                                          trkDetItr->pointECAL, a_coneR,
00268                                          trkDetItr->directionHCAL,
00269                                          multiplicity);
00270 
00271         t_hsimInfoMatched   ->push_back(hsimInfo["eMatched"   ]);
00272         t_hsimInfoRest      ->push_back(hsimInfo["eRest"      ]);
00273         t_hsimInfoPhoton    ->push_back(hsimInfo["eGamma"     ]);
00274         t_hsimInfoNeutHad   ->push_back(hsimInfo["eNeutralHad"]);
00275         t_hsimInfoCharHad   ->push_back(hsimInfo["eChargedHad"]);
00276         t_hsimInfoPdgMatched->push_back(hsimInfo["pdgMatched" ]);
00277         t_hsimInfoTotal     ->push_back(hsimInfo["eTotal"     ]);
00278 
00279         t_hsimInfoNMatched  ->push_back(multiplicity.at(0));
00280         t_hsimInfoNTotal    ->push_back(multiplicity.at(1));
00281         t_hsimInfoNNeutHad  ->push_back(multiplicity.at(2));
00282         t_hsimInfoNCharHad  ->push_back(multiplicity.at(3));
00283         t_hsimInfoNPhoton   ->push_back(multiplicity.at(4));
00284         t_hsimInfoNRest     ->push_back(multiplicity.at(5));
00285 
00286         t_hsim              ->push_back(hsim                   );
00287         t_nSimHits          ->push_back(nSimHits               );
00288 
00289         if (myverbose > 0) {
00290           std::cout << "Matched (E) " << hsimInfo["eMatched"] << " (N) "
00291                     << multiplicity.at(0) << " Rest (E) " << hsimInfo["eRest"] 
00292                     << " (N) " << multiplicity.at(1) << " Gamma (E) "
00293                     << hsimInfo["eGamma"] << " (N) "  << multiplicity.at(2) 
00294                     << " Neutral Had (E) " << hsimInfo["eNeutralHad"]  
00295                     << " (N) "  << multiplicity.at(3) << " Charged Had (E) "
00296                     << hsimInfo["eChargedHad"] << " (N) " << multiplicity.at(4)
00297                     << " Total (E) " << hsimInfo["eTotal"] << " (N) "
00298                     << multiplicity.at(5) << " PDG " << hsimInfo["pdgMatched"] 
00299                     << " Total E " << hsim << " NHit " << nSimHits <<std::endl;
00300         }
00301       }
00302     }
00303   }
00304 
00305   //  delete associate;
00306   if (associate) delete associate;
00307   tree->Fill();
00308 }
00309 
00310 void IsolatedTracksHcalScale::beginJob() {
00311 
00312   nEventProc=0;
00313 
00314 
00316   tree = fs->make<TTree>("tree", "tree");
00317   tree->SetAutoSave(10000);
00318 
00319   tree->Branch("t_RunNo"              ,&t_RunNo               ,"t_RunNo/I");
00320   tree->Branch("t_Lumi"               ,&t_Lumi                ,"t_Lumi/I");
00321   tree->Branch("t_Bunch"              ,&t_Bunch               ,"t_Bunch/I");
00322 
00323   t_trackP              = new std::vector<double>();
00324   t_trackPt             = new std::vector<double>();
00325   t_trackEta            = new std::vector<double>();
00326   t_trackPhi            = new std::vector<double>();
00327   t_trackHcalEta        = new std::vector<double>();
00328   t_trackHcalPhi        = new std::vector<double>();
00329   t_hCone               = new std::vector<double>();
00330   t_conehmaxNearP       = new std::vector<double>();
00331   t_eMipDR              = new std::vector<double>();
00332   t_eECALDR             = new std::vector<double>();
00333   t_eHCALDR             = new std::vector<double>();
00334   t_e11x11_20Sig        = new std::vector<double>();
00335   t_e15x15_20Sig        = new std::vector<double>();
00336   t_eMipDR_1            = new std::vector<double>();
00337   t_eECALDR_1           = new std::vector<double>();
00338   t_eMipDR_2            = new std::vector<double>();
00339   t_eECALDR_2           = new std::vector<double>();
00340 
00341   tree->Branch("t_trackP",            "vector<double>", &t_trackP           );
00342   tree->Branch("t_trackPt",           "vector<double>", &t_trackPt          );
00343   tree->Branch("t_trackEta",          "vector<double>", &t_trackEta         );
00344   tree->Branch("t_trackPhi",          "vector<double>", &t_trackPhi         );
00345   tree->Branch("t_trackHcalEta",      "vector<double>", &t_trackHcalEta     );
00346   tree->Branch("t_trackHcalPhi",      "vector<double>", &t_trackHcalPhi     );
00347   tree->Branch("t_hCone",             "vector<double>", &t_hCone            );
00348   tree->Branch("t_conehmaxNearP",     "vector<double>", &t_conehmaxNearP    );
00349   tree->Branch("t_eMipDR",            "vector<double>", &t_eMipDR           );
00350   tree->Branch("t_eECALDR",           "vector<double>", &t_eECALDR          );
00351   tree->Branch("t_eHCALDR",           "vector<double>", &t_eHCALDR          );
00352   tree->Branch("t_e11x11_20Sig",      "vector<double>", &t_e11x11_20Sig     );
00353   tree->Branch("t_e15x15_20Sig",      "vector<double>", &t_e15x15_20Sig     );
00354   tree->Branch("t_eMipDR_1",          "vector<double>", &t_eMipDR_1         );
00355   tree->Branch("t_eECALDR_1",         "vector<double>", &t_eECALDR_1        );
00356   tree->Branch("t_eMipDR_2",          "vector<double>", &t_eMipDR_2         );
00357   tree->Branch("t_eECALDR_2",         "vector<double>", &t_eECALDR_2        );
00358 
00359   if (doMC) {
00360     t_hsimInfoMatched    = new std::vector<double>();
00361     t_hsimInfoRest       = new std::vector<double>();
00362     t_hsimInfoPhoton     = new std::vector<double>();
00363     t_hsimInfoNeutHad    = new std::vector<double>();
00364     t_hsimInfoCharHad    = new std::vector<double>();
00365     t_hsimInfoPdgMatched = new std::vector<double>();
00366     t_hsimInfoTotal      = new std::vector<double>();
00367     t_hsimInfoNMatched   = new std::vector<int>();
00368     t_hsimInfoNTotal     = new std::vector<int>();
00369     t_hsimInfoNNeutHad   = new std::vector<int>();
00370     t_hsimInfoNCharHad   = new std::vector<int>();
00371     t_hsimInfoNPhoton    = new std::vector<int>();
00372     t_hsimInfoNRest      = new std::vector<int>();
00373     t_hsim               = new std::vector<double>();
00374     t_nSimHits           = new std::vector<int>();
00375 
00376     tree->Branch("t_hsimInfoMatched",    "vector<double>", &t_hsimInfoMatched    );
00377     tree->Branch("t_hsimInfoRest",       "vector<double>", &t_hsimInfoRest    );
00378     tree->Branch("t_hsimInfoPhoton",     "vector<double>", &t_hsimInfoPhoton    );
00379     tree->Branch("t_hsimInfoNeutHad",    "vector<double>", &t_hsimInfoNeutHad    );
00380     tree->Branch("t_hsimInfoCharHad",    "vector<double>", &t_hsimInfoCharHad    );
00381     tree->Branch("t_hsimInfoPdgMatched", "vector<double>", &t_hsimInfoPdgMatched );
00382     tree->Branch("t_hsimInfoTotal",      "vector<double>", &t_hsimInfoTotal    );
00383     tree->Branch("t_hsimInfoNMatched",   "vector<int>",    &t_hsimInfoNMatched    );
00384     tree->Branch("t_hsimInfoNTotal",     "vector<int>",    &t_hsimInfoNTotal    );
00385     tree->Branch("t_hsimInfoNNeutHad",   "vector<int>",    &t_hsimInfoNNeutHad    );
00386     tree->Branch("t_hsimInfoNCharHad",   "vector<int>",    &t_hsimInfoNCharHad    );
00387     tree->Branch("t_hsimInfoNPhoton",    "vector<int>",    &t_hsimInfoNPhoton    );
00388     tree->Branch("t_hsimInfoNRest",      "vector<int>",    &t_hsimInfoNRest    );
00389     tree->Branch("t_hsim",               "vector<double>", &t_hsim    );
00390     tree->Branch("t_nSimHits",           "vector<int>",    &t_nSimHits    );
00391   }
00392 }
00393 
00394 void IsolatedTracksHcalScale::endJob() {
00395 
00396   std::cout << "Number of Events Processed " << nEventProc << std::endl;
00397 }
00398 
00399 void IsolatedTracksHcalScale::clearTreeVectors() {
00400 
00401   t_trackP            ->clear();
00402   t_trackPt           ->clear();
00403   t_trackEta          ->clear();
00404   t_trackPhi          ->clear();
00405   t_trackHcalEta      ->clear();
00406   t_trackHcalPhi      ->clear();
00407   t_hCone             ->clear();
00408   t_conehmaxNearP     ->clear();
00409   t_eMipDR            ->clear();
00410   t_eECALDR           ->clear();
00411   t_eHCALDR           ->clear();
00412   t_e11x11_20Sig      ->clear();
00413   t_e15x15_20Sig      ->clear();
00414   t_eMipDR_1          ->clear();
00415   t_eECALDR_1         ->clear();
00416   t_eMipDR_2          ->clear();
00417   t_eECALDR_2         ->clear();
00418 
00419   if (doMC) {
00420     t_hsimInfoMatched    ->clear();
00421     t_hsimInfoRest       ->clear();
00422     t_hsimInfoPhoton     ->clear();
00423     t_hsimInfoNeutHad    ->clear();
00424     t_hsimInfoCharHad    ->clear();
00425     t_hsimInfoPdgMatched ->clear();
00426     t_hsimInfoTotal      ->clear();
00427     t_hsimInfoNMatched   ->clear();
00428     t_hsimInfoNTotal     ->clear();
00429     t_hsimInfoNNeutHad   ->clear();
00430     t_hsimInfoNCharHad   ->clear();
00431     t_hsimInfoNPhoton    ->clear();
00432     t_hsimInfoNRest      ->clear();
00433     t_hsim               ->clear();
00434     t_nSimHits           ->clear();
00435   }
00436 }
00437 
00438 
00439 //define this as a plug-in
00440 DEFINE_FWK_MODULE(IsolatedTracksHcalScale);