Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 11:59:16 

0001 // system include files
0002 #include <memory>
0003 #include <string>
0004 #include <vector>
0005 
0006 // Root objects
0007 #include "TROOT.h"
0008 #include "TSystem.h"
0009 #include "TFile.h"
0010 #include "TProfile.h"
0011 #include "TDirectory.h"
0012 #include "TTree.h"
0013 #include "TH1F.h"
0014 #include "TLorentzVector.h"
0015 #include "TInterpreter.h"
0016 
0017 //Tracks
0018 #include "DataFormats/TrackReco/interface/Track.h"
0019 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0020 #include "DataFormats/TrackReco/interface/HitPattern.h"
0021 #include "DataFormats/TrackReco/interface/TrackBase.h"
0022 
0023 // Vertices
0024 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0025 #include "DataFormats/VertexReco/interface/Vertex.h"
0026 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0027 
0028 //Triggers
0029 #include "DataFormats/Common/interface/TriggerResults.h"
0030 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0031 #include "DataFormats/HLTReco/interface/TriggerObject.h"
0032 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
0033 #include "DataFormats/Luminosity/interface/LumiDetails.h"
0034 #include "DataFormats/Math/interface/LorentzVector.h"
0035 #include "DataFormats/Math/interface/deltaR.h"
0036 #include "DataFormats/Math/interface/deltaPhi.h"
0037 
0038 #include "FWCore/Framework/interface/Frameworkfwd.h"
0039 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0040 #include "FWCore/Framework/interface/Event.h"
0041 #include "FWCore/Framework/interface/MakerMacros.h"
0042 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0043 #include "FWCore/ServiceRegistry/interface/Service.h"
0044 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0045 #include "FWCore/Common/interface/TriggerNames.h"
0046 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0047 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0048 
0049 #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h"
0050 #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h"
0051 #include "Calibration/IsolatedParticles/interface/eCone.h"
0052 #include "Calibration/IsolatedParticles/interface/eECALMatrix.h"
0053 #include "Calibration/IsolatedParticles/interface/TrackSelection.h"
0054 
0055 #include "MagneticField/Engine/interface/MagneticField.h"
0056 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0057 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0058 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0059 #include "Geometry/Records/interface/CaloTopologyRecord.h"
0060 #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
0061 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0062 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0063 #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h"
0064 #include "Geometry/CaloTopology/interface/EcalTrigTowerConstituentsMap.h"
0065 #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
0066 #include "SimDataFormats/PileupSummaryInfo/interface/PileupSummaryInfo.h"
0067 #include "DataFormats/JetReco/interface/PFJet.h"
0068 
0069 class IsoTrackCalibration : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> {
0070 public:
0071   explicit IsoTrackCalibration(const edm::ParameterSet &);
0072   ~IsoTrackCalibration() override;
0073 
0074   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0075 
0076 private:
0077   void analyze(edm::Event const &, edm::EventSetup const &) override;
0078   void beginJob() override;
0079   void beginRun(edm::Run const &, edm::EventSetup const &) override;
0080   void endRun(edm::Run const &, edm::EventSetup const &) override;
0081   virtual void beginLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &);
0082   virtual void endLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &);
0083 
0084   edm::Service<TFileService> fs_;
0085   HLTConfigProvider hltConfig_;
0086   std::vector<std::string> trigNames_, HLTNames_;
0087   int verbosity_;
0088   spr::trackSelectionParameters selectionParameters_;
0089   std::string theTrackQuality_;
0090   double a_mipR_, a_coneR_, a_charIsoR_;
0091   bool isMC_, isQCD_, isAOD_;
0092   double constTrackPt_, slopeTrackPt_;
0093   double maxEcalEnr_;
0094   double maxNeighborTrackEnr_;
0095   int nRun_;
0096   edm::InputTag triggerEvent_, theTriggerResultsLabel_;
0097   edm::EDGetTokenT<trigger::TriggerEvent> tok_trigEvt_;
0098   edm::EDGetTokenT<edm::TriggerResults> tok_trigRes_;
0099 
0100   edm::EDGetTokenT<reco::TrackCollection> tok_genTrack_;
0101   edm::EDGetTokenT<reco::VertexCollection> tok_recVtx_;
0102   edm::EDGetTokenT<reco::BeamSpot> tok_bs_;
0103   edm::EDGetTokenT<EcalRecHitCollection> tok_EB_;
0104   edm::EDGetTokenT<EcalRecHitCollection> tok_EE_;
0105   edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_;
0106   edm::EDGetTokenT<GenEventInfoProduct> tok_ew_;
0107 
0108   edm::ESGetToken<CaloGeometry, CaloGeometryRecord> tok_geom_;
0109   edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> tok_magField_;
0110 
0111   TTree *tree;
0112   int t_Run, t_Event, t_nVtx, t_nTrk, t_ieta;
0113   double t_EventWeight, t_p, t_pt, t_phi;
0114   double t_eHcal, t_eHcal10, t_eHcal30;
0115   double t_eMipDR, t_hmaxNearP;
0116   bool t_selectTk, t_qltyMissFlag, t_qltyPVFlag;
0117   std::vector<unsigned int> *t_DetIds, *t_DetIds1, *t_DetIds3;
0118   std::vector<double> *t_HitEnergies, *t_HitEnergies1, *t_HitEnergies3;
0119 
0120   TH1F *h_nTrk, *h_nVtx;
0121   TProfile *h_RecHit_iEta, *h_RecHit_num;
0122   TH1I *h_iEta, *h_tketa0[5], *h_tketa1[5], *h_tketa2[5];
0123   TH1I *h_tketa3[5], *h_tketa4[5], *h_tketa5[5];
0124   TH1F *h_Rechit_E, *h_jetp;
0125   TH1F *h_jetpt[4];
0126   TH1I *h_tketav1[5][6], *h_tketav2[5][6];
0127 };
0128 
0129 IsoTrackCalibration::IsoTrackCalibration(const edm::ParameterSet &iConfig) : nRun_(0) {
0130   usesResource("TFileService");
0131   //now do whatever initialization is needed
0132   const double isolationRadius(28.9);
0133   verbosity_ = iConfig.getUntrackedParameter<int>("Verbosity", 0);
0134   trigNames_ = iConfig.getUntrackedParameter<std::vector<std::string> >("Triggers");
0135   theTrackQuality_ = iConfig.getUntrackedParameter<std::string>("TrackQuality", "highPurity");
0136   reco::TrackBase::TrackQuality trackQuality_ = reco::TrackBase::qualityByName(theTrackQuality_);
0137   constTrackPt_ = iConfig.getUntrackedParameter<double>("ConstTrackPt", 10.0);
0138   slopeTrackPt_ = iConfig.getUntrackedParameter<double>("SlopeTrackPt", 0.16);
0139   selectionParameters_.minPt = constTrackPt_;
0140   selectionParameters_.minQuality = trackQuality_;
0141   selectionParameters_.maxDxyPV = iConfig.getUntrackedParameter<double>("MaxDxyPV", 0.02);
0142   selectionParameters_.maxDzPV = iConfig.getUntrackedParameter<double>("MaxDzPV", 0.02);
0143   selectionParameters_.maxChi2 = iConfig.getUntrackedParameter<double>("MaxChi2", 5.0);
0144   selectionParameters_.maxDpOverP = iConfig.getUntrackedParameter<double>("MaxDpOverP", 0.1);
0145   selectionParameters_.minOuterHit = iConfig.getUntrackedParameter<int>("MinOuterHit", 4);
0146   selectionParameters_.minLayerCrossed = iConfig.getUntrackedParameter<int>("MinLayerCrossed", 8);
0147   selectionParameters_.maxInMiss = iConfig.getUntrackedParameter<int>("MaxInMiss", 0);
0148   selectionParameters_.maxOutMiss = iConfig.getUntrackedParameter<int>("MaxOutMiss", 0);
0149   a_coneR_ = iConfig.getUntrackedParameter<double>("ConeRadius", 34.98);
0150   a_charIsoR_ = a_coneR_ + isolationRadius;
0151   a_mipR_ = iConfig.getUntrackedParameter<double>("ConeRadiusMIP", 14.0);
0152   maxEcalEnr_ = iConfig.getUntrackedParameter<double>("MaxEcalEnergyInCone", 2.5);
0153   maxNeighborTrackEnr_ = iConfig.getUntrackedParameter<double>("MaxNeighborTrackEnergy", 40.0);
0154   isMC_ = iConfig.getUntrackedParameter<bool>("IsMC", false);
0155   isQCD_ = iConfig.getUntrackedParameter<bool>("IsQCD", false);
0156   isAOD_ = iConfig.getUntrackedParameter<bool>("IsAOD", true);
0157   triggerEvent_ = edm::InputTag("hltTriggerSummaryAOD", "", "HLT");
0158   theTriggerResultsLabel_ = edm::InputTag("TriggerResults", "", "HLT");
0159 
0160   // define tokens for access
0161   tok_trigEvt_ = consumes<trigger::TriggerEvent>(triggerEvent_);
0162   tok_trigRes_ = consumes<edm::TriggerResults>(theTriggerResultsLabel_);
0163   tok_genTrack_ = consumes<reco::TrackCollection>(edm::InputTag("generalTracks"));
0164   tok_recVtx_ = consumes<reco::VertexCollection>(edm::InputTag("offlinePrimaryVertices"));
0165   tok_bs_ = consumes<reco::BeamSpot>(edm::InputTag("offlineBeamSpot"));
0166   tok_ew_ = consumes<GenEventInfoProduct>(edm::InputTag("generator"));
0167 
0168   if (isAOD_) {
0169     tok_EB_ = consumes<EcalRecHitCollection>(edm::InputTag("reducedEcalRecHitsEB"));
0170     tok_EE_ = consumes<EcalRecHitCollection>(edm::InputTag("reducedEcalRecHitsEE"));
0171     tok_hbhe_ = consumes<HBHERecHitCollection>(edm::InputTag("reducedHcalRecHits", "hbhereco"));
0172   } else {
0173     tok_EB_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit", "EcalRecHitsEB"));
0174     tok_EE_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit", "EcalRecHitsEE"));
0175     tok_hbhe_ = consumes<HBHERecHitCollection>(edm::InputTag("hbhereco"));
0176   }
0177 
0178   if (verbosity_ >= 0) {
0179     edm::LogVerbatim("IsoTrack") << "Parameters read from config file \n"
0180                                  << "\t minPt " << selectionParameters_.minPt << "\t theTrackQuality "
0181                                  << theTrackQuality_ << "\t minQuality " << selectionParameters_.minQuality
0182                                  << "\t maxDxyPV " << selectionParameters_.maxDxyPV << "\t maxDzPV "
0183                                  << selectionParameters_.maxDzPV << "\t maxChi2 " << selectionParameters_.maxChi2
0184                                  << "\t maxDpOverP " << selectionParameters_.maxDpOverP << "\t minOuterHit "
0185                                  << selectionParameters_.minOuterHit << "\t minLayerCrossed "
0186                                  << selectionParameters_.minLayerCrossed << "\t maxInMiss "
0187                                  << selectionParameters_.maxInMiss << "\t maxOutMiss "
0188                                  << selectionParameters_.maxOutMiss << "\t a_coneR " << a_coneR_ << "\t a_charIsoR "
0189                                  << a_charIsoR_ << "\t a_mipR " << a_mipR_ << "\t isMC " << isMC_ << "\t isQCD "
0190                                  << isQCD_ << "\t isAOD " << isAOD_;
0191     edm::LogVerbatim("IsoTrack") << trigNames_.size() << " triggers to be studied:";
0192     for (unsigned int k = 0; k < trigNames_.size(); ++k)
0193       edm::LogVerbatim("IsoTrack") << "[" << k << "] " << trigNames_[k];
0194   }
0195 
0196   tok_geom_ = esConsumes<CaloGeometry, CaloGeometryRecord>();
0197   tok_magField_ = esConsumes<MagneticField, IdealMagneticFieldRecord>();
0198 }
0199 
0200 IsoTrackCalibration::~IsoTrackCalibration() {
0201   // do anything here that needs to be done at desctruction time
0202   // (e.g. close files, deallocate resources etc.)
0203 }
0204 
0205 void IsoTrackCalibration::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0206   t_Run = iEvent.id().run();
0207   t_Event = iEvent.id().event();
0208   if (verbosity_ % 10 > 0)
0209     edm::LogVerbatim("IsoTrack") << "Run " << t_Run << " Event " << t_Event << " Luminosity "
0210                                  << iEvent.luminosityBlock() << " Bunch " << iEvent.bunchCrossing()
0211                                  << " starts ==========";
0212   //Get magnetic field
0213   const MagneticField *bField = &iSetup.getData(tok_magField_);
0214 
0215   // get handles to calogeometry
0216   const CaloGeometry *geo = &iSetup.getData(tok_geom_);
0217 
0218   //Get track collection
0219   edm::Handle<reco::TrackCollection> trkCollection;
0220   iEvent.getByToken(tok_genTrack_, trkCollection);
0221 
0222   //event weight for FLAT sample and PU information
0223   t_EventWeight = 1.0;
0224   edm::Handle<GenEventInfoProduct> genEventInfo;
0225   iEvent.getByToken(tok_ew_, genEventInfo);
0226   if (genEventInfo.isValid())
0227     t_EventWeight = genEventInfo->weight();
0228 
0229   //Define the best vertex and the beamspot
0230   edm::Handle<reco::VertexCollection> recVtxs;
0231   iEvent.getByToken(tok_recVtx_, recVtxs);
0232   edm::Handle<reco::BeamSpot> beamSpotH;
0233   iEvent.getByToken(tok_bs_, beamSpotH);
0234   math::XYZPoint leadPV(0, 0, 0);
0235 
0236   t_nVtx = recVtxs->size();
0237   h_nVtx->Fill(t_nVtx);
0238 
0239   if (!recVtxs->empty() && !((*recVtxs)[0].isFake())) {
0240     leadPV = math::XYZPoint((*recVtxs)[0].x(), (*recVtxs)[0].y(), (*recVtxs)[0].z());
0241   } else if (beamSpotH.isValid()) {
0242     leadPV = beamSpotH->position();
0243   }
0244   if ((verbosity_ / 100) % 10 > 0) {
0245     edm::LogVerbatim("IsoTrack") << "Primary Vertex " << leadPV;
0246     if (beamSpotH.isValid())
0247       edm::LogVerbatim("IsoTrack") << " Beam Spot " << beamSpotH->position();
0248   }
0249 
0250   // RecHits
0251   edm::Handle<EcalRecHitCollection> barrelRecHitsHandle;
0252   edm::Handle<EcalRecHitCollection> endcapRecHitsHandle;
0253   iEvent.getByToken(tok_EB_, barrelRecHitsHandle);
0254   iEvent.getByToken(tok_EE_, endcapRecHitsHandle);
0255   edm::Handle<HBHERecHitCollection> hbhe;
0256   iEvent.getByToken(tok_hbhe_, hbhe);
0257 
0258   //Trigger
0259   bool triggerOK = false;
0260   if (isMC_ && !isQCD_) {
0261     triggerOK = true;  // ignore HLT for single pion MC
0262   } else {
0263     trigger::TriggerEvent triggerEvent;
0264     edm::Handle<trigger::TriggerEvent> triggerEventHandle;
0265     iEvent.getByToken(tok_trigEvt_, triggerEventHandle);
0266     if (!triggerEventHandle.isValid()) {
0267       edm::LogVerbatim("IsoTrack") << "Error! Can't get the product " << triggerEvent_.label();
0268     } else {
0269       /////////////////////////////TriggerResults
0270       edm::Handle<edm::TriggerResults> triggerResults;
0271       iEvent.getByToken(tok_trigRes_, triggerResults);
0272       if (triggerResults.isValid()) {
0273         const edm::TriggerNames &triggerNames = iEvent.triggerNames(*triggerResults);
0274         const std::vector<std::string> &triggerNames_ = triggerNames.triggerNames();
0275         for (unsigned int iHLT = 0; iHLT < triggerResults->size(); iHLT++) {
0276           int hlt = triggerResults->accept(iHLT);
0277           if (hlt > 0) {
0278             for (unsigned int i = 0; i < trigNames_.size(); ++i) {
0279               if (triggerNames_[iHLT].find(trigNames_[i]) != std::string::npos) {
0280                 triggerOK = true;
0281                 if (verbosity_ % 10 > 0)
0282                   edm::LogVerbatim("IsoTrack")
0283                       << "This is the trigger we are looking for " << triggerNames_[iHLT] << " Flag " << hlt;
0284               }
0285             }
0286           }
0287         }
0288       }
0289     }
0290   }
0291 
0292   if (triggerOK) {
0293     //Propagate tracks to calorimeter surface)
0294     std::vector<spr::propagatedTrackDirection> trkCaloDirections;
0295     spr::propagateCALO(trkCollection, geo, bField, theTrackQuality_, trkCaloDirections, ((verbosity_ / 100) % 10 > 2));
0296     //Loop over tracks
0297     std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
0298     unsigned int nTracks(0);
0299     t_nTrk = trkCaloDirections.size();
0300 
0301     for (trkDetItr = trkCaloDirections.begin(), nTracks = 0; trkDetItr != trkCaloDirections.end();
0302          trkDetItr++, nTracks++) {
0303       const reco::Track *pTrack = &(*(trkDetItr->trkItr));
0304       if (verbosity_ % 10 > 0)
0305         edm::LogVerbatim("IsoTrack") << "This track : " << nTracks << " (pt/eta/phi/p) :" << pTrack->pt() << "/"
0306                                      << pTrack->eta() << "/" << pTrack->phi() << "/" << pTrack->p();
0307 
0308       t_ieta = 0;
0309       if (trkDetItr->okHCAL) {
0310         HcalDetId detId = (HcalDetId)(trkDetItr->detIdHCAL);
0311         t_ieta = detId.ieta();
0312       }
0313       // ---------- eta-dependent restriction on Pt ----------------------------
0314       selectionParameters_.minPt = constTrackPt_ - std::abs((double)t_ieta) * slopeTrackPt_;
0315       // -----------------------------------------------------------------------
0316 
0317       //Selection of good track
0318       t_selectTk = spr::goodTrack(pTrack, leadPV, selectionParameters_, ((verbosity_ / 100) % 10 > 2));
0319       spr::trackSelectionParameters oneCutParameters = selectionParameters_;
0320       oneCutParameters.maxDxyPV = 10;
0321       oneCutParameters.maxDzPV = 100;
0322       oneCutParameters.maxInMiss = 2;
0323       oneCutParameters.maxOutMiss = 2;
0324       bool qltyFlag = spr::goodTrack(pTrack, leadPV, oneCutParameters, ((verbosity_ / 100) % 10 > 2));
0325       oneCutParameters = selectionParameters_;
0326       oneCutParameters.maxDxyPV = 10;
0327       oneCutParameters.maxDzPV = 100;
0328       t_qltyMissFlag = spr::goodTrack(pTrack, leadPV, oneCutParameters, ((verbosity_ / 100) % 10 > 2));
0329       oneCutParameters = selectionParameters_;
0330       oneCutParameters.maxInMiss = 2;
0331       oneCutParameters.maxOutMiss = 2;
0332       t_qltyPVFlag = spr::goodTrack(pTrack, leadPV, oneCutParameters, ((verbosity_ / 100) % 10 > 2));
0333 
0334       if (verbosity_ % 10 > 0)
0335         edm::LogVerbatim("IsoTrack") << "qltyFlag|okECAL|okHCAL : " << qltyFlag << "|" << trkDetItr->okECAL << "/"
0336                                      << trkDetItr->okHCAL;
0337       if (qltyFlag && trkDetItr->okECAL && trkDetItr->okHCAL) {
0338         int nRH_eMipDR(0), nNearTRKs(0);
0339         //------ ecal energy around track -------------------------------
0340         t_eMipDR = spr::eCone_ecal(geo,
0341                                    barrelRecHitsHandle,
0342                                    endcapRecHitsHandle,
0343                                    trkDetItr->pointHCAL,
0344                                    trkDetItr->pointECAL,
0345                                    a_mipR_,
0346                                    trkDetItr->directionECAL,
0347                                    nRH_eMipDR);
0348         //---- isolation criteria ----------------------------------------------
0349         t_hmaxNearP =
0350             spr::chargeIsolationCone(nTracks, trkCaloDirections, a_charIsoR_, nNearTRKs, ((verbosity_ / 100) % 10 > 2));
0351 
0352         if (t_eMipDR < maxEcalEnr_ && t_hmaxNearP < maxNeighborTrackEnr_) {
0353           //------------ HCAL --------------------------------------------------
0354           //------ initialize arrays of DetID and hit energies -----------------
0355           t_DetIds->clear();
0356           t_DetIds1->clear();
0357           t_DetIds3->clear();
0358           t_HitEnergies->clear();
0359           t_HitEnergies1->clear();
0360           t_HitEnergies3->clear();
0361           int nRecHits(-999), nRecHits1(-999), nRecHits3(-999);
0362           std::vector<DetId> ids, ids1, ids3;
0363           //------ hcal energy in the main cone -------------------------------
0364           t_eHcal = spr::eCone_hcal(geo,
0365                                     hbhe,
0366                                     trkDetItr->pointHCAL,
0367                                     trkDetItr->pointECAL,
0368                                     a_coneR_,
0369                                     trkDetItr->directionHCAL,
0370                                     nRecHits,
0371                                     ids,
0372                                     *t_HitEnergies);
0373           t_DetIds->reserve(ids.size());
0374           for (unsigned int k = 0; k < ids.size(); ++k) {
0375             t_DetIds->push_back(ids[k].rawId());
0376           }
0377           //----- hcal energy in the extended cone 1 (a_coneR+10) --------------
0378           t_eHcal10 = spr::eCone_hcal(geo,
0379                                       hbhe,
0380                                       trkDetItr->pointHCAL,
0381                                       trkDetItr->pointECAL,
0382                                       a_coneR_ + 10,
0383                                       trkDetItr->directionHCAL,
0384                                       nRecHits1,
0385                                       ids1,
0386                                       *t_HitEnergies1);
0387           t_DetIds1->reserve(ids1.size());
0388           for (unsigned int k = 0; k < ids1.size(); ++k) {
0389             t_DetIds1->push_back(ids1[k].rawId());
0390           }
0391           //----- hcal energy in the extended cone 3 (a_coneR+30) --------------
0392           t_eHcal30 = spr::eCone_hcal(geo,
0393                                       hbhe,
0394                                       trkDetItr->pointHCAL,
0395                                       trkDetItr->pointECAL,
0396                                       a_coneR_ + 30,
0397                                       trkDetItr->directionHCAL,
0398                                       nRecHits3,
0399                                       ids3,
0400                                       *t_HitEnergies3);
0401           t_DetIds3->reserve(ids3.size());
0402           for (unsigned int k = 0; k < ids3.size(); ++k) {
0403             t_DetIds3->push_back(ids3[k].rawId());
0404           }
0405 
0406           t_p = pTrack->p();
0407           t_pt = pTrack->pt();
0408           t_phi = pTrack->phi();
0409           if (verbosity_ % 10 > 0) {
0410             edm::LogVerbatim("IsoTrack") << "This track : " << nTracks << " (pt/eta/phi/p) :" << pTrack->pt() << "/"
0411                                          << pTrack->eta() << "/" << pTrack->phi() << "/" << t_p;
0412             edm::LogVerbatim("IsoTrack") << "e_MIP " << t_eMipDR << " Chg Isolation " << t_hmaxNearP << " eHcal"
0413                                          << t_eHcal << " ieta " << t_ieta << " Quality " << t_qltyMissFlag << ":"
0414                                          << t_qltyPVFlag << ":" << t_selectTk;
0415             for (unsigned int lll = 0; lll < t_DetIds->size(); lll++) {
0416               edm::LogVerbatim("IsoTrack")
0417                   << "det id is = " << t_DetIds->at(lll) << "   hit enery is  = " << t_HitEnergies->at(lll);
0418             }
0419           }
0420           tree->Fill();
0421         }  // end of conditions on t_eMipDR and t_hmaxNearP
0422       }    // end of loose check of track quality
0423     }      // end of loop over tracks
0424 
0425     h_nTrk->Fill(nTracks);
0426   }  // end of triggerOK
0427 }
0428 
0429 void IsoTrackCalibration::beginJob() {
0430   h_nVtx = fs_->make<TH1F>("h_nVtx", "h_nVtx", 100, 0, 100);
0431   h_nTrk = fs_->make<TH1F>("h_nTrk", "h_nTrk", 100, 0, 2000);
0432 
0433   tree = fs_->make<TTree>("CalibTree", "CalibTree");
0434 
0435   tree->Branch("t_Run", &t_Run, "t_Run/I");
0436   tree->Branch("t_Event", &t_Event, "t_Event/I");
0437   tree->Branch("t_nVtx", &t_nVtx, "t_nVtx/I");
0438   tree->Branch("t_nTrk", &t_nTrk, "t_nTrk/I");
0439   tree->Branch("t_EventWeight", &t_EventWeight, "t_EventWeight/D");
0440   tree->Branch("t_p", &t_p, "t_p/D");
0441   tree->Branch("t_pt", &t_pt, "t_pt/D");
0442   tree->Branch("t_ieta", &t_ieta, "t_ieta/I");
0443   tree->Branch("t_phi", &t_phi, "t_phi/D");
0444   tree->Branch("t_eMipDR", &t_eMipDR, "t_eMipDR/D");
0445   tree->Branch("t_eHcal", &t_eHcal, "t_eHcal/D");
0446   tree->Branch("t_eHcal10", &t_eHcal10, "t_eHcal10/D");
0447   tree->Branch("t_eHcal30", &t_eHcal30, "t_eHcal30/D");
0448   tree->Branch("t_hmaxNearP", &t_hmaxNearP, "t_hmaxNearP/D");
0449   tree->Branch("t_selectTk", &t_selectTk, "t_selectTk/O");
0450   tree->Branch("t_qltyMissFlag", &t_qltyMissFlag, "t_qltyMissFlag/O");
0451   tree->Branch("t_qltyPVFlag", &t_qltyPVFlag, "t_qltyPVFlag/O)");
0452 
0453   t_DetIds = new std::vector<unsigned int>();
0454   t_DetIds1 = new std::vector<unsigned int>();
0455   t_DetIds3 = new std::vector<unsigned int>();
0456   t_HitEnergies = new std::vector<double>();
0457   t_HitEnergies1 = new std::vector<double>();
0458   t_HitEnergies3 = new std::vector<double>();
0459 
0460   tree->Branch("t_DetIds", "std::vector<unsigned int>", &t_DetIds);
0461   //tree->Branch("t_DetIds1",      "std::vector<unsigned int>", &t_DetIds1);
0462   //tree->Branch("t_DetIds3",      "std::vector<unsigned int>", &t_DetIds3);
0463   tree->Branch("t_HitEnergies", "std::vector<double>", &t_HitEnergies);
0464   //tree->Branch("t_HitEnergies1", "std::vector<double>",       &t_HitEnergies1);
0465   //tree->Branch("t_HitEnergies3", "std::vector<double>",       &t_HitEnergies3);
0466 }
0467 
0468 // ------------ method called when starting to processes a run  ------------
0469 void IsoTrackCalibration::beginRun(edm::Run const &iRun, edm::EventSetup const &iSetup) {
0470   bool changed_(true);
0471   bool flag = hltConfig_.init(iRun, iSetup, "HLT", changed_);
0472   edm::LogInfo("HcalIsoTrack") << "Run[" << nRun_ << "] " << iRun.run() << " process HLT init flag " << flag
0473                                << " change flag " << changed_;
0474 
0475   // check if trigger names in (new) config
0476   if (changed_) {
0477     edm::LogInfo("HcalIsoTrack") << "New trigger menu found !!!";
0478     const unsigned int n(hltConfig_.size());
0479     for (unsigned itrig = 0; itrig < trigNames_.size(); itrig++) {
0480       unsigned int triggerindx = hltConfig_.triggerIndex(trigNames_[itrig]);
0481       if (triggerindx >= n) {
0482         edm::LogWarning("HcalIsoTrack") << trigNames_[itrig] << " " << triggerindx << " does not exist in "
0483                                         << "the current menu";
0484       } else {
0485         edm::LogInfo("HcalIsoTrack") << trigNames_[itrig] << " " << triggerindx << " exists";
0486       }
0487     }
0488   }
0489 }
0490 
0491 // ------------ method called when ending the processing of a run  ------------
0492 void IsoTrackCalibration::endRun(edm::Run const &iRun, edm::EventSetup const &) {
0493   ++nRun_;
0494   edm::LogWarning("HcalIsoTrack") << "endRun[" << nRun_ << "] " << iRun.run();
0495 }
0496 
0497 // ------------ method called when starting to processes a luminosity block  ------------
0498 void IsoTrackCalibration::beginLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &) {}
0499 // ------------ method called when ending the processing of a luminosity block  ------------
0500 void IsoTrackCalibration::endLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &) {}
0501 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0502 void IsoTrackCalibration::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0503   //The following says we do not know what parameters are allowed so do no validation
0504   // Please change this to state exactly what you do use, even if it is no parameters
0505   edm::ParameterSetDescription desc;
0506   desc.setUnknown();
0507   descriptions.addDefault(desc);
0508 }
0509 
0510 //define this as a plug-in
0511 DEFINE_FWK_MODULE(IsoTrackCalibration);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!