/
githubmirror
/
cmssw
Обзор
Документация
Войти
/
githubmirror
/
cmssw
Код
Запросы
0
Пакеты
0
Релизы
0
Аналитика
Безопасность
master
Alignment/OfflineValidation/plugins/EopElecTreeWriter.cc
915 строк
32 KB
Shahzad Malik Muzaffar
[ALCA] Changes suggested by new llvm18 clang-format
05 сен 2024, 17:03
05 сен 2024, 17:03
dfd2283
Код
Авторство
О чём код?
// -*- C++ -*- // // Package: EopTreeWriter // Class: EopTreeWriter // /**\class EopTreeWriter EopTreeWriter.cc Alignment/OfflineValidation/plugins/EopTreeWriter.cc // Description: <one line class summary> Implementation: <Notes on implementation> */ // // Original Author: Holger Enderle // Created: Thu Dec 4 11:22:48 CET 2008 // $Id: EopElecTreeWriter.cc,v 1.8 2012/08/30 08:40:51 cgoetzma Exp $ // // // framework include files #include "Alignment/OfflineValidation/interface/EopElecVariables.h" #include "CommonTools/UtilAlgos/interface/TFileService.h" #include "CommonTools/Utils/interface/TFileDirectory.h" #include "DataFormats/CaloTowers/interface/CaloTower.h" #include "DataFormats/Common/interface/Handle.h" #include "DataFormats/Common/interface/TriggerResults.h" #include "DataFormats/DTRecHit/interface/DTRecHitCollection.h" #include "DataFormats/EcalDetId/interface/EBDetId.h" #include "DataFormats/EcalDetId/interface/EEDetId.h" #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h" #include "DataFormats/EgammaCandidates/interface/GsfElectron.h" #include "DataFormats/EgammaCandidates/interface/GsfElectronCore.h" #include "DataFormats/EgammaReco/interface/BasicCluster.h" #include "DataFormats/EgammaReco/interface/ElectronSeed.h" #include "DataFormats/EgammaReco/interface/ElectronSeedFwd.h" #include "DataFormats/EgammaReco/interface/SuperCluster.h" #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h" #include "DataFormats/GsfTrackReco/interface/GsfTrack.h" #include "DataFormats/HLTReco/interface/TriggerEvent.h" #include "DataFormats/HcalIsolatedTrack/interface/IsolatedPixelTrackCandidate.h" #include "DataFormats/HcalIsolatedTrack/interface/IsolatedPixelTrackCandidateFwd.h" #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h" #include "DataFormats/JetReco/interface/CaloJet.h" #include "DataFormats/Math/interface/deltaR.h" #include "DataFormats/ParticleFlowReco/interface/GsfPFRecTrack.h" #include "DataFormats/RecoCandidate/interface/RecoCaloTowerCandidate.h" #include "DataFormats/TrackReco/interface/Track.h" #include "DataFormats/TrackReco/interface/TrackFwd.h" #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h" #include "DataFormats/VertexReco/interface/Vertex.h" #include "DataFormats/VertexReco/interface/VertexFwd.h" #include "FWCore/Framework/interface/ESHandle.h" #include "FWCore/Framework/interface/Event.h" #include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/Framework/interface/Frameworkfwd.h" #include "FWCore/Framework/interface/MakerMacros.h" #include "FWCore/Framework/interface/one/EDAnalyzer.h" #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/ServiceRegistry/interface/Service.h" #include "FWCore/Utilities/interface/EDGetToken.h" #include "Geometry/CaloGeometry/interface/CaloGeometry.h" #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h" #include "Geometry/Records/interface/CaloGeometryRecord.h" #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h" #include "MagneticField/Engine/interface/MagneticField.h" #include "RecoEcal/EgammaCoreTools/interface/SuperClusterShapeAlgo.h" #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h" #include "TrackPropagation/SteppingHelixPropagator/interface/SteppingHelixPropagator.h" #include "TrackingTools/GeomPropagators/interface/PropagationExceptions.h" #include "TrackingTools/GsfTools/interface/MultiTrajectoryStateMode.h" #include "TrackingTools/GsfTools/interface/MultiTrajectoryStateTransform.h" #include "TrackingTools/TrackAssociator/interface/TrackAssociatorParameters.h" #include "TrackingTools/TrackAssociator/interface/TrackDetectorAssociator.h" #include "TCanvas.h" #include "TH1.h" #include "TH2D.h" #include "TLorentzVector.h" #include "TMath.h" #include "TTree.h" struct EopTriggerType { bool fired; double prescale; int index; EopTriggerType() { fired = false; prescale = 0; index = 0; } }; // // class decleration // using namespace std; class EopElecTreeWriter : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> { public: explicit EopElecTreeWriter(const edm::ParameterSet&); ~EopElecTreeWriter() override; static void fillDescriptions(edm::ConfigurationDescriptions&); private: // methods void beginRun(const edm::Run&, const edm::EventSetup&) override; void endRun(const edm::Run&, const edm::EventSetup&) override {} void analyze(const edm::Event&, const edm::EventSetup&) override; void endJob() override; // ----------member data --------------------------- // ES tokens const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magFieldToken_; const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkGeomToken_; const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeomToken_; // EDM tokens const edm::EDGetTokenT<reco::VertexCollection> theVertexCollectionToken_; const edm::EDGetTokenT<HBHERecHitCollection> theHBHERecHitCollectionToken_; const edm::EDGetTokenT<EcalRecHitCollection> theEcalRecHitCollectionToken_; const edm::EDGetTokenT<reco::SuperClusterCollection> theBarrelSupClusCollectionToken_; const edm::EDGetTokenT<reco::SuperClusterCollection> theEndCapSupClusCollectionToken_; const edm::EDGetTokenT<edm::TriggerResults> theTriggerResultsToken_; const edm::EDGetTokenT<trigger::TriggerEvent> theTriggerEventToken_; const edm::EDGetTokenT<reco::GsfTrackCollection> theGsfTrackCollectionToken_; const edm::EDGetTokenT<reco::GsfElectronCoreCollection> theGsfElectronCoreCollectionToken_; // some static constants static constexpr float k_etaBarrel = 1.55; static constexpr float k_etaEndcap = 1.44; // other member data std::string theTrigger_; std::string theFilter_; bool debugTriggerSelection_; edm::Service<TFileService> fs_; TTree* tree_; EopElecVariables* treeMemPtr_; HLTConfigProvider hltConfig_; std::vector<std::string> triggerNames_; MultiTrajectoryStateTransform* mtsTransform_; // histograms // Cut flow (events number) TH1D* h_nEvents; TH1D* h_nEventsWithVertex; TH1D* h_nEventsTriggered; TH1D* h_nEventsHLTFilter; TH1D* h_nEventsHLTelectron; TH1D* h_nEventsHLTrejected; TH1D* h_nEvents2Elec; TH1D* h_nHLTelectrons; TH1D* h_nTrkRejectedPerEvt; TH1D* h_nTrkSelectedPerEvt; // Cut flow (tracks number) TH1D* h_nTracks; TH1D* h_nTracksFiltered; TH1D* h_cut_Ptmin; TH1D* h_cut_OneSCmatch; TH1D* h_counter1; TH1D* h_counter2; TH1D* h_distToClosestSCgsf; TH1D* h_distToClosestSC; TH1D* h_EcalEnergy; TH1D* h_Momentum; TH1D* HcalEnergy; TH1D* h_fBREM; TH1D* h_Eop_InnerNegative; TH1D* h_Eop_InnerPositive; TH2D* HcalVSEcal; }; namespace eopUtils { static constexpr float R_ECAL = 136.5; static constexpr float Z_Endcap = 328.0; static constexpr float etaBarrelEndcap = 1.479; // Function to convert the eta of the track to a detector eta (for matching with SC) float ecalEta(float EtaParticle, float Zvertex, float RhoVertex) { if (EtaParticle != 0.) { float Theta = 0.0; float ZEcal = (R_ECAL - RhoVertex) * sinh(EtaParticle) + Zvertex; if (ZEcal != 0.0) Theta = atan(R_ECAL / ZEcal); if (Theta < 0.0) Theta = Theta + M_PI; float ETA = -log(tan(0.5 * Theta)); if (fabs(ETA) > etaBarrelEndcap) { float Zend = Z_Endcap; if (EtaParticle < 0.0) Zend = -Zend; float Zlen = Zend - Zvertex; float RR = Zlen / sinh(EtaParticle); Theta = atan((RR + RhoVertex) / Zend); if (Theta < 0.0) Theta = Theta + M_PI; ETA = -log(tan(0.5 * Theta)); } return ETA; } else { edm::LogWarning("") << "[EcalPositionFromTrack::etaTransformation] Warning: Eta equals to zero, not correcting"; return EtaParticle; } } } // namespace eopUtils // constructors and destructor EopElecTreeWriter::EopElecTreeWriter(const edm::ParameterSet& iConfig) : magFieldToken_(esConsumes()), tkGeomToken_(esConsumes()), caloGeomToken_(esConsumes()), theVertexCollectionToken_(consumes<reco::VertexCollection>(edm::InputTag("offlinePrimaryVertices"))), theHBHERecHitCollectionToken_(consumes<HBHERecHitCollection>(edm::InputTag("hbhereco", ""))), theEcalRecHitCollectionToken_(consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit", "EcalRecHitsEB"))), theBarrelSupClusCollectionToken_( consumes<reco::SuperClusterCollection>(edm::InputTag("hybridSuperClusters", ""))), theEndCapSupClusCollectionToken_(consumes<reco::SuperClusterCollection>( edm::InputTag("multi5x5SuperClusters", "multi5x5EndcapSuperClusters"))), theTriggerResultsToken_(consumes<edm::TriggerResults>(edm::InputTag("TriggerResults", "", "HLT"))), theTriggerEventToken_(consumes<trigger::TriggerEvent>(edm::InputTag("hltTriggerSummaryAOD"))), theGsfTrackCollectionToken_(consumes<reco::GsfTrackCollection>(iConfig.getParameter<edm::InputTag>("src"))), theGsfElectronCoreCollectionToken_( consumes<reco::GsfElectronCoreCollection>(edm::InputTag("gedGsfElectronCores"))), theTrigger_(iConfig.getParameter<std::string>("triggerPath")), theFilter_(iConfig.getParameter<std::string>("hltFilter")), debugTriggerSelection_(iConfig.getParameter<bool>("debugTriggerSelection")) { usesResource(TFileService::kSharedResource); // TTree creation tree_ = fs_->make<TTree>("EopTree", "EopTree"); treeMemPtr_ = new EopElecVariables; tree_->Branch("EopElecVariables", &treeMemPtr_); // address of pointer! // Control histograms declaration h_distToClosestSC = fs_->make<TH1D>("distToClosestSC", "distToClosestSC", 100, 0, 0.1); h_distToClosestSCgsf = fs_->make<TH1D>("distToClosestSCgsf", "distToClosestSCgsf", 100, 0, 0.1); h_EcalEnergy = fs_->make<TH1D>("EcalEnergy", "EcalEnergy", 100, 0, 200); h_Momentum = fs_->make<TH1D>("Momentum", "Momentum", 100, 0, 200); HcalEnergy = fs_->make<TH1D>("HcalEnergy", "HcalEnergy", 100, 0, 40); h_fBREM = fs_->make<TH1D>("fBREM", "fBREM", 100, -0.2, 1); h_Eop_InnerNegative = fs_->make<TH1D>("Eop_InnerNegative", "Eop_InnerNegative", 100, 0, 3); h_Eop_InnerPositive = fs_->make<TH1D>("Eop_InnerPositive", "Eop_InnerPositive", 100, 0, 3); HcalVSEcal = fs_->make<TH2D>("HcalVSEcal", "HcalVSEcal", 100, 0, 160, 100, 0, 10); h_nEvents = fs_->make<TH1D>("nEvents", "nEvents", 1, 0, 1); h_nEventsWithVertex = fs_->make<TH1D>("nEventsWithVertex", "nEventsWithVertex", 1, 0, 1); h_nEventsTriggered = fs_->make<TH1D>("nEventsTriggered", "nEventsTriggered", 1, 0, 1); h_nEventsHLTFilter = fs_->make<TH1D>("nEventsHLTFilter", "nEventsHLTFilter", 1, 0, 1); h_nEventsHLTelectron = fs_->make<TH1D>("nEventsHLTelectron", "nEventsHLTelectron", 1, 0, 1); h_nEventsHLTrejected = fs_->make<TH1D>("nEventsHLTrejected", "nEventsHLTrejected", 1, 0, 1); h_nEvents2Elec = fs_->make<TH1D>("nEvents2Elec", "nEvents2Elec", 1, 0, 1); h_nHLTelectrons = fs_->make<TH1D>("nHLTelectrons", "nHLTelectrons", 20, 0, 20); h_nTrkRejectedPerEvt = fs_->make<TH1D>("nTrkRejectedPerEvt", "nTrkRejectedPerEvt", 20, 0, 20); h_nTrkSelectedPerEvt = fs_->make<TH1D>("nTrkSelectedPerEvt", "nTrkSelectedPerEvt", 20, 0, 20); h_nTracks = fs_->make<TH1D>("nTracks", "nTracks", 1, 0, 1); h_nTracksFiltered = fs_->make<TH1D>("nTracksFiltered", "nTracksFiltered", 1, 0, 1); h_cut_Ptmin = fs_->make<TH1D>("cut_Ptmin", "cut_Ptmin", 1, 0, 1); h_cut_OneSCmatch = fs_->make<TH1D>("cut_OneSCmatch", "cut_OneSCmatch", 1, 0, 1); h_counter1 = fs_->make<TH1D>("counter1", "counter1", 1, 0, 1); h_counter2 = fs_->make<TH1D>("counter2", "counter2", 1, 0, 1); } EopElecTreeWriter::~EopElecTreeWriter() { // control histograms h_distToClosestSC->SetXTitle("distance from track to closest SuperCluster in eta-phi plan (weighted matching)"); h_distToClosestSC->SetYTitle("# Tracks"); h_distToClosestSCgsf->SetXTitle("distance from track to closest SuperCluster in eta-phi plan (gsfElectronCore)"); h_distToClosestSCgsf->SetYTitle("# Tracks"); h_EcalEnergy->SetXTitle("Ecal energy deposit (GeV)"); h_EcalEnergy->SetYTitle("# tracks"); HcalEnergy->SetXTitle("Hcal energy deposit (GeV)"); HcalEnergy->SetYTitle("# tracks"); h_Momentum->SetXTitle("Momentum magnitude (GeV)"); h_Momentum->SetYTitle("# tracks"); h_Eop_InnerNegative->SetXTitle("E/p"); h_Eop_InnerNegative->SetYTitle("# tracks"); h_Eop_InnerPositive->SetXTitle("E/p"); h_Eop_InnerPositive->SetYTitle("# tracks"); HcalVSEcal->SetXTitle("Ecal energy (GeV)"); HcalVSEcal->SetYTitle("Hcal energy (GeV)"); } //########################################### //# method called to for each event # //########################################### void EopElecTreeWriter::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) { using namespace edm; h_nEvents->Fill(0.5); Double_t EnergyHcalIn01; Double_t EnergyHcalIn02; Double_t EnergyHcalIn03; Double_t EnergyHcalIn04; Double_t EnergyHcalIn05; Double_t etaWidth; Double_t phiWidth; Int_t algo_ID; Double_t EnergyEcal; Double_t fbrem; Double_t pin; Double_t etaIn; Double_t phiIn; Double_t pout; Double_t etaOut; Double_t phiOut; Double_t MaxPtIn01; Double_t SumPtIn01; Bool_t NoTrackIn0015; Double_t MaxPtIn02; Double_t SumPtIn02; Bool_t NoTrackIn0020; Double_t MaxPtIn03; Double_t SumPtIn03; Bool_t NoTrackIn0025; Double_t MaxPtIn04; Double_t SumPtIn04; Bool_t NoTrackIn0030; Double_t MaxPtIn05; Double_t SumPtIn05; Bool_t NoTrackIn0035; Bool_t NoTrackIn0040; Double_t dRSC_first; Double_t dRSC_second; Double_t etaSC; Double_t phiSC; Int_t nbSC; //to count the nb of SuperCluster matching with a given track Int_t nBasicClus; //to count the nb of basic cluster in a given superCluster Bool_t isEcalDriven; Bool_t isTrackerDriven; Bool_t isBarrel; Bool_t isEndcap; Bool_t TrigTag; //--------------- for GsfTrack propagation through tracker --------------- const MagneticField* magField_ = &iSetup.getData(magFieldToken_); const TrackerGeometry* trackerGeom_ = &iSetup.getData(tkGeomToken_); MultiTrajectoryStateTransform mtsTransform(trackerGeom_, magField_); //--------------- Super Cluster ----------------- // getting primary vertex (necessary to convert eta track to eta detector const reco::VertexCollection& vertex = iEvent.get(theVertexCollectionToken_); if (vertex.empty()) { edm::LogError("EopElecTreeWriter") << "Error: no primary vertex found!"; return; } const reco::Vertex& vert = vertex.front(); h_nEventsWithVertex->Fill(0.5); // getting calorimeter geometry const CaloGeometry* geo = &iSetup.getData(caloGeomToken_); const CaloSubdetectorGeometry* subGeo = geo->getSubdetectorGeometry(DetId::Ecal, EcalBarrel); if (subGeo == nullptr) edm::LogError("EopElecTreeWriter") << "ERROR: unable to find SubDetector geometry!!!"; // getting Hcal rechits const HBHERecHitCollection* HcalHits = &iEvent.get(theHBHERecHitCollectionToken_); // getting Ecal rechits const EcalRecHitCollection* rhc = &iEvent.get(theEcalRecHitCollectionToken_); if (rhc == nullptr) edm::LogError("EopElecTreeWriter") << "ERROR: unable to find the EcalRecHit collection !!!"; // getting SuperCluster const reco::SuperClusterCollection* BarrelSupClusCollection = &iEvent.get(theBarrelSupClusCollectionToken_); const reco::SuperClusterCollection* EndcapSupClusCollection = &iEvent.get(theEndCapSupClusCollectionToken_); // necessary to re-calculate phi and eta width of SuperClusters SuperClusterShapeAlgo SCShape(rhc, subGeo); //--------------- Trigger ----------------- TrigTag = false; const edm::TriggerResults* trigRes = &iEvent.get(theTriggerResultsToken_); // trigger event const trigger::TriggerEvent* triggerEvent = &iEvent.get(theTriggerEventToken_); // our trigger table std::map<std::string, EopTriggerType> HLTpaths; for (const auto& triggerName : triggerNames_) { if (triggerName.find(theTrigger_) != 0) continue; EopTriggerType myTrigger; const unsigned int prescaleSize = hltConfig_.prescaleSize(); for (unsigned int ps = 0; ps < prescaleSize; ps++) { auto const prescaleValue = hltConfig_.prescaleValue<double>(ps, triggerName); if (prescaleValue != 1) { myTrigger.prescale = prescaleValue; } } myTrigger.index = hltConfig_.triggerIndex(triggerName); if (myTrigger.index == -1) continue; myTrigger.fired = trigRes->wasrun(myTrigger.index) && trigRes->accept(myTrigger.index) && !trigRes->error(myTrigger.index); HLTpaths[triggerName] = myTrigger; } // First cut : trigger cut std::string firstFiredPath = ""; for (const auto& it : HLTpaths) { if (it.second.fired) { TrigTag = true; firstFiredPath = it.first; break; } } if (!TrigTag) return; h_nEventsTriggered->Fill(0.5); // Displaying filters label from the first fired trigger // Useful for finding the good filter label std::vector<std::string> filters = hltConfig_.moduleLabels(firstFiredPath); if (debugTriggerSelection_) { edm::LogInfo("EopElecTreeWriter") << "filters : "; for (unsigned int i = 0; i < filters.size(); i++) { edm::LogInfo("EopElecTreeWriter") << filters[i] << " "; } } // Getting HLT electrons edm::InputTag testTag(theFilter_, "", "HLT"); int testindex = triggerEvent->filterIndex(testTag); if (testindex >= triggerEvent->sizeFilters()) return; h_nEventsHLTFilter->Fill(0.5); const trigger::Keys& KEYS_el(triggerEvent->filterKeys(testindex)); std::vector<const trigger::TriggerObject*> HLTelectrons; for (unsigned int i = 0; i != KEYS_el.size(); ++i) { const trigger::TriggerObject* triggerObject_el = &(triggerEvent->getObjects().at(KEYS_el[i])); HLTelectrons.push_back(triggerObject_el); } h_nHLTelectrons->Fill(HLTelectrons.size()); if (HLTelectrons.empty()) return; h_nEventsHLTelectron->Fill(0.5); // finding the HLT electron with highest pt and saving the corresponding index unsigned int HighPtIndex = 0; double maxPtHLT = -5.; for (unsigned int j = 0; j < HLTelectrons.size(); j++) { if (HLTelectrons[j]->pt() > maxPtHLT) { maxPtHLT = HLTelectrons[j]->pt(); HighPtIndex = j; } } //----------------- Tracks ------------------- // getting GsfTrack const reco::GsfTrackCollection& tracks = iEvent.get(theGsfTrackCollectionToken_); // filtering track int nRejected = 0; int nSelected = 0; std::vector<const reco::GsfTrack*> filterTracks; for (const auto& track : tracks) { h_nTracks->Fill(0.5); double deltar = reco::deltaR(track.eta(), track.phi(), HLTelectrons[HighPtIndex]->eta(), HLTelectrons[HighPtIndex]->phi()); // remove the triggered electron with highest pt if (deltar < 0.025) { treeMemPtr_->px_rejected_track = track.px(); treeMemPtr_->py_rejected_track = track.py(); treeMemPtr_->pz_rejected_track = track.pz(); treeMemPtr_->p_rejected_track = track.p(); nRejected++; continue; } filterTracks.push_back(&track); // we use all the others nSelected++; h_nTracksFiltered->Fill(0.5); } h_nTrkRejectedPerEvt->Fill(nRejected); h_nTrkSelectedPerEvt->Fill(nSelected); if (nRejected == 0) return; h_nEventsHLTrejected->Fill(0.5); if (filterTracks.empty()) return; h_nEvents2Elec->Fill(0.5); //-------- test:Matching SC/track using gsfElectonCore collection -------- const reco::GsfElectronCoreCollection* electrons = &iEvent.get(theGsfElectronCoreCollectionToken_); for (const auto& elec : *electrons) { double etaGSF = eopUtils::ecalEta((elec.gsfTrack())->eta(), vert.z(), (vert.position()).rho()); if ((elec.gsfTrack())->pt() < 10.) continue; double DELTAR = 0; DELTAR = reco::deltaR((elec.superCluster())->eta(), (elec.superCluster())->phi(), etaGSF, (elec.gsfTrack())->phi()); if (DELTAR < 0.1) h_distToClosestSCgsf->Fill(DELTAR); } //------------------------------------------------------------ //-------------- Loop on tracks ------------------------- for (const auto& track : filterTracks) { // initializing variables isEcalDriven = false; isTrackerDriven = false; isBarrel = false; isEndcap = false; etaWidth = 0; phiWidth = 0; etaSC = 0; phiSC = 0; fbrem = 0; pin = 0; etaIn = 0; phiIn = 0; pout = 0; etaOut = 0; phiOut = 0; algo_ID = 0; EnergyEcal = 0; dRSC_first = 999; dRSC_second = 9999; nbSC = 0; nBasicClus = 0; // First cut on momentum magnitude h_Momentum->Fill(track->p()); if (track->pt() < 10.) continue; h_cut_Ptmin->Fill(0.5); // calculating track parameters at innermost and outermost for Gsf tracks TrajectoryStateOnSurface inTSOS = mtsTransform.innerStateOnSurface(*track); TrajectoryStateOnSurface outTSOS = mtsTransform.outerStateOnSurface(*track); if (inTSOS.isValid() && outTSOS.isValid()) { GlobalVector inMom; //multiTrajectoryStateMode::momentumFromModeCartesian(inTSOS, inMom_); multiTrajectoryStateMode::momentumFromModePPhiEta(inTSOS, inMom); pin = inMom.mag(); etaIn = inMom.eta(); phiIn = inMom.phi(); GlobalVector outMom; //multiTrajectoryStateMode::momentumFromModeCartesian(outTSOS, outMom_); multiTrajectoryStateMode::momentumFromModePPhiEta(outTSOS, outMom); pout = outMom.mag(); etaOut = outMom.eta(); phiOut = outMom.phi(); fbrem = (pin - pout) / pin; h_fBREM->Fill(fbrem); } // Matching track with Hcal rec hits EnergyHcalIn01 = 0; EnergyHcalIn02 = 0; EnergyHcalIn03 = 0; EnergyHcalIn04 = 0; EnergyHcalIn05 = 0; //for (std::vector<HBHERecHit>::const_iterator hcal = (*HcalHits).begin(); hcal != (*HcalHits).end(); hcal++) { for (const auto& hcal : *HcalHits) { GlobalPoint posH = geo->getPosition(hcal.detid()); double phihit = posH.phi(); double etahit = posH.eta(); double dR = reco::deltaR(etahit, phihit, etaOut, phiOut); // saving Hcal energy deposit measured for different eta-phi radius if (dR < 0.1) EnergyHcalIn01 += hcal.energy(); if (dR < 0.2) EnergyHcalIn02 += hcal.energy(); if (dR < 0.3) EnergyHcalIn03 += hcal.energy(); if (dR < 0.4) EnergyHcalIn04 += hcal.energy(); if (dR < 0.5) EnergyHcalIn05 += hcal.energy(); } HcalEnergy->Fill(EnergyHcalIn02); //Isolation against charged particles MaxPtIn01 = 0.; SumPtIn01 = 0.; NoTrackIn0015 = true; MaxPtIn02 = 0.; SumPtIn02 = 0.; NoTrackIn0020 = true; MaxPtIn03 = 0.; SumPtIn03 = 0.; NoTrackIn0025 = true; MaxPtIn04 = 0.; SumPtIn04 = 0.; NoTrackIn0030 = true; MaxPtIn05 = 0.; SumPtIn05 = 0.; NoTrackIn0035 = true; NoTrackIn0040 = true; for (const auto& track1 : filterTracks) { if (track == track1) continue; double etaIn1 = 0.; double phiIn1 = 0.; TrajectoryStateOnSurface inTSOS1 = mtsTransform.innerStateOnSurface(*track1); if (inTSOS1.isValid()) { GlobalVector inMom1; multiTrajectoryStateMode::momentumFromModePPhiEta(inTSOS1, inMom1); etaIn1 = inMom1.eta(); phiIn1 = inMom1.phi(); } if (etaIn1 == 0 && phiIn1 == 0) continue; double dR = reco::deltaR(etaIn1, phiIn1, etaIn, phiIn); // different radius of inner isolation cone if (dR < 0.015) NoTrackIn0015 = false; if (dR < 0.020) NoTrackIn0020 = false; if (dR < 0.025) NoTrackIn0025 = false; if (dR < 0.030) NoTrackIn0030 = false; if (dR < 0.035) NoTrackIn0035 = false; if (dR < 0.040) NoTrackIn0040 = false; //calculate maximum Pt and sum Pt inside cones of different radius if (dR < 0.1) { SumPtIn01 += track1->pt(); if (track1->pt() > MaxPtIn01) { MaxPtIn01 = track1->pt(); } } if (dR < 0.2) { SumPtIn02 += track1->pt(); if (track1->pt() > MaxPtIn02) { MaxPtIn02 = track1->pt(); } } if (dR < 0.3) { SumPtIn03 += track1->pt(); if (track1->pt() > MaxPtIn03) { MaxPtIn03 = track1->pt(); } } if (dR < 0.4) { SumPtIn04 += track1->pt(); if (track1->pt() > MaxPtIn04) { MaxPtIn04 = track1->pt(); } } if (dR < 0.5) { SumPtIn05 += track1->pt(); if (track1->pt() > MaxPtIn05) { MaxPtIn05 = track1->pt(); } } } // Track-SuperCluster matching double dRSC; double etaAtEcal = 0; // to convert eta from track to detector frame etaAtEcal = eopUtils::ecalEta(etaIn, vert.z(), (vert.position()).rho()); //Barrel if (std::abs(track->eta()) < k_etaBarrel) { for (const auto& SC : *BarrelSupClusCollection) { dRSC = reco::deltaR(SC.eta(), SC.phi(), etaAtEcal, phiIn); if (dRSC < dRSC_first) { dRSC_first = dRSC; } // distance in eta-phi plan to closest SC else if (dRSC < dRSC_second) { dRSC_second = dRSC; } // to next closest SC if (dRSC < 0.09) { //Calculate phiWidth & etaWidth for associated SuperClusters SCShape.Calculate_Covariances(SC); phiWidth = SCShape.phiWidth(); etaWidth = SCShape.etaWidth(); etaSC = SC.eta(); phiSC = SC.phi(); algo_ID = SC.algoID(); EnergyEcal = SC.energy(); nBasicClus = SC.clustersSize(); nbSC++; isBarrel = true; } } } // Endcap if (std::abs(track->eta()) > k_etaEndcap) { for (const auto& SC : *EndcapSupClusCollection) { dRSC = reco::deltaR(SC.eta(), SC.phi(), etaAtEcal, phiIn); if (dRSC < dRSC_first) { dRSC_first = dRSC; } // distance in eta-phi plan to closest SC else if (dRSC < dRSC_second) { dRSC_second = dRSC; } // to next closest SC if (dRSC < 0.09) { //Calculate phiWidth & etaWidth for associated SuperClusters SCShape.Calculate_Covariances(SC); phiWidth = SCShape.phiWidth(); etaWidth = SCShape.etaWidth(); etaSC = SC.eta(); phiSC = SC.phi(); algo_ID = SC.algoID(); EnergyEcal = SC.energy(); nBasicClus = SC.clustersSize(); nbSC++; isEndcap = true; } } } if (dRSC_first < 0.1) h_distToClosestSC->Fill(dRSC_first); if (nbSC == 1) h_counter1->Fill(0.5); if (nbSC == 0) h_counter2->Fill(0.5); if (nbSC > 1 || nbSC == 0) continue; h_cut_OneSCmatch->Fill(0.5); if (isBarrel && isEndcap) { edm::LogError("EopElecTreeWriter") << "Error: Super Cluster double matching!"; return; } h_EcalEnergy->Fill(EnergyEcal); HcalVSEcal->Fill(EnergyEcal, EnergyHcalIn02); // E over p plots if (track->charge() < 0) h_Eop_InnerNegative->Fill(EnergyEcal / pin); if (track->charge() > 0) h_Eop_InnerPositive->Fill(EnergyEcal / pin); //Check if track-SuperCluster matching is Ecal or Tracker driven edm::RefToBase<TrajectorySeed> seed = track->extra()->seedRef(); if (seed.isNull()) { edm::LogError("GsfElectronCore") << "The GsfTrack has no seed ?!"; } else { reco::ElectronSeedRef elseed = seed.castTo<reco::ElectronSeedRef>(); if (elseed.isNull()) { edm::LogError("GsfElectronCore") << "The GsfTrack seed is not an ElectronSeed ?!"; } else { isEcalDriven = elseed->isEcalDriven(); isTrackerDriven = elseed->isTrackerDriven(); } } treeMemPtr_->charge = track->charge(); treeMemPtr_->nHits = track->numberOfValidHits(); treeMemPtr_->nLostHits = track->numberOfLostHits(); treeMemPtr_->innerOk = track->innerOk(); treeMemPtr_->outerRadius = track->outerRadius(); treeMemPtr_->chi2 = track->chi2(); treeMemPtr_->normalizedChi2 = track->normalizedChi2(); treeMemPtr_->px = track->px(); treeMemPtr_->py = track->py(); treeMemPtr_->pz = track->pz(); treeMemPtr_->p = track->p(); treeMemPtr_->pIn = pin; treeMemPtr_->etaIn = etaIn; treeMemPtr_->phiIn = phiIn; treeMemPtr_->pOut = pout; treeMemPtr_->etaOut = etaOut; treeMemPtr_->phiOut = phiOut; treeMemPtr_->pt = track->pt(); treeMemPtr_->ptError = track->ptError(); treeMemPtr_->theta = track->theta(); treeMemPtr_->eta = track->eta(); treeMemPtr_->phi = track->phi(); treeMemPtr_->fbrem = fbrem; treeMemPtr_->MaxPtIn01 = MaxPtIn01; treeMemPtr_->SumPtIn01 = SumPtIn01; treeMemPtr_->NoTrackIn0015 = NoTrackIn0015; treeMemPtr_->MaxPtIn02 = MaxPtIn02; treeMemPtr_->SumPtIn02 = SumPtIn02; treeMemPtr_->NoTrackIn0020 = NoTrackIn0020; treeMemPtr_->MaxPtIn03 = MaxPtIn03; treeMemPtr_->SumPtIn03 = SumPtIn03; treeMemPtr_->NoTrackIn0025 = NoTrackIn0025; treeMemPtr_->MaxPtIn04 = MaxPtIn04; treeMemPtr_->SumPtIn04 = SumPtIn04; treeMemPtr_->NoTrackIn0030 = NoTrackIn0030; treeMemPtr_->MaxPtIn05 = MaxPtIn05; treeMemPtr_->SumPtIn05 = SumPtIn05; treeMemPtr_->NoTrackIn0035 = NoTrackIn0035; treeMemPtr_->NoTrackIn0040 = NoTrackIn0040; treeMemPtr_->SC_algoID = algo_ID; treeMemPtr_->SC_energy = EnergyEcal; treeMemPtr_->SC_nBasicClus = nBasicClus; treeMemPtr_->SC_etaWidth = etaWidth; treeMemPtr_->SC_phiWidth = phiWidth; treeMemPtr_->SC_eta = etaSC; treeMemPtr_->SC_phi = phiSC; treeMemPtr_->SC_isBarrel = isBarrel; treeMemPtr_->SC_isEndcap = isEndcap; treeMemPtr_->dRto1stSC = dRSC_first; treeMemPtr_->dRto2ndSC = dRSC_second; treeMemPtr_->HcalEnergyIn01 = EnergyHcalIn01; treeMemPtr_->HcalEnergyIn02 = EnergyHcalIn02; treeMemPtr_->HcalEnergyIn03 = EnergyHcalIn03; treeMemPtr_->HcalEnergyIn04 = EnergyHcalIn04; treeMemPtr_->HcalEnergyIn05 = EnergyHcalIn05; treeMemPtr_->isEcalDriven = isEcalDriven; treeMemPtr_->isTrackerDriven = isTrackerDriven; treeMemPtr_->RunNumber = iEvent.id().run(); treeMemPtr_->EvtNumber = iEvent.id().event(); tree_->Fill(); } // loop on tracks } // analyze // ------------ method called once each job just after ending the event loop ------------ void EopElecTreeWriter::endJob() { delete treeMemPtr_; treeMemPtr_ = nullptr; } // ------------ method called once each job just before starting event loop ------------ void EopElecTreeWriter::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) { bool changed = true; // Load trigger configuration at each begin of run // Fill the trigger names if (hltConfig_.init(iRun, iSetup, "HLT", changed)) { if (changed) { triggerNames_ = hltConfig_.triggerNames(); } } // Displaying the trigger names if (debugTriggerSelection_) { unsigned int i = 0; for (const auto& it : triggerNames_) { edm::LogInfo("EopElecTreeWriter") << "HLTpath: " << (i++) << " = " << it; } } } //************************************************************* void EopElecTreeWriter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) //************************************************************* { edm::ParameterSetDescription desc; desc.setComment("Generate tree for Tracker Alignment E/p validation"); desc.add<edm::InputTag>("src", edm::InputTag("electronGsfTracks")); desc.add<std::string>("triggerPath", "HLT_Ele"); desc.add<std::string>("hltFilter", "hltDiEle27L1DoubleEGWPTightHcalIsoFilter"); desc.add<bool>("debugTriggerSelection", false); descriptions.addWithDefaultLabel(desc); } //define this as a plug-in DEFINE_FWK_MODULE(EopElecTreeWriter);