Reference Manual

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    IsolatedTracksCone
00004 // Class:      IsolatedTracksCone
00005 // 
00006 
00015 //
00016 // Original Author: Jim Hirschauer (adaptation of Seema Sharma's
00017 // IsolatedTracksNew)
00018 //         Created:  Thu Nov  6 15:30:40 CST 2008
00019 //
00020 //
00021 #define _DEBUG_QUIET
00022 
00023 #include "Calibration/IsolatedParticles/plugins/IsolatedTracksCone.h"
00024 #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h"
00025 
00026 #include "DataFormats/HLTReco/interface/TriggerObject.h"
00027 #include "FWCore/Common/interface/TriggerNames.h"
00028 #include "FWCore/Framework/interface/LuminosityBlock.h"
00029 #include "DataFormats/Common/interface/TriggerResults.h"
00030 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
00031 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
00032 
00033 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
00034 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
00035 
00036 IsolatedTracksCone::IsolatedTracksCone(const edm::ParameterSet& iConfig) {
00037 
00038   //now do what ever initialization is needed
00039   doMC            = iConfig.getUntrackedParameter<bool>  ("DoMC", false); 
00040   myverbose_      = 
00041     iConfig.getUntrackedParameter<int>( "Verbosity", 5 );
00042   useJetTrigger_  = 
00043     iConfig.getUntrackedParameter<bool>( "useJetTrigger", false);
00044   drLeadJetVeto_  = 
00045     iConfig.getUntrackedParameter<double>( "drLeadJetVeto",  1.2 );
00046   ptMinLeadJet_   = 
00047     iConfig.getUntrackedParameter<double>( "ptMinLeadJet",  15.0 );
00048 
00049   debugTrks_          = 
00050     iConfig.getUntrackedParameter<int>("DebugTracks");
00051   printTrkHitPattern_ = 
00052     iConfig.getUntrackedParameter<bool>("PrintTrkHitPattern");
00053   
00054   minTrackP_     = 
00055     iConfig.getUntrackedParameter<double>( "minTrackP", 10.0);
00056   maxTrackEta_    = 
00057     iConfig.getUntrackedParameter<double>( "maxTrackEta", 5.0);
00058   maxNearTrackP_  = 
00059     iConfig.getUntrackedParameter<double>( "maxNearTrackP", 1.0);
00060 
00061   debugEcalSimInfo_   = 
00062     iConfig.getUntrackedParameter<int>("DebugEcalSimInfo");
00063 
00064   applyEcalIsolation_ = 
00065     iConfig.getUntrackedParameter<bool>("ApplyEcalIsolation");
00066 
00067   _L1extraTauJetSource = 
00068     iConfig.getParameter<edm::InputTag>("L1extraTauJetSource");
00069   _L1extraCenJetSource = 
00070     iConfig.getParameter<edm::InputTag>("L1extraCenJetSource");
00071   _L1extraFwdJetSource = 
00072     iConfig.getParameter<edm::InputTag>("L1extraFwdJetSource");
00073 
00074   edm::ParameterSet parameters = 
00075     iConfig.getParameter<edm::ParameterSet>("TrackAssociatorParameters");
00076   parameters_.loadParameters( parameters );
00077   trackAssociator_ =  new TrackDetectorAssociator();
00078   trackAssociator_->useDefaultPropagator();
00079 
00080   if(myverbose_>=0) {
00081     std::cout <<"Parameters read from config file \n" 
00082               << "myverbose_ "          << myverbose_     << "\n"     
00083               << "useJetTrigger_ "      << useJetTrigger_ << "\n"
00084               << "drLeadJetVeto_ "      << drLeadJetVeto_ << "\n"
00085               << "minTrackP_ "         << minTrackP_    << "\n"
00086               << "maxTrackEta_ "        << maxTrackEta_   << "\n"
00087               << "maxNearTrackP_ "      << maxNearTrackP_ 
00088               << std::endl;
00089   }
00090 }
00091 
00092 IsolatedTracksCone::~IsolatedTracksCone() {
00093   delete  trackAssociator_;
00094 }
00095 
00096 void IsolatedTracksCone::analyze(const edm::Event& iEvent, 
00097                                  const edm::EventSetup& iSetup) {
00098   
00099   unsigned int irun = (unsigned int)iEvent.id().run();
00100   unsigned int ilum = (unsigned int)iEvent.getLuminosityBlock().luminosityBlock();
00101   unsigned int ievt = (unsigned int)iEvent.id().event();
00102 
00103   
00104   
00105   clearTrackVectors();
00106   
00107   // check the L1 objects
00108   bool   L1Pass = false;
00109   leadL1JetPT=-999, leadL1JetEta=-999,  leadL1JetPhi=-999;
00110   if( !useJetTrigger_ ) {
00111     L1Pass = true;
00112   } else {
00113     edm::Handle<l1extra::L1JetParticleCollection> l1TauHandle;
00114     iEvent.getByLabel(_L1extraTauJetSource,l1TauHandle);
00115     l1extra::L1JetParticleCollection::const_iterator itr;
00116     for(itr = l1TauHandle->begin(); itr != l1TauHandle->end(); ++itr ) 
00117     {
00118       if( itr->pt()>leadL1JetPT ) {
00119         leadL1JetPT  = itr->pt();
00120         leadL1JetEta = itr->eta();
00121         leadL1JetPhi = itr->phi();
00122       }
00123     }
00124     edm::Handle<l1extra::L1JetParticleCollection> l1CenJetHandle;
00125     iEvent.getByLabel(_L1extraCenJetSource,l1CenJetHandle);
00126     for( itr = l1CenJetHandle->begin();  itr != l1CenJetHandle->end(); ++itr ) 
00127     {
00128       if( itr->pt()>leadL1JetPT ) {
00129         leadL1JetPT  = itr->pt();
00130         leadL1JetEta = itr->eta();
00131         leadL1JetPhi = itr->phi();
00132       }
00133     }
00134     edm::Handle<l1extra::L1JetParticleCollection> l1FwdJetHandle;
00135     iEvent.getByLabel(_L1extraFwdJetSource,l1FwdJetHandle);
00136     for( itr = l1FwdJetHandle->begin();  itr != l1FwdJetHandle->end(); ++itr ) 
00137     {
00138       if( itr->pt()>leadL1JetPT ) {
00139         leadL1JetPT  = itr->pt();
00140         leadL1JetEta = itr->eta();
00141         leadL1JetPhi = itr->phi();
00142       }
00143     }
00144     if(leadL1JetPT>ptMinLeadJet_) 
00145     {
00146       L1Pass = true;
00147     }
00148     
00149   }
00150   
00151 
00153   // Break now if L1Pass is false
00155   //   if (!L1Pass) {
00156   //     nEVT_failL1++;
00157   //     //    std::cout << "L1Pass is false : " << L1Pass << std::endl;
00158   //     return;  
00159   //   }
00160   
00162   // Get the collection handles
00164   
00165 
00166   edm::ESHandle<CaloGeometry> pG;
00167   iSetup.get<CaloGeometryRecord>().get(pG);
00168   const CaloGeometry* geo = pG.product();
00169   const CaloSubdetectorGeometry* gEB = 
00170     geo->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
00171   const CaloSubdetectorGeometry* gEE = 
00172     geo->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
00173   const CaloSubdetectorGeometry* gHB = 
00174     geo->getSubdetectorGeometry(DetId::Hcal,HcalBarrel);
00175   const CaloSubdetectorGeometry* gHE = 
00176     geo->getSubdetectorGeometry(DetId::Hcal,HcalEndcap);
00177   
00178   edm::ESHandle<CaloTopology> theCaloTopology;
00179   iSetup.get<CaloTopologyRecord>().get(theCaloTopology); 
00180   const CaloTopology *caloTopology = theCaloTopology.product();
00181   //  const CaloSubdetectorTopology* theEBTopology   = theCaloTopology->getSubdetectorTopology(DetId::Ecal,EcalBarrel);
00182   //  const CaloSubdetectorTopology* theEETopology   = theCaloTopology->getSubdetectorTopology(DetId::Ecal,EcalEndcap);
00183   
00184   edm::ESHandle<HcalTopology> htopo;
00185   iSetup.get<IdealGeometryRecord>().get(htopo);
00186   const HcalTopology* theHBHETopology = htopo.product();
00187   
00188   edm::Handle<EcalRecHitCollection> barrelRecHitsHandle;
00189   edm::Handle<EcalRecHitCollection> endcapRecHitsHandle;
00190   iEvent.getByLabel("ecalRecHit","EcalRecHitsEB",barrelRecHitsHandle);
00191   iEvent.getByLabel("ecalRecHit","EcalRecHitsEE",endcapRecHitsHandle);
00192   
00193   edm::Handle<HBHERecHitCollection> hbhe;
00194   iEvent.getByLabel("hbhereco",hbhe);
00195   const HBHERecHitCollection Hithbhe = *(hbhe.product());
00196 
00197   edm::Handle<reco::TrackCollection> trkCollection;
00198   iEvent.getByLabel("generalTracks", trkCollection);
00199   reco::TrackCollection::const_iterator trkItr;
00200   if(debugTrks_>1){
00201     std::cout << "Track Collection: " << std::endl;
00202     std::cout << "Number of Tracks " << trkCollection->size() << std::endl;
00203   }
00204   std::string theTrackQuality = "highPurity";
00205   reco::TrackBase::TrackQuality trackQuality_=
00206     reco::TrackBase::qualityByName(theTrackQuality);
00207   
00208   //get Handles to SimTracks and SimHits
00209   edm::Handle<edm::SimTrackContainer> SimTk;
00210   if (doMC) iEvent.getByLabel("g4SimHits",SimTk);
00211   edm::SimTrackContainer::const_iterator simTrkItr;
00212 
00213   edm::Handle<edm::SimVertexContainer> SimVtx;
00214   if (doMC) iEvent.getByLabel("g4SimHits",SimVtx);
00215   edm::SimVertexContainer::const_iterator vtxItr = SimVtx->begin();
00216 
00217   //get Handles to PCaloHitContainers of eb/ee/hbhe
00218   edm::Handle<edm::PCaloHitContainer> pcaloeb;
00219   if (doMC) iEvent.getByLabel("g4SimHits", "EcalHitsEB", pcaloeb);
00220 
00221   edm::Handle<edm::PCaloHitContainer> pcaloee;
00222   if (doMC) iEvent.getByLabel("g4SimHits", "EcalHitsEE", pcaloee);
00223 
00224   edm::Handle<edm::PCaloHitContainer> pcalohh;
00225   if (doMC) iEvent.getByLabel("g4SimHits", "HcalHits", pcalohh);
00226   
00227   
00228   
00230   // Get HLT_IsoTrackHB/HE Information
00232     
00233   edm::InputTag theTriggerResultsLabel ("TriggerResults","","HLT");
00234   edm::Handle<edm::TriggerResults> triggerResults;
00235   iEvent.getByLabel( theTriggerResultsLabel, triggerResults);
00236 
00237   
00238 
00239   std::vector<int> v_hlTriggers;
00240   int hltHB(-99);
00241   int hltHE(-99);
00242   int hltL1Jet15                        (-99);
00243   int hltJet30                  (-99);
00244   int hltJet50                  (-99);
00245   int hltJet80                  (-99);
00246   int hltJet110                 (-99);
00247   int hltJet140                 (-99);
00248   int hltJet180                 (-99);
00249   int hltL1SingleEG5            (-99);
00250   int hltZeroBias               (-99);
00251   int hltMinBiasHcal            (-99);
00252   int hltMinBiasEcal            (-99);
00253   int hltMinBiasPixel           (-99);
00254   int hltSingleIsoTau30_Trk5    (-99);
00255   int hltDoubleLooseIsoTau15_Trk5(-99);
00256 
00257   if (triggerResults.isValid()) {
00258     
00259     const edm::TriggerNames & triggerNames = iEvent.triggerNames(*triggerResults);
00260     // TriggerNames class  triggerNames.init(*triggerResults);
00261     
00262 
00263     for (unsigned int i=0; i<triggerResults->size(); i++){
00264       //      std::cout << "triggerNames.triggerName(" << i << ") = " << triggerNames.triggerName(i) << std::endl;
00265       if (triggerNames.triggerName(i) == "HLT_IsoTrackHE_1E31") hltHE = triggerResults->accept(i);
00266       if (triggerNames.triggerName(i) == "HLT_IsoTrackHB_1E31") hltHB = triggerResults->accept(i);
00267       if (triggerNames.triggerName(i) == "HLT_L1Jet15"                  ) hltL1Jet15                     = triggerResults->accept(i);
00268       if (triggerNames.triggerName(i) == "HLT_Jet30"                    ) hltJet30                       = triggerResults->accept(i);
00269       if (triggerNames.triggerName(i) == "HLT_Jet50"                    ) hltJet50                       = triggerResults->accept(i);
00270       if (triggerNames.triggerName(i) == "HLT_Jet80"                    ) hltJet80                       = triggerResults->accept(i);
00271       if (triggerNames.triggerName(i) == "HLT_Jet110"                   ) hltJet110                      = triggerResults->accept(i);
00272       if (triggerNames.triggerName(i) == "HLT_Jet140"                   ) hltJet140                      = triggerResults->accept(i);
00273       if (triggerNames.triggerName(i) == "HLT_Jet180"                   ) hltJet180                      = triggerResults->accept(i);
00274       if (triggerNames.triggerName(i) == "HLT_L1SingleEG5"              ) hltL1SingleEG5                 = triggerResults->accept(i);
00275       if (triggerNames.triggerName(i) == "HLT_ZeroBias"                 ) hltZeroBias                    = triggerResults->accept(i);
00276       if (triggerNames.triggerName(i) == "HLT_MinBiasHcal"              ) hltMinBiasHcal                 = triggerResults->accept(i);
00277       if (triggerNames.triggerName(i) == "HLT_MinBiasEcal"              ) hltMinBiasEcal                 = triggerResults->accept(i);
00278       if (triggerNames.triggerName(i) == "HLT_MinBiasPixel"             ) hltMinBiasPixel                = triggerResults->accept(i);
00279       if (triggerNames.triggerName(i) == "HLT_SingleIsoTau30_Trk5"      ) hltSingleIsoTau30_Trk5         = triggerResults->accept(i);
00280       if (triggerNames.triggerName(i) == "HLT_DoubleLooseIsoTau15_Trk5" ) hltDoubleLooseIsoTau15_Trk5    = triggerResults->accept(i);
00281     }
00282   }
00283   
00284     
00285 
00286   
00288   // Primary loop over tracks
00290   TrackerHitAssociator* associate=0;
00291   if (doMC) associate = new TrackerHitAssociator(iEvent);
00292   
00293 
00294   nTRK      = 0;
00295   nRawTRK   = 0;
00296   nFailPt   = 0;
00297   nFailEta  = 0;
00298   nFailHighPurityQaul = 0;
00299   nMissEcal = 0;
00300   nMissHcal = 0;
00301 
00302   for( trkItr = trkCollection->begin(); 
00303        trkItr != trkCollection->end(); ++trkItr)
00304   {
00305 
00306     nRawTRK++;
00307 
00308     const reco::Track* pTrack = &(*trkItr);
00309 
00311     // Check for min Pt and max Eta P
00313 
00314     bool trkQual  = pTrack->quality(trackQuality_);
00315     bool goodPt   = pTrack->p()>minTrackP_;
00316     bool goodEta  = std::abs(pTrack->momentum().eta())<maxTrackEta_;
00317 
00318     double eta1       = pTrack->momentum().eta();
00319     double phi1       = pTrack->momentum().phi();
00320     double pt1        = pTrack->pt();
00321     double p1         = pTrack->p();
00322 
00323 
00324     if (!goodEta){
00325       nFailEta++;
00326     }
00327     if (!goodPt){
00328       nFailPt++;
00329     }
00330     if (!trkQual){
00331       nFailHighPurityQaul++;
00332     }
00333     
00334     hRawPt ->Fill(pt1 );
00335     hRawP  ->Fill(p1  );
00336     hRawEta->Fill(eta1);
00337     hRawPhi->Fill(phi1);
00338       
00339     if( !goodEta || !goodPt || !trkQual ) continue; // Skip to next track
00340     
00342     // Find track trajectory
00344     
00345       
00346     const FreeTrajectoryState fts1 = 
00347       trackAssociator_->getFreeTrajectoryState(iSetup, *pTrack);
00348     
00349       
00350     TrackDetMatchInfo info1 = 
00351       trackAssociator_->associate(iEvent, iSetup, fts1, parameters_);
00352     
00353       
00354 
00356     // First confirm track makes it to Hcal
00358 
00359     if (info1.trkGlobPosAtHcal.x()==0 && 
00360         info1.trkGlobPosAtHcal.y()==0 && 
00361         info1.trkGlobPosAtHcal.z()==0) 
00362     {
00363       nMissHcal++;
00364       continue;      
00365     }
00366     
00367     const GlobalPoint hpoint1(info1.trkGlobPosAtHcal.x(),
00368                               info1.trkGlobPosAtHcal.y(),
00369                               info1.trkGlobPosAtHcal.z());
00370 
00371       
00372 
00374     // Get basic quantities
00376 
00377     const reco::HitPattern& hitp = pTrack->hitPattern();
00378     int nLayersCrossed = hitp.trackerLayersWithMeasurement();        
00379     int nOuterHits     = hitp.stripTOBLayersWithMeasurement()
00380       +hitp.stripTECLayersWithMeasurement() ;
00381 
00382     
00383     double simP = 0;
00384     if (doMC) {
00385       edm::SimTrackContainer::const_iterator matchedSimTrk = 
00386         spr::matchedSimTrack(iEvent, SimTk, SimVtx, pTrack, *associate, false);
00387       simP = matchedSimTrk->momentum().P();
00388     }
00390     // Get Ecal Point
00392 
00393     const GlobalPoint point1(info1.trkGlobPosAtEcal.x(),
00394                              info1.trkGlobPosAtEcal.y(),
00395                              info1.trkGlobPosAtEcal.z());
00396     
00397 
00398     // Sanity check that track hits Ecal
00399 
00400     if (info1.trkGlobPosAtEcal.x()==0 && 
00401         info1.trkGlobPosAtEcal.y()==0 && 
00402         info1.trkGlobPosAtEcal.z()==0) 
00403     {
00404       std::cout << "Track doesn't reach Ecal." << std::endl;
00405       nMissEcal++;
00406       continue;
00407     }
00408 
00409     // Get Track Momentum - make sure you have latest version of TrackDetMatchInfo
00410     
00411     GlobalVector trackMomAtEcal = info1.trkMomAtEcal;
00412     GlobalVector trackMomAtHcal = info1.trkMomAtHcal;
00413 
00415     // If using Jet trigger, get distance from leading jet
00417 
00418     double drFromLeadJet = 999.0;
00419     if( useJetTrigger_ ) {
00420       double dphi = DeltaPhi(phi1, leadL1JetPhi);
00421       double deta = eta1 - leadL1JetEta;
00422       drFromLeadJet =  sqrt(dphi*dphi + deta*deta);
00423     }
00424     
00425 
00427     // Define Arrays for sizes of Charge, Neut Iso Radii and
00428     // Clustering Cone Radius.
00430 
00431     const int a_size = 7;
00432     double a_coneR[a_size];
00433     double a_charIsoR[a_size];
00434     double a_neutIsoR[a_size];
00435 
00436     a_coneR[0] = 17.49; // = area of 2x2
00437     a_coneR[1] = 26.23; // = area of 3x3
00438     a_coneR[2] = 30.61;
00439     a_coneR[3] = 34.98; // = area of 4x4
00440     a_coneR[4] = 39.35;
00441     a_coneR[5] = 43.72; // = area of 5x5
00442     a_coneR[6] = 52.46; // = area of 6x6
00443     
00444     for (int i=0; i<a_size; i++){
00445       a_charIsoR[i] = a_coneR[i]+28.9; // 28.9 gives 55.1 for 3x3 benchmark 
00446       a_neutIsoR[i] = a_charIsoR[i]*0.726; // Ecal radius = 0.726*Hcal radius
00447     }
00448     
00450     // Do Neutral Iso in radius on Ecal surface.
00452   
00453     // NxN cluster
00454     double e3x3=-999.0;
00455     double trkEcalEne =-999.0;
00456     edm::ESHandle<EcalSeverityLevelAlgo> sevlv;
00457     iSetup.get<EcalSeverityLevelAlgoRcd>().get(sevlv);
00458 
00459     if(std::abs(point1.eta())<1.479) {
00460       const DetId isoCell = gEB->getClosestCell(point1);
00461       e3x3   = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, geo, caloTopology, sevlv.product(),1,1);  
00462       trkEcalEne   = spr::eCaloSimInfo(iEvent, geo, pcaloeb, pcaloee, SimTk, SimVtx, pTrack, *associate);
00463     } else {
00464       const DetId isoCell = gEE->getClosestCell(point1);
00465       e3x3   = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, geo, caloTopology, sevlv.product(),1,1);
00466       trkEcalEne   = spr::eCaloSimInfo(iEvent, geo, pcaloeb, pcaloee, SimTk, SimVtx, pTrack, *associate);
00467     }
00468 
00469     // Cone cluster
00470 
00471     // Set up array of cone sizes for MIP cut
00472     const int a_mip_size = 5;
00473     double a_mipR[a_mip_size];
00474     a_mipR[0] = 3.84; // = area of 3x3 ecal
00475     a_mipR[1] = 14.0;
00476     a_mipR[2] = 19.0;
00477     a_mipR[3] = 24.0;
00478     a_mipR[4] = 9.0;  // = from standard analyzer
00479 
00480     std::vector<double> v_eDR;
00481     for (int i = 0 ; i < a_size ; i++){
00482       int nRH_eDR = 0;
00483 
00484       // Cone in ecal
00485       double eDR = spr::eCone_ecal(geo, 
00486                                    barrelRecHitsHandle, 
00487                                    endcapRecHitsHandle, 
00488                                    hpoint1, point1, 
00489                                    a_neutIsoR[i],  
00490                                    trackMomAtEcal, nRH_eDR);
00491       v_eDR.push_back(eDR);
00492       
00493     }
00494 
00495     std::vector<double> v_eMipDR;
00496     for (int i = 0 ; i < a_mip_size ; i++){
00497       int nRH_eMipDR = 0;
00498       double eMipDR = spr::eCone_ecal(geo, barrelRecHitsHandle, 
00499                                       endcapRecHitsHandle, 
00500                                       hpoint1, point1, 
00501                                       a_mipR[i], trackMomAtEcal, nRH_eMipDR);
00502       
00503       v_eMipDR.push_back(eMipDR);
00504     }
00505     
00506             
00508     // Do charge isolation in radius at Hcal surface for 5 different
00509     // radii defined above in a_charIso
00511     
00512     std::vector<double> v_hmaxNearP_goodTrk;
00513     std::vector<double> v_hmaxNearP        ;
00514     std::vector<int>    v_hnNearTRKs       ;
00515     std::vector<int>    v_hnLayers_maxNearP;
00516     std::vector<int>    v_htrkQual_maxNearP;
00517 
00518     std::vector<double> v_cone_hmaxNearP_goodTrk;
00519     std::vector<double> v_cone_hmaxNearP        ;
00520     std::vector<int>    v_cone_hnNearTRKs       ;
00521     std::vector<int>    v_cone_hnLayers_maxNearP;
00522     std::vector<int>    v_cone_htrkQual_maxNearP;
00523 
00524     for (int i = 0 ; i < a_size ; i++){
00525 
00526       double hmaxNearP         = -999.0;
00527       int    hnNearTRKs        =  0;
00528       int    hnLayers_maxNearP =  0;
00529       int    htrkQual_maxNearP = -1; 
00530       double hmaxNearP_goodTrk = -999.0;
00531 
00532       double conehmaxNearP         = -999.0;
00533       int    conehnNearTRKs        =  0;
00534       int    conehnLayers_maxNearP =  0;
00535       int    conehtrkQual_maxNearP = -1; 
00536       double conehmaxNearP_goodTrk = -999.0;
00537 
00538       conehmaxNearP = spr::coneChargeIsolation(iEvent, iSetup, 
00539                                                trkItr, trkCollection, 
00540                                                *trackAssociator_, parameters_, 
00541                                                theTrackQuality, conehnNearTRKs,
00542                                                conehnLayers_maxNearP,
00543                                                conehtrkQual_maxNearP, 
00544                                                conehmaxNearP_goodTrk, 
00545                                                hpoint1, trackMomAtHcal, 
00546                                                a_charIsoR[i]); 
00547 
00548       v_hmaxNearP_goodTrk.push_back(hmaxNearP_goodTrk);
00549       v_hmaxNearP        .push_back(hmaxNearP        );
00550       v_hnNearTRKs       .push_back(hnNearTRKs       );
00551       v_hnLayers_maxNearP.push_back(hnLayers_maxNearP);
00552       v_htrkQual_maxNearP.push_back(htrkQual_maxNearP);
00553 
00554       v_cone_hmaxNearP_goodTrk.push_back(conehmaxNearP_goodTrk);
00555       v_cone_hmaxNearP        .push_back(conehmaxNearP        );
00556       v_cone_hnNearTRKs       .push_back(conehnNearTRKs       );
00557       v_cone_hnLayers_maxNearP.push_back(conehnLayers_maxNearP);
00558       v_cone_htrkQual_maxNearP.push_back(conehtrkQual_maxNearP);
00559       
00560     }
00561     
00562     double h3x3=-999.0, h5x5=-999.0;
00563     double hsim3x3=-999.0, hsim5x5=-999.0, trkHcalEne=-999.0;
00564     std::map<std::string, double> hsimInfo3x3, hsimInfo5x5;
00565     double distFromHotCell_h3x3 = -99.;
00566     int ietaFromHotCell_h3x3    = -99;
00567     int iphiFromHotCell_h3x3    = -99;
00568     double distFromHotCell_h5x5 = -99.;
00569     int ietaFromHotCell_h5x5    = -99;
00570     int iphiFromHotCell_h5x5    = -99;
00571 
00572     GlobalPoint gPosHotCell_h3x3(0.,0.,0.);
00573     GlobalPoint gPosHotCell_h5x5(0.,0.,0.);
00574 
00575     int nRH_h3x3(0), nRH_h5x5(0);
00576 
00577     // Hcal Energy Clustering
00578     
00579     // Get closetcell for ietaFromHotCell and iphiFromHotCell
00580     DetId ClosestCell;
00581     if( std::abs(pTrack->eta())<1.392 ) {
00582       ClosestCell = gHB->getClosestCell(hpoint1);
00583     } else {
00584       ClosestCell = gHE->getClosestCell(hpoint1);
00585     }
00586     // Transform into HcalDetId so that I can get ieta, iphi later.
00587     HcalDetId ClosestCell_HcalDetId(ClosestCell.rawId());
00588 
00589     // Using NxN Cluster
00590     std::vector<int>    v_RH_h3x3_ieta;
00591     std::vector<int>    v_RH_h3x3_iphi;
00592     std::vector<double> v_RH_h3x3_ene;
00593     std::vector<int>    v_RH_h5x5_ieta;
00594     std::vector<int>    v_RH_h5x5_iphi;
00595     std::vector<double> v_RH_h5x5_ene;
00596     
00597     
00598     h3x3 = spr::eHCALmatrix(geo, theHBHETopology, ClosestCell, hbhe,1,1,
00599                             nRH_h3x3, v_RH_h3x3_ieta, v_RH_h3x3_iphi, v_RH_h3x3_ene, 
00600                             gPosHotCell_h3x3);
00601     distFromHotCell_h3x3 = spr::getDistInPlaneTrackDir(hpoint1, trackMomAtHcal, gPosHotCell_h3x3);
00602     
00603     h5x5 = spr::eHCALmatrix(geo, theHBHETopology, ClosestCell, hbhe,2,2,
00604                             nRH_h5x5, v_RH_h5x5_ieta, v_RH_h5x5_iphi,  v_RH_h5x5_ene, 
00605                             gPosHotCell_h5x5);
00606     distFromHotCell_h5x5 = spr::getDistInPlaneTrackDir(hpoint1, trackMomAtHcal, gPosHotCell_h5x5);
00607     
00608 
00609     //     double heta = (double)hpoint1.eta();
00610     //     double hphi = (double)hpoint1.phi();
00611     std::vector<int> multiplicity3x3;
00612     std::vector<int> multiplicity5x5;
00613     if (doMC) {
00614       hsim3x3   = spr::eHCALmatrix(theHBHETopology, ClosestCell, 
00615                                    pcalohh,1,1);
00616       hsim5x5   = spr::eHCALmatrix(theHBHETopology, ClosestCell,
00617                                    pcalohh,2,2);
00618 
00619       hsimInfo3x3 = spr::eHCALSimInfo(iEvent, theHBHETopology, ClosestCell, pcalohh, SimTk, SimVtx, pTrack, *associate, 1,1, multiplicity3x3);
00620       hsimInfo5x5 = spr::eHCALSimInfo(iEvent, theHBHETopology, ClosestCell, pcalohh, SimTk, SimVtx, pTrack, *associate, 2,2, multiplicity5x5);
00621     
00622       // Get energy from all simhits in hcal associated with iso track
00623       trkHcalEne   = spr::eCaloSimInfo(iEvent, geo, pcalohh, SimTk, SimVtx, pTrack, *associate);
00624     }
00625 
00626     // Finally for cones of varying radii.
00627     std::vector<double> v_hsimInfoConeMatched;
00628     std::vector<double> v_hsimInfoConeRest   ;
00629     std::vector<double> v_hsimInfoConePhoton  ;
00630     std::vector<double> v_hsimInfoConeNeutHad;
00631     std::vector<double> v_hsimInfoConeCharHad;
00632     std::vector<double> v_hsimInfoConePdgMatched;
00633     std::vector<double> v_hsimInfoConeTotal  ;
00634 
00635     std::vector<int> v_hsimInfoConeNMatched;
00636     std::vector<int> v_hsimInfoConeNTotal  ;
00637     std::vector<int> v_hsimInfoConeNNeutHad;
00638     std::vector<int> v_hsimInfoConeNCharHad;
00639     std::vector<int> v_hsimInfoConeNPhoton ;
00640     std::vector<int> v_hsimInfoConeNRest   ;
00641 
00642     std::vector<double> v_hsimCone           ;
00643     std::vector<double> v_hCone              ;
00644     
00645     std::vector<int>    v_nRecHitsCone       ;
00646     std::vector<int>    v_nSimHitsCone       ;
00647 
00648     std::vector<double> v_distFromHotCell;
00649     std::vector<int> v_ietaFromHotCell;
00650     std::vector<int> v_iphiFromHotCell;
00651     GlobalPoint gposHotCell(0.,0.,0.);
00652     
00653     
00654     std::vector<int>    v_RH_r26_ieta;
00655     std::vector<int>    v_RH_r26_iphi;
00656     std::vector<double> v_RH_r26_ene;
00657     std::vector<int>    v_RH_r44_ieta;
00658     std::vector<int>    v_RH_r44_iphi;
00659     std::vector<double> v_RH_r44_ene;
00660     
00661     
00662         
00663     for (int i = 0 ; i < a_size ; i++){
00664 
00665       
00666       std::map<std::string, double> hsimInfoCone;
00667       double hsimCone = -999.0, hCone = -999.0;
00668       double distFromHotCell = -99.0;
00669       int ietaFromHotCell = -99;
00670       int iphiFromHotCell = -99;
00671       int ietaHotCell = -99;
00672       int iphiHotCell = -99;
00673       int nRecHitsCone = -999;
00674       int nSimHitsCone = -999;
00675 
00676       std::vector<int> multiplicityCone;
00677       std::vector<DetId> coneRecHitDetIds;
00678       if (doMC) 
00679         hsimCone = spr::eCone_hcal(geo, pcalohh, hpoint1, point1, 
00680                                    a_coneR[i], trackMomAtHcal, nSimHitsCone);
00681       
00682       // If needed, get ieta and iphi of rechits for cones of 23.25
00683       // and for hitmap for debugging
00684       bool makeHitmaps = false;
00685       if (a_coneR[i] == 26.23 && makeHitmaps)
00686       {
00687         
00688         hCone = spr::eCone_hcal(geo, hbhe, hpoint1, point1, 
00689                                 a_coneR[i], trackMomAtHcal,nRecHitsCone,
00690                                 v_RH_r26_ieta, v_RH_r26_iphi,  v_RH_r26_ene, 
00691                                 coneRecHitDetIds, distFromHotCell, 
00692                                 ietaHotCell, iphiHotCell, gposHotCell);
00693       } 
00694       else if (a_coneR[i] == 43.72 && makeHitmaps)
00695       {
00696         
00697         hCone = spr::eCone_hcal(geo, hbhe, hpoint1, point1, 
00698                                 a_coneR[i], trackMomAtHcal,nRecHitsCone,
00699                                 v_RH_r44_ieta, v_RH_r44_iphi,  v_RH_r44_ene, 
00700                                 coneRecHitDetIds, distFromHotCell, 
00701                                 ietaHotCell, iphiHotCell, gposHotCell);
00702       } 
00703       else 
00704       {
00705         
00706         hCone = spr::eCone_hcal(geo, hbhe, hpoint1, point1, 
00707                                 a_coneR[i], trackMomAtHcal, nRecHitsCone, 
00708                                 coneRecHitDetIds, distFromHotCell, 
00709                                 ietaHotCell, iphiHotCell, gposHotCell);
00710       }
00711 
00712       
00713       
00714       if (ietaHotCell != 99){
00715         ietaFromHotCell = ietaHotCell-ClosestCell_HcalDetId.ieta();
00716         iphiFromHotCell = iphiHotCell-ClosestCell_HcalDetId.iphi();
00717       }
00718       
00719       // SimHits NOT matched to RecHits
00720       if (doMC) {
00721         hsimInfoCone = spr::eHCALSimInfoCone(iEvent,pcalohh, SimTk, SimVtx, pTrack, *associate, geo, hpoint1, point1, a_coneR[i], trackMomAtHcal, multiplicityCone);      
00722       
00723       
00724           
00725         // SimHits matched to RecHits
00726         //       hsimInfoCone = spr::eHCALSimInfoCone(iEvent,pcalohh, SimTk, SimVtx, 
00727         //                                                 pTrack, *associate, 
00728         //                                                 geo, hpoint1, point1, 
00729         //                                                 a_coneR[i], trackMomAtHcal, 
00730         //                                                 multiplicityCone,
00731         //                                                 coneRecHitDetIds);      
00732       
00733         v_hsimInfoConeMatched   .push_back(hsimInfoCone["eMatched"   ]);
00734         v_hsimInfoConeRest      .push_back(hsimInfoCone["eRest"      ]);
00735         v_hsimInfoConePhoton    .push_back(hsimInfoCone["eGamma"     ]);
00736         v_hsimInfoConeNeutHad   .push_back(hsimInfoCone["eNeutralHad"]);
00737         v_hsimInfoConeCharHad   .push_back(hsimInfoCone["eChargedHad"]);
00738         v_hsimInfoConePdgMatched.push_back(hsimInfoCone["pdgMatched" ]);
00739         v_hsimInfoConeTotal     .push_back(hsimInfoCone["eTotal"     ]);
00740 
00741         v_hsimInfoConeNMatched  .push_back(multiplicityCone.at(0));
00742 
00743         v_hsimInfoConeNTotal      .push_back(multiplicityCone.at(1));
00744         v_hsimInfoConeNNeutHad    .push_back(multiplicityCone.at(2));
00745         v_hsimInfoConeNCharHad    .push_back(multiplicityCone.at(3));
00746         v_hsimInfoConeNPhoton     .push_back(multiplicityCone.at(4));
00747         v_hsimInfoConeNRest       .push_back(multiplicityCone.at(5));
00748 
00749         v_hsimCone                .push_back(hsimCone                   );
00750         v_nSimHitsCone            .push_back(nSimHitsCone               );
00751       }
00752       v_hCone                   .push_back(hCone                      );
00753       v_nRecHitsCone            .push_back(nRecHitsCone               );
00754       
00755       v_distFromHotCell         .push_back(distFromHotCell            );
00756       v_ietaFromHotCell         .push_back(ietaFromHotCell            );
00757       v_iphiFromHotCell         .push_back(iphiFromHotCell            );
00758       
00759           
00760     }
00761     
00762  
00764     // Fill Vectors that go into root file
00766 
00767     t_v_hnNearTRKs        ->push_back(v_hnNearTRKs                 );
00768     t_v_hnLayers_maxNearP ->push_back(v_hnLayers_maxNearP          );
00769     t_v_htrkQual_maxNearP ->push_back(v_htrkQual_maxNearP          );
00770     t_v_hmaxNearP_goodTrk ->push_back(v_hmaxNearP_goodTrk          ); 
00771     t_v_hmaxNearP         ->push_back(v_hmaxNearP                  );    
00772 
00773     t_v_cone_hnNearTRKs        ->push_back(v_cone_hnNearTRKs       );
00774     t_v_cone_hnLayers_maxNearP ->push_back(v_cone_hnLayers_maxNearP);
00775     t_v_cone_htrkQual_maxNearP ->push_back(v_cone_htrkQual_maxNearP);
00776     t_v_cone_hmaxNearP_goodTrk ->push_back(v_cone_hmaxNearP_goodTrk); 
00777     t_v_cone_hmaxNearP         ->push_back(v_cone_hmaxNearP        );    
00778 
00779     //    t_hScale            ->push_back(hScale                     );
00780     t_trkNOuterHits     ->push_back(nOuterHits                 );
00781     t_trkNLayersCrossed ->push_back(nLayersCrossed             );
00782     t_dtFromLeadJet     ->push_back(drFromLeadJet              );
00783     t_trkP              ->push_back(p1                         );
00784     t_trkPt             ->push_back(pt1                        );
00785     t_trkEta            ->push_back(eta1                       );
00786     t_trkPhi            ->push_back(phi1                       );
00787 
00788     t_e3x3              ->push_back(e3x3                       );
00789     t_v_eDR             ->push_back(v_eDR                      );
00790     t_v_eMipDR          ->push_back(v_eMipDR                   );
00791 
00792     t_h3x3              ->push_back(h3x3                       );
00793     t_h5x5              ->push_back(h5x5                       );
00794     t_nRH_h3x3          ->push_back(nRH_h3x3                   );
00795     t_nRH_h5x5          ->push_back(nRH_h5x5                   );
00796     
00797     t_v_RH_h3x3_ieta    ->push_back(v_RH_h3x3_ieta);
00798     t_v_RH_h3x3_iphi    ->push_back(v_RH_h3x3_iphi);
00799     t_v_RH_h3x3_ene     ->push_back(v_RH_h3x3_ene);
00800     t_v_RH_h5x5_ieta    ->push_back(v_RH_h5x5_ieta);
00801     t_v_RH_h5x5_iphi    ->push_back(v_RH_h5x5_iphi);
00802     t_v_RH_h5x5_ene     ->push_back(v_RH_h5x5_ene);
00803 
00804     if (doMC) {
00805       t_simP              ->push_back(simP                       );
00806       t_hsim3x3           ->push_back(hsim3x3                    );
00807       t_hsim5x5           ->push_back(hsim5x5                    );
00808 
00809       t_hsim3x3Matched    ->push_back(hsimInfo3x3["eMatched"]    );
00810       t_hsim5x5Matched    ->push_back(hsimInfo5x5["eMatched"]    );
00811       t_hsim3x3Rest       ->push_back(hsimInfo3x3["eRest"]       );
00812       t_hsim5x5Rest       ->push_back(hsimInfo5x5["eRest"]       );
00813       t_hsim3x3Photon     ->push_back(hsimInfo3x3["eGamma"]      );
00814       t_hsim5x5Photon     ->push_back(hsimInfo5x5["eGamma"]      );
00815       t_hsim3x3NeutHad    ->push_back(hsimInfo3x3["eNeutralHad"] );
00816       t_hsim5x5NeutHad    ->push_back(hsimInfo5x5["eNeutralHad"] );
00817       t_hsim3x3CharHad    ->push_back(hsimInfo3x3["eChargedHad"] );
00818       t_hsim5x5CharHad    ->push_back(hsimInfo5x5["eChargedHad"] );
00819       t_hsim3x3Total      ->push_back(hsimInfo3x3["eTotal"]      );
00820       t_hsim5x5Total      ->push_back(hsimInfo5x5["eTotal"]      );
00821       t_hsim3x3PdgMatched ->push_back(hsimInfo3x3["pdgMatched"]  );
00822       t_hsim5x5PdgMatched ->push_back(hsimInfo5x5["pdgMatched"]  );
00823 
00824       t_hsim3x3NMatched   ->push_back(multiplicity3x3.at(0));
00825       t_hsim3x3NTotal     ->push_back(multiplicity3x3.at(1));
00826       t_hsim3x3NNeutHad   ->push_back(multiplicity3x3.at(2));
00827       t_hsim3x3NCharHad   ->push_back(multiplicity3x3.at(3));
00828       t_hsim3x3NPhoton    ->push_back(multiplicity3x3.at(4));
00829       t_hsim3x3NRest      ->push_back(multiplicity3x3.at(5));
00830 
00831       t_hsim5x5NMatched   ->push_back(multiplicity5x5.at(0));
00832       t_hsim5x5NTotal     ->push_back(multiplicity5x5.at(1));
00833       t_hsim5x5NNeutHad   ->push_back(multiplicity5x5.at(2));
00834       t_hsim5x5NCharHad   ->push_back(multiplicity5x5.at(3));
00835       t_hsim5x5NPhoton    ->push_back(multiplicity5x5.at(4));
00836       t_hsim5x5NRest      ->push_back(multiplicity5x5.at(5));
00837     }
00838     
00839     t_distFromHotCell_h3x3->push_back(distFromHotCell_h3x3);
00840     t_ietaFromHotCell_h3x3->push_back(ietaFromHotCell_h3x3);
00841     t_iphiFromHotCell_h3x3->push_back(iphiFromHotCell_h3x3);
00842     t_distFromHotCell_h5x5->push_back(distFromHotCell_h5x5);
00843     t_ietaFromHotCell_h5x5->push_back(ietaFromHotCell_h5x5);
00844     t_iphiFromHotCell_h5x5->push_back(iphiFromHotCell_h5x5);
00845     
00846     if (doMC) {
00847       t_trkHcalEne                ->push_back(trkHcalEne                 );
00848       t_trkEcalEne                ->push_back(trkEcalEne                 );
00849 
00850       t_v_hsimInfoConeMatched     ->push_back(v_hsimInfoConeMatched   );
00851       t_v_hsimInfoConeRest        ->push_back(v_hsimInfoConeRest      );
00852       t_v_hsimInfoConePhoton      ->push_back(v_hsimInfoConePhoton   );
00853       t_v_hsimInfoConeNeutHad     ->push_back(v_hsimInfoConeNeutHad   );
00854       t_v_hsimInfoConeCharHad     ->push_back(v_hsimInfoConeCharHad   );
00855       t_v_hsimInfoConePdgMatched  ->push_back(v_hsimInfoConePdgMatched);  
00856       t_v_hsimInfoConeTotal       ->push_back(v_hsimInfoConeTotal     );
00857 
00858       t_v_hsimInfoConeNMatched    ->push_back(v_hsimInfoConeNMatched  );
00859       t_v_hsimInfoConeNTotal      ->push_back(v_hsimInfoConeNTotal    );
00860       t_v_hsimInfoConeNNeutHad    ->push_back(v_hsimInfoConeNNeutHad  );
00861       t_v_hsimInfoConeNCharHad    ->push_back(v_hsimInfoConeNCharHad  );
00862       t_v_hsimInfoConeNPhoton     ->push_back(v_hsimInfoConeNPhoton    );
00863       t_v_hsimInfoConeNRest       ->push_back(v_hsimInfoConeNRest     );  
00864 
00865       t_v_hsimCone    ->push_back(v_hsimCone    );
00866       t_v_hCone       ->push_back(v_hCone       );
00867       t_v_nRecHitsCone->push_back(v_nRecHitsCone);
00868       t_v_nSimHitsCone->push_back(v_nSimHitsCone);
00869     }
00870 
00871     
00872     t_v_distFromHotCell->push_back(v_distFromHotCell);
00873     t_v_ietaFromHotCell->push_back(v_ietaFromHotCell);
00874     t_v_iphiFromHotCell->push_back(v_iphiFromHotCell);
00875     
00876     t_v_RH_r26_ieta    ->push_back(v_RH_r26_ieta);
00877     t_v_RH_r26_iphi    ->push_back(v_RH_r26_iphi);
00878     t_v_RH_r26_ene     ->push_back(v_RH_r26_ene);
00879     t_v_RH_r44_ieta    ->push_back(v_RH_r44_ieta);
00880     t_v_RH_r44_iphi    ->push_back(v_RH_r44_iphi);
00881     t_v_RH_r44_ene     ->push_back(v_RH_r44_ene);
00882 
00883     
00884 
00885     t_v_hlTriggers    ->push_back(v_hlTriggers);
00886     t_hltHB ->push_back(hltHB);
00887     t_hltHE ->push_back(hltHE);
00888     t_hltL1Jet15                 ->push_back(hltL1Jet15                 );
00889     t_hltJet30                   ->push_back(hltJet30                   );
00890     t_hltJet50                   ->push_back(hltJet50                   );
00891     t_hltJet80                   ->push_back(hltJet80                   );
00892     t_hltJet110                  ->push_back(hltJet110                  );
00893     t_hltJet140                  ->push_back(hltJet140                  );
00894     t_hltJet180                  ->push_back(hltJet180                  );
00895     t_hltL1SingleEG5             ->push_back(hltL1SingleEG5             );
00896     t_hltZeroBias                ->push_back(hltZeroBias                );
00897     t_hltMinBiasHcal             ->push_back(hltMinBiasHcal             );
00898     t_hltMinBiasEcal             ->push_back(hltMinBiasEcal             );
00899     t_hltMinBiasPixel            ->push_back(hltMinBiasPixel            );
00900     t_hltSingleIsoTau30_Trk5     ->push_back(hltSingleIsoTau30_Trk5     );
00901     t_hltDoubleLooseIsoTau15_Trk5->push_back(hltDoubleLooseIsoTau15_Trk5);
00902 
00903     t_irun->push_back(irun);
00904     t_ievt->push_back(ievt);
00905     t_ilum->push_back(ilum);
00906     
00907     nTRK++;
00908     
00909     
00910   } // Loop over track collection
00911   
00912     //  std::cout << "nEVT = " << nEVT << std::endl;
00913   
00914   ntp->Fill();
00915   nEVT++;
00916   
00917   
00918   delete associate;
00919 }
00920   
00921 
00922 
00923 void IsolatedTracksCone::beginJob(const edm::EventSetup&) {
00924 
00925   //   hbScale = 120.0;
00926   //   heScale = 135.0;
00927   nEVT=0;
00928   nEVT_failL1=0;
00929   nTRK=0;
00930   
00931   double tempgen_TH[22] = { 0.0,  1.0,  2.0,  3.0,  4.0,  
00932                             5.0,  6.0,  7.0,  8.0,  9.0, 
00933                             10.0, 12.0, 15.0, 20.0, 25.0, 
00934                             30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 100};
00935   for(int i=0; i<22; i++)  genPartPBins[i]  = tempgen_TH[i];
00936 
00937 
00938   double tempgen_Eta[5] = {0.0, 0.5, 1.1, 1.7, 2.0};
00939   for(int i=0; i<5; i++) genPartEtaBins[i] = tempgen_Eta[i];
00940 
00941   t_v_hnNearTRKs           = new std::vector<std::vector<int> >   (); 
00942   t_v_hnLayers_maxNearP    = new std::vector<std::vector<int> >   (); 
00943   t_v_htrkQual_maxNearP    = new std::vector<std::vector<int> >  (); 
00944   t_v_hmaxNearP_goodTrk    = new std::vector<std::vector<double> >();
00945   t_v_hmaxNearP            = new std::vector<std::vector<double> >();
00946                                                             
00947   t_v_cone_hnNearTRKs           = new std::vector<std::vector<int> >   (); 
00948   t_v_cone_hnLayers_maxNearP    = new std::vector<std::vector<int> >   (); 
00949   t_v_cone_htrkQual_maxNearP    = new std::vector<std::vector<int> >  (); 
00950   t_v_cone_hmaxNearP_goodTrk    = new std::vector<std::vector<double> >();
00951   t_v_cone_hmaxNearP            = new std::vector<std::vector<double> >();
00952                                                             
00953   //  t_hScale             = new std::vector<double>();
00954   t_trkNOuterHits      = new std::vector<double>();
00955   t_trkNLayersCrossed  = new std::vector<double>();
00956   t_dtFromLeadJet      = new std::vector<double>();
00957   t_trkP               = new std::vector<double>();
00958   t_trkPt              = new std::vector<double>();
00959   t_trkEta             = new std::vector<double>();
00960   t_trkPhi             = new std::vector<double>();
00961 
00962   t_e3x3               = new std::vector<double>();
00963   t_v_eDR              = new std::vector<std::vector<double> >();
00964   t_v_eMipDR           = new std::vector<std::vector<double> >();
00965 
00966   t_h3x3               = new std::vector<double>();
00967   t_h5x5               = new std::vector<double>();
00968 
00969   t_nRH_h3x3           = new std::vector<double>();
00970   t_nRH_h5x5           = new std::vector<double>();
00971 
00972   if (doMC) {
00973     t_simP               = new std::vector<double>();
00974     t_hsim3x3            = new std::vector<double>();
00975     t_hsim5x5            = new std::vector<double>();
00976 
00977     t_hsim3x3Matched     = new std::vector<double>();
00978     t_hsim5x5Matched     = new std::vector<double>();
00979     t_hsim3x3Rest        = new std::vector<double>();
00980     t_hsim5x5Rest        = new std::vector<double>();
00981     t_hsim3x3Photon      = new std::vector<double>();
00982     t_hsim5x5Photon      = new std::vector<double>();
00983     t_hsim3x3NeutHad     = new std::vector<double>();
00984     t_hsim5x5NeutHad     = new std::vector<double>();
00985     t_hsim3x3CharHad     = new std::vector<double>();
00986     t_hsim5x5CharHad     = new std::vector<double>();
00987     t_hsim3x3PdgMatched  = new std::vector<double>();
00988     t_hsim5x5PdgMatched  = new std::vector<double>();
00989     t_hsim3x3Total       = new std::vector<double>();
00990     t_hsim5x5Total       = new std::vector<double>();
00991 
00992     t_hsim3x3NMatched  = new std::vector<int>();
00993     t_hsim3x3NTotal    = new std::vector<int>();
00994     t_hsim3x3NNeutHad  = new std::vector<int>();
00995     t_hsim3x3NCharHad  = new std::vector<int>();
00996     t_hsim3x3NPhoton   = new std::vector<int>();
00997     t_hsim3x3NRest     = new std::vector<int>();
00998 
00999     t_hsim5x5NMatched  = new std::vector<int>();
01000     t_hsim5x5NTotal    = new std::vector<int>();
01001     t_hsim5x5NNeutHad  = new std::vector<int>();
01002     t_hsim5x5NCharHad  = new std::vector<int>();
01003     t_hsim5x5NPhoton   = new std::vector<int>();
01004     t_hsim5x5NRest     = new std::vector<int>();
01005 
01006     t_trkHcalEne         = new std::vector<double>();
01007     t_trkEcalEne         = new std::vector<double>();
01008   }
01009 
01010   t_distFromHotCell_h3x3  = new std::vector<double>();
01011   t_ietaFromHotCell_h3x3  = new std::vector<int>();
01012   t_iphiFromHotCell_h3x3  = new std::vector<int>();
01013   t_distFromHotCell_h5x5  = new std::vector<double>();
01014   t_ietaFromHotCell_h5x5  = new std::vector<int>();
01015   t_iphiFromHotCell_h5x5  = new std::vector<int>();
01016 
01017   if (doMC) {
01018     t_v_hsimInfoConeMatched   = new std::vector<std::vector<double> >();
01019     t_v_hsimInfoConeRest      = new std::vector<std::vector<double> >();
01020     t_v_hsimInfoConePhoton     = new std::vector<std::vector<double> >();
01021     t_v_hsimInfoConeNeutHad   = new std::vector<std::vector<double> >();
01022     t_v_hsimInfoConeCharHad   = new std::vector<std::vector<double> >();
01023     t_v_hsimInfoConePdgMatched= new std::vector<std::vector<double> >();
01024     t_v_hsimInfoConeTotal     = new std::vector<std::vector<double> >();
01025 
01026     t_v_hsimInfoConeNMatched  = new std::vector<std::vector<int> >();
01027     t_v_hsimInfoConeNTotal    = new std::vector<std::vector<int> >();
01028     t_v_hsimInfoConeNNeutHad  = new std::vector<std::vector<int> >();
01029     t_v_hsimInfoConeNCharHad  = new std::vector<std::vector<int> >();
01030     t_v_hsimInfoConeNPhoton    = new std::vector<std::vector<int> >();
01031     t_v_hsimInfoConeNRest     = new std::vector<std::vector<int> >();
01032 
01033     t_v_hsimCone    = new std::vector<std::vector<double> >();
01034   }
01035 
01036   t_v_hCone       = new std::vector<std::vector<double> >();
01037   t_v_nRecHitsCone= new std::vector<std::vector<int> >();
01038   t_v_nSimHitsCone= new std::vector<std::vector<int> >();
01039 
01040   t_v_distFromHotCell= new std::vector<std::vector<double> >();
01041   t_v_ietaFromHotCell= new std::vector<std::vector<int> >();
01042   t_v_iphiFromHotCell= new std::vector<std::vector<int> >();
01043 
01044   t_v_RH_h3x3_ieta = new std::vector<std::vector<int> >();
01045   t_v_RH_h3x3_iphi = new std::vector<std::vector<int> >();
01046   t_v_RH_h3x3_ene  = new std::vector<std::vector<double> >();
01047   t_v_RH_h5x5_ieta = new std::vector<std::vector<int> >();
01048   t_v_RH_h5x5_iphi = new std::vector<std::vector<int> >();
01049   t_v_RH_h5x5_ene  = new std::vector<std::vector<double> >();
01050   t_v_RH_r26_ieta  = new std::vector<std::vector<int> >();
01051   t_v_RH_r26_iphi  = new std::vector<std::vector<int> >();
01052   t_v_RH_r26_ene   = new std::vector<std::vector<double> >();
01053   t_v_RH_r44_ieta  = new std::vector<std::vector<int> >();
01054   t_v_RH_r44_iphi  = new std::vector<std::vector<int> >();
01055   t_v_RH_r44_ene   = new std::vector<std::vector<double> >();
01056 
01057 
01058   t_v_hlTriggers    = new std::vector<std::vector<int> >();
01059 
01060   t_hltHE                         = new std::vector<int>();
01061   t_hltHB                         = new std::vector<int>();
01062   t_hltL1Jet15                    = new std::vector<int>();
01063   t_hltJet30                      = new std::vector<int>();
01064   t_hltJet50                      = new std::vector<int>();
01065   t_hltJet80                      = new std::vector<int>();
01066   t_hltJet110                     = new std::vector<int>();
01067   t_hltJet140                     = new std::vector<int>();
01068   t_hltJet180                     = new std::vector<int>();
01069   t_hltL1SingleEG5                = new std::vector<int>();
01070   t_hltZeroBias                   = new std::vector<int>();
01071   t_hltMinBiasHcal                = new std::vector<int>();
01072   t_hltMinBiasEcal                = new std::vector<int>();
01073   t_hltMinBiasPixel               = new std::vector<int>();
01074   t_hltSingleIsoTau30_Trk5        = new std::vector<int>();
01075   t_hltDoubleLooseIsoTau15_Trk5   = new std::vector<int>();
01076   
01077 
01078   t_irun = new std::vector<unsigned int>();
01079   t_ievt = new std::vector<unsigned int>();
01080   t_ilum = new std::vector<unsigned int>();
01081 
01082   BuildTree();
01083 }
01084 
01085 void IsolatedTracksCone::clearTrackVectors() {
01086 
01087   t_v_hnNearTRKs          ->clear();   
01088   t_v_hnLayers_maxNearP   ->clear();   
01089   t_v_htrkQual_maxNearP   ->clear();   
01090   t_v_hmaxNearP_goodTrk   ->clear();
01091   t_v_hmaxNearP           ->clear();
01092 
01093   t_v_cone_hnNearTRKs          ->clear();   
01094   t_v_cone_hnLayers_maxNearP   ->clear();   
01095   t_v_cone_htrkQual_maxNearP   ->clear();   
01096   t_v_cone_hmaxNearP_goodTrk   ->clear();
01097   t_v_cone_hmaxNearP           ->clear();
01098 
01099   //  t_hScale             ->clear();
01100   t_trkNOuterHits      ->clear();
01101   t_trkNLayersCrossed  ->clear();
01102   t_dtFromLeadJet      ->clear();
01103   t_trkP               ->clear();
01104   t_trkPt              ->clear();
01105   t_trkEta             ->clear();
01106   t_trkPhi             ->clear();
01107   t_e3x3               ->clear();
01108   t_v_eDR              ->clear();
01109   t_v_eMipDR           ->clear();
01110   t_h3x3               ->clear();
01111   t_h5x5               ->clear();
01112   t_nRH_h3x3           ->clear();
01113   t_nRH_h5x5           ->clear();
01114 
01115   if (doMC) {
01116     t_simP               ->clear();
01117     t_hsim3x3            ->clear();
01118     t_hsim5x5            ->clear();
01119     t_hsim3x3Matched     ->clear();
01120     t_hsim5x5Matched     ->clear();
01121     t_hsim3x3Rest        ->clear();
01122     t_hsim5x5Rest        ->clear();
01123     t_hsim3x3Photon      ->clear();
01124     t_hsim5x5Photon      ->clear();
01125     t_hsim3x3NeutHad     ->clear();
01126     t_hsim5x5NeutHad     ->clear();
01127     t_hsim3x3CharHad     ->clear();
01128     t_hsim5x5CharHad     ->clear();
01129     t_hsim3x3PdgMatched  ->clear();
01130     t_hsim5x5PdgMatched  ->clear();
01131     t_hsim3x3Total       ->clear();
01132     t_hsim5x5Total       ->clear();
01133 
01134     t_hsim3x3NMatched    ->clear();
01135     t_hsim3x3NTotal      ->clear();
01136     t_hsim3x3NNeutHad    ->clear();
01137     t_hsim3x3NCharHad    ->clear();
01138     t_hsim3x3NPhoton      ->clear();
01139     t_hsim3x3NRest       ->clear();
01140 
01141     t_hsim5x5NMatched    ->clear();
01142     t_hsim5x5NTotal      ->clear();
01143     t_hsim5x5NNeutHad    ->clear();
01144     t_hsim5x5NCharHad    ->clear();
01145     t_hsim5x5NPhoton      ->clear();
01146     t_hsim5x5NRest       ->clear();
01147 
01148     t_trkHcalEne         ->clear();
01149     t_trkEcalEne         ->clear();
01150   }
01151 
01152   t_distFromHotCell_h3x3  ->clear();
01153   t_ietaFromHotCell_h3x3  ->clear();
01154   t_iphiFromHotCell_h3x3  ->clear();
01155   t_distFromHotCell_h5x5  ->clear();
01156   t_ietaFromHotCell_h5x5  ->clear();
01157   t_iphiFromHotCell_h5x5  ->clear();
01158 
01159   if (doMC) {
01160     t_v_hsimInfoConeMatched   ->clear();
01161     t_v_hsimInfoConeRest      ->clear();
01162     t_v_hsimInfoConePhoton     ->clear();
01163     t_v_hsimInfoConeNeutHad   ->clear();
01164     t_v_hsimInfoConeCharHad   ->clear();
01165     t_v_hsimInfoConePdgMatched->clear();  
01166     t_v_hsimInfoConeTotal     ->clear();
01167 
01168     t_v_hsimInfoConeNMatched   ->clear();
01169     t_v_hsimInfoConeNRest      ->clear();
01170     t_v_hsimInfoConeNPhoton     ->clear();
01171     t_v_hsimInfoConeNNeutHad   ->clear();
01172     t_v_hsimInfoConeNCharHad   ->clear();
01173     t_v_hsimInfoConeNTotal     ->clear();
01174 
01175     t_v_hsimCone    ->clear();
01176   }
01177 
01178   t_v_hCone       ->clear();
01179   t_v_nRecHitsCone->clear();
01180   t_v_nSimHitsCone->clear();
01181 
01182   t_v_distFromHotCell->clear();
01183   t_v_ietaFromHotCell->clear();
01184   t_v_iphiFromHotCell->clear();
01185   
01186   t_v_RH_h3x3_ieta ->clear();
01187   t_v_RH_h3x3_iphi ->clear();
01188   t_v_RH_h3x3_ene ->clear();
01189   t_v_RH_h5x5_ieta ->clear();
01190   t_v_RH_h5x5_iphi ->clear();
01191   t_v_RH_h5x5_ene ->clear();
01192   t_v_RH_r26_ieta  ->clear();
01193   t_v_RH_r26_iphi  ->clear();
01194   t_v_RH_r26_ene  ->clear();
01195   t_v_RH_r44_ieta  ->clear();
01196   t_v_RH_r44_iphi  ->clear();
01197   t_v_RH_r44_ene  ->clear();
01198 
01199   t_v_hlTriggers  ->clear();
01200   t_hltHB->clear();
01201   t_hltHE->clear();
01202   t_hltL1Jet15               ->clear();
01203   t_hltJet30                 ->clear();
01204   t_hltJet50                 ->clear();
01205   t_hltJet80                 ->clear();
01206   t_hltJet110                ->clear();
01207   t_hltJet140                ->clear();
01208   t_hltJet180                ->clear();
01209   t_hltL1SingleEG5           ->clear();
01210   t_hltZeroBias              ->clear();
01211   t_hltMinBiasHcal           ->clear();
01212   t_hltMinBiasEcal           ->clear();
01213   t_hltMinBiasPixel          ->clear();
01214   t_hltSingleIsoTau30_Trk5     ->clear();
01215   t_hltDoubleLooseIsoTau15_Trk5->clear();
01216   
01217   t_irun->clear();
01218   t_ievt->clear();
01219   t_ilum->clear();
01220 
01221 
01222 
01223 }
01224 
01225 // ----- method called once each job just after ending the event loop ----
01226 void IsolatedTracksCone::endJob() {
01227 
01228   std::cout << "Number of Events Processed " << nEVT << std::endl;
01229   
01230 }
01231 
01232 
01233 void IsolatedTracksCone::BuildTree(){
01234 
01235 
01236   hRawPt  = fs->make<TH1F>("hRawPt ", "hRawPt ",  100,  0.0, 100.0);
01237   hRawP   = fs->make<TH1F>("hRawP  ", "hRawP  ",  100,  0.0, 100.0);
01238   hRawEta = fs->make<TH1F>("hRawEta", "hRawEta",   15,  0.0,   3.0);
01239   hRawPhi = fs->make<TH1F>("hRawPhi", "hRawPhi",  100, -3.2,   3.2);
01240 
01241   ntp = fs->make<TTree>("ntp", "ntp");
01242 
01243   
01244   // Counters
01245   ntp->Branch("nEVT"        , &nEVT        , "nEVT/I"        );
01246   ntp->Branch("leadL1JetPT" , &leadL1JetPT , "leadL1JetPT/D" );
01247   ntp->Branch("leadL1JetEta", &leadL1JetEta, "leadL1JetEta/D");
01248   ntp->Branch("leadL1JetPhi", &leadL1JetPhi, "leadL1JetPhi/D");
01249   ntp->Branch("nTRK",         &nTRK,         "nTRK/I");
01250   ntp->Branch("nRawTRK"            , &nRawTRK            ,"nRawTRK/I"            );
01251   ntp->Branch("nFailHighPurityQaul", &nFailHighPurityQaul,"nFailHighPurityQaul/I");
01252   ntp->Branch("nFailPt"            , &nFailPt            ,"nFailPt/I"            );
01253   ntp->Branch("nFailEta"           , &nFailEta           ,"nFailEta/I"           );
01254   ntp->Branch("nMissEcal"          , &nMissEcal          ,"nMissEcal/I"          );
01255   ntp->Branch("nMissHcal"          , &nMissHcal          ,"nMissHcal/I"          );
01256 
01257   ntp->Branch("hnNearTRKs"          ,"vector<vector<int> >   ",&t_v_hnNearTRKs          );
01258   ntp->Branch("hnLayers_maxNearP"   ,"vector<vector<int> >   ",&t_v_hnLayers_maxNearP   );
01259   ntp->Branch("htrkQual_maxNearP"   ,"vector<vector<int> >   ",&t_v_htrkQual_maxNearP   );
01260   ntp->Branch("hmaxNearP_goodTrk"   ,"vector<vector<double> >",&t_v_hmaxNearP_goodTrk   );
01261   ntp->Branch("hmaxNearP"           ,"vector<vector<double> >",&t_v_hmaxNearP           );
01262 
01263   ntp->Branch("cone_hnNearTRKs"       ,"vector<vector<int> >   ",&t_v_cone_hnNearTRKs       );
01264   ntp->Branch("cone_hnLayers_maxNearP","vector<vector<int> >   ",&t_v_cone_hnLayers_maxNearP);
01265   ntp->Branch("cone_htrkQual_maxNearP","vector<vector<int> >   ",&t_v_cone_htrkQual_maxNearP);
01266   ntp->Branch("cone_hmaxNearP_goodTrk","vector<vector<double> >",&t_v_cone_hmaxNearP_goodTrk);
01267   ntp->Branch("cone_hmaxNearP"        ,"vector<vector<double> >",&t_v_cone_hmaxNearP        );
01268                                                                     
01269   //  ntp->Branch("hScale"           , "vector<double>", &t_hScale           );
01270   ntp->Branch("trkNOuterHits"    , "vector<double>", &t_trkNOuterHits    );
01271   ntp->Branch("trkNLayersCrossed", "vector<double>", &t_trkNLayersCrossed);
01272   ntp->Branch("drFromLeadJet"    , "vector<double>", &t_dtFromLeadJet    );
01273   ntp->Branch("trkP"             , "vector<double>", &t_trkP             );
01274   ntp->Branch("trkPt"            , "vector<double>", &t_trkPt            );
01275   ntp->Branch("trkEta"           , "vector<double>", &t_trkEta           );
01276   ntp->Branch("trkPhi"           , "vector<double>", &t_trkPhi           );
01277   ntp->Branch("e3x3"             , "vector<double>", &t_e3x3             );
01278 
01279   ntp->Branch("e3x3"          , "vector<double>"         , &t_e3x3 );
01280   ntp->Branch("v_eDR"         , "vector<vector<double> >", &t_v_eDR);
01281   ntp->Branch("v_eMipDR"      , "vector<vector<double> >", &t_v_eMipDR);
01282 
01283   ntp->Branch("h3x3"             , "vector<double>", &t_h3x3             );
01284   ntp->Branch("h5x5"             , "vector<double>", &t_h5x5             );
01285   ntp->Branch("nRH_h3x3"         , "vector<double>", &t_nRH_h3x3         );
01286   ntp->Branch("nRH_h5x5"         , "vector<double>", &t_nRH_h5x5         );
01287 
01288   if (doMC) {
01289     ntp->Branch("simP"             , "vector<double>", &t_simP             );
01290     ntp->Branch("hsim3x3"          , "vector<double>", &t_hsim3x3          );
01291     ntp->Branch("hsim5x5"          , "vector<double>", &t_hsim5x5          );
01292 
01293     ntp->Branch("hsim3x3Matched"   , "vector<double>", &t_hsim3x3Matched   );
01294     ntp->Branch("hsim5x5Matched"   , "vector<double>", &t_hsim5x5Matched   );
01295     ntp->Branch("hsim3x3Rest"      , "vector<double>", &t_hsim3x3Rest      );
01296     ntp->Branch("hsim5x5Rest"      , "vector<double>", &t_hsim5x5Rest      );
01297     ntp->Branch("hsim3x3Photon"    , "vector<double>", &t_hsim3x3Photon    );
01298     ntp->Branch("hsim5x5Photon"    , "vector<double>", &t_hsim5x5Photon    );
01299     ntp->Branch("hsim3x3NeutHad"   , "vector<double>", &t_hsim3x3NeutHad   );
01300     ntp->Branch("hsim5x5NeutHad"   , "vector<double>", &t_hsim5x5NeutHad   );
01301     ntp->Branch("hsim3x3CharHad"   , "vector<double>", &t_hsim3x3CharHad   );
01302     ntp->Branch("hsim5x5CharHad"   , "vector<double>", &t_hsim5x5CharHad   );
01303     ntp->Branch("hsim3x3PdgMatched", "vector<double>", &t_hsim3x3PdgMatched);
01304     ntp->Branch("hsim5x5PdgMatched", "vector<double>", &t_hsim5x5PdgMatched);
01305     ntp->Branch("hsim3x3Total"     , "vector<double>", &t_hsim3x3Total     );
01306     ntp->Branch("hsim5x5Total"     , "vector<double>", &t_hsim5x5Total     );
01307 
01308     ntp->Branch("hsim3x3NMatched"   , "vector<int>", &t_hsim3x3NMatched   );
01309     ntp->Branch("hsim3x3NRest"      , "vector<int>", &t_hsim3x3NRest      );
01310     ntp->Branch("hsim3x3NPhoton"    , "vector<int>", &t_hsim3x3NPhoton    );
01311     ntp->Branch("hsim3x3NNeutHad"   , "vector<int>", &t_hsim3x3NNeutHad   );
01312     ntp->Branch("hsim3x3NCharHad"   , "vector<int>", &t_hsim3x3NCharHad   );
01313     ntp->Branch("hsim3x3NTotal"     , "vector<int>", &t_hsim3x3NTotal     );
01314 
01315     ntp->Branch("hsim5x5NMatched"   , "vector<int>", &t_hsim5x5NMatched   );
01316     ntp->Branch("hsim5x5NRest"      , "vector<int>", &t_hsim5x5NRest      );
01317     ntp->Branch("hsim5x5NPhoton"    , "vector<int>", &t_hsim5x5NPhoton    );
01318     ntp->Branch("hsim5x5NNeutHad"   , "vector<int>", &t_hsim5x5NNeutHad   );
01319     ntp->Branch("hsim5x5NCharHad"   , "vector<int>", &t_hsim5x5NCharHad   );
01320     ntp->Branch("hsim5x5NTotal"     , "vector<int>", &t_hsim5x5NTotal     );
01321 
01322     ntp->Branch("trkHcalEne"       , "vector<double>", &t_trkHcalEne       );
01323     ntp->Branch("trkEcalEne"       , "vector<double>", &t_trkEcalEne       );
01324   }
01325 
01326   ntp->Branch("distFromHotCell_h3x3", "vector<double>", &t_distFromHotCell_h3x3);
01327   ntp->Branch("ietaFromHotCell_h3x3", "vector<int>", &t_ietaFromHotCell_h3x3);
01328   ntp->Branch("iphiFromHotCell_h3x3", "vector<int>", &t_iphiFromHotCell_h3x3);
01329   ntp->Branch("distFromHotCell_h5x5", "vector<double>", &t_distFromHotCell_h5x5);
01330   ntp->Branch("ietaFromHotCell_h5x5", "vector<int>", &t_ietaFromHotCell_h5x5);
01331   ntp->Branch("iphiFromHotCell_h5x5", "vector<int>", &t_iphiFromHotCell_h5x5);
01332 
01333   if (doMC) {
01334     ntp->Branch("v_hsimInfoConeMatched"   ,"vector<vector<double> >",&t_v_hsimInfoConeMatched   );
01335     ntp->Branch("v_hsimInfoConeRest"      ,"vector<vector<double> >",&t_v_hsimInfoConeRest      );
01336     ntp->Branch("v_hsimInfoConePhoton"     ,"vector<vector<double> >",&t_v_hsimInfoConePhoton     );
01337     ntp->Branch("v_hsimInfoConeNeutHad"   ,"vector<vector<double> >",&t_v_hsimInfoConeNeutHad   );
01338     ntp->Branch("v_hsimInfoConeCharHad"   ,"vector<vector<double> >",&t_v_hsimInfoConeCharHad   );
01339     ntp->Branch("v_hsimInfoConePdgMatched","vector<vector<double> >",&t_v_hsimInfoConePdgMatched);
01340     ntp->Branch("v_hsimInfoConeTotal"     ,"vector<vector<double> >",&t_v_hsimInfoConeTotal     );
01341 
01342     ntp->Branch("v_hsimInfoConeNMatched"  ,"vector<vector<int> >"   ,&t_v_hsimInfoConeNMatched   );
01343     ntp->Branch("v_hsimInfoConeNRest"     ,"vector<vector<int> >"   ,&t_v_hsimInfoConeNRest      );
01344     ntp->Branch("v_hsimInfoConeNPhoton"    ,"vector<vector<int> >"   ,&t_v_hsimInfoConeNPhoton     );
01345     ntp->Branch("v_hsimInfoConeNNeutHad"  ,"vector<vector<int> >"   ,&t_v_hsimInfoConeNNeutHad   );
01346     ntp->Branch("v_hsimInfoConeNCharHad"  ,"vector<vector<int> >"   ,&t_v_hsimInfoConeNCharHad   );
01347     ntp->Branch("v_hsimInfoConeNTotal"    ,"vector<vector<int> >"   ,&t_v_hsimInfoConeNTotal     );
01348 
01349     ntp->Branch("v_hsimCone"     ,"vector<vector<double> >",&t_v_hsimCone    );
01350   }
01351 
01352   ntp->Branch("v_hCone"        ,"vector<vector<double> >",&t_v_hCone       );
01353   ntp->Branch("v_nRecHitsCone" ,"vector<vector<int> >"   ,&t_v_nRecHitsCone);
01354   ntp->Branch("v_nSimHitsCone" ,"vector<vector<int> >"   ,&t_v_nSimHitsCone);
01355 
01356   ntp->Branch("v_distFromHotCell"        ,"vector<vector<double> >",&t_v_distFromHotCell       );
01357   ntp->Branch("v_ietaFromHotCell"        ,"vector<vector<int> >",&t_v_ietaFromHotCell       );
01358   ntp->Branch("v_iphiFromHotCell"        ,"vector<vector<int> >",&t_v_iphiFromHotCell       );
01359 
01360   ntp->Branch("v_RH_h3x3_ieta" ,"vector<vector<int> >"   ,&t_v_RH_h3x3_ieta);
01361   ntp->Branch("v_RH_h3x3_iphi" ,"vector<vector<int> >"   ,&t_v_RH_h3x3_iphi);
01362   ntp->Branch("v_RH_h3x3_ene"  ,"vector<vector<double> >",&t_v_RH_h3x3_ene );
01363   ntp->Branch("v_RH_h5x5_ieta" ,"vector<vector<int> >"   ,&t_v_RH_h5x5_ieta);
01364   ntp->Branch("v_RH_h5x5_iphi" ,"vector<vector<int> >"   ,&t_v_RH_h5x5_iphi);
01365   ntp->Branch("v_RH_h5x5_ene"  ,"vector<vector<double> >",&t_v_RH_h5x5_ene );
01366   ntp->Branch("v_RH_r26_ieta"  ,"vector<vector<int> >"   ,&t_v_RH_r26_ieta );
01367   ntp->Branch("v_RH_r26_iphi"  ,"vector<vector<int> >"   ,&t_v_RH_r26_iphi );
01368   ntp->Branch("v_RH_r26_ene"   ,"vector<vector<double> >",&t_v_RH_r26_ene  );
01369   ntp->Branch("v_RH_r44_ieta"  ,"vector<vector<int> >"   ,&t_v_RH_r44_ieta );
01370   ntp->Branch("v_RH_r44_iphi"  ,"vector<vector<int> >"   ,&t_v_RH_r44_iphi );
01371   ntp->Branch("v_RH_r44_ene"   ,"vector<vector<double> >",&t_v_RH_r44_ene  );
01372 
01373   ntp->Branch("v_hlTriggers" ,"vector<vector<int> >",&t_v_hlTriggers);
01374   ntp->Branch("v_hltHB" ,"vector<int>",&t_hltHB);
01375   ntp->Branch("v_hltHE" ,"vector<int>",&t_hltHE);
01376   ntp->Branch("v_hltL1Jet15" ,"vector<int>",&t_hltL1Jet15                    );
01377   ntp->Branch("v_hltJet30" ,"vector<int>",&t_hltJet30                );
01378   ntp->Branch("v_hltJet50" ,"vector<int>",&t_hltJet50                );
01379   ntp->Branch("v_hltJet80" ,"vector<int>",&t_hltJet80                );
01380   ntp->Branch("v_hltJet110" ,"vector<int>",&t_hltJet110              );
01381   ntp->Branch("v_hltJet140" ,"vector<int>",&t_hltJet140              );
01382   ntp->Branch("v_hltJet180" ,"vector<int>",&t_hltJet180              );
01383   ntp->Branch("v_hltL1SingleEG5" ,"vector<int>",&t_hltL1SingleEG5            );
01384   ntp->Branch("v_hltZeroBias" ,"vector<int>",&t_hltZeroBias                  );
01385   ntp->Branch("v_hltMinBiasHcal" ,"vector<int>",&t_hltMinBiasHcal            );
01386   ntp->Branch("v_hltMinBiasEcal" ,"vector<int>",&t_hltMinBiasEcal            );
01387   ntp->Branch("v_hltMinBiasPixel" ,"vector<int>",&t_hltMinBiasPixel          );
01388   ntp->Branch("v_hltSingleIsoTau30_Trk5" ,"vector<int>",&t_hltSingleIsoTau30_Trk5     );
01389   ntp->Branch("v_hltDoubleLooseIsoTau15_Trk5" ,"vector<int>",&t_hltDoubleLooseIsoTau15_Trk5);
01390 
01391   ntp->Branch("irun" ,"vector<unsigned int>", &t_irun);
01392   ntp->Branch("ievt" ,"vector<unsigned int>", &t_ievt);
01393   ntp->Branch("ilum" ,"vector<unsigned int>", &t_ilum);
01394    
01395 }
01396 
01397 
01398 void IsolatedTracksCone::printTrack(const reco::Track* pTrack) {
01399   
01400   std::string theTrackQuality = "highPurity";
01401   reco::TrackBase::TrackQuality trackQuality_ = reco::TrackBase::qualityByName(theTrackQuality);
01402 
01403   std::cout << " Reference Point " << pTrack->referencePoint() <<"\n"
01404             << " TrackMmentum " << pTrack->momentum()
01405             << " (pt,eta,phi)(" << pTrack->pt()<<","<<pTrack->eta()<<","<<pTrack->phi()<<")"
01406             << " p " << pTrack->p() << "\n"
01407             << " Normalized chi2 " << pTrack->normalizedChi2() <<"  charge " << pTrack->charge()
01408             << " qoverp() " << pTrack->qoverp() <<"+-" << pTrack->qoverpError()
01409             << " d0 " << pTrack->d0() << "\n"
01410             << " NValidHits " << pTrack->numberOfValidHits() << "  NLostHits " << pTrack->numberOfLostHits()
01411             << " TrackQuality " << pTrack->qualityName(trackQuality_) << " " << pTrack->quality(trackQuality_) 
01412             << std::endl;
01413   
01414   if( printTrkHitPattern_ ) {
01415     const reco::HitPattern& p = pTrack->hitPattern();
01416     
01417     for (int i=0; i<p.numberOfHits(); i++) {
01418       p.printHitPattern(i, std::cout);
01419     }
01420   }
01421 
01422 }
01423 
01424 double IsolatedTracksCone::DeltaPhi(double v1, double v2) {
01425   // Computes the correctly normalized phi difference
01426   // v1, v2 = phi of object 1 and 2
01427   
01428   double pi    = 3.141592654;
01429   double twopi = 6.283185307;
01430   
01431   double diff = std::abs(v2 - v1);
01432   double corr = twopi - diff;
01433   if (diff < pi){ return diff;} else { return corr;} 
01434 }
01435 
01436 
01437 double IsolatedTracksCone::DeltaR(double eta1, double phi1, 
01438                                   double eta2, double phi2) {
01439   double deta = eta1 - eta2;
01440   double dphi = DeltaPhi(phi1, phi2);
01441   return std::sqrt(deta*deta + dphi*dphi);
01442 }
01443 
01444 
01445 //define this as a plug-in
01446 DEFINE_FWK_MODULE(IsolatedTracksCone);
01447           
01448           
01449           
01450           
01451           
01452 
01453           
01454           
01455 
01456 
01457 
01458 
01459 
01460 
01461
CMSSW_4_4_3_patch1/src/Calibration/IsolatedParticles/plugins/IsolatedTracksCone.cc
