/
githubmirror
/
cmssw
Обзор
Документация
Войти
/
githubmirror
/
cmssw
Код
Запросы
0
Пакеты
0
Релизы
0
Аналитика
Безопасность
master
Calibration/IsolatedParticles/plugins/IsolatedTracksHcalScale.cc
691 строка
32 KB
Matti Kortelainen
Remove TrackAssociatorRecord as obsolete
19 авг 2021, 03:57
19 авг 2021, 03:57
b19bdaa
Код
Авторство
О чём код?
// System include files #include <cmath> #include <memory> #include <map> #include <string> #include <vector> // root objects #include "TROOT.h" #include "TSystem.h" #include "TFile.h" #include "TH1F.h" #include "TH2F.h" #include "TProfile.h" #include "TDirectory.h" #include "TTree.h" #include <Math/GenVector/VectorUtil.h> #include "Calibration/IsolatedParticles/interface/FindCaloHit.h" #include "Calibration/IsolatedParticles/interface/eECALMatrix.h" #include "Calibration/IsolatedParticles/interface/eHCALMatrix.h" #include "Calibration/IsolatedParticles/interface/MatchingSimTrack.h" #include "Calibration/IsolatedParticles/interface/CaloSimInfo.h" #include "Calibration/IsolatedParticles/interface/TrackSelection.h" #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h" #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h" #include "Calibration/IsolatedParticles/interface/eCone.h" #include "Calibration/IsolatedParticles/interface/eHCALMatrixExtra.h" #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h" #include "DataFormats/Common/interface/Ref.h" #include "DataFormats/Common/interface/Handle.h" #include "DataFormats/Math/interface/Point3D.h" #include "DataFormats/Candidate/interface/Candidate.h" // muons and tracks #include "DataFormats/TrackReco/interface/Track.h" #include "DataFormats/TrackReco/interface/TrackFwd.h" #include "DataFormats/TrackReco/interface/HitPattern.h" #include "DataFormats/TrackReco/interface/TrackBase.h" // Vertices #include "DataFormats/VertexReco/interface/VertexFwd.h" #include "DataFormats/VertexReco/interface/Vertex.h" #include "DataFormats/BeamSpot/interface/BeamSpot.h" // Calorimeters #include "DataFormats/DetId/interface/DetId.h" #include "DataFormats/EcalDetId/interface/EBDetId.h" #include "DataFormats/EcalDetId/interface/EEDetId.h" #include "DataFormats/HcalDetId/interface/HcalDetId.h" #include "DataFormats/HcalDetId/interface/HcalSubdetector.h" #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h" #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h" #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetupFwd.h" // Jets in the event #include "DataFormats/JetReco/interface/CaloJet.h" #include "DataFormats/JetReco/interface/CaloJetCollection.h" #include "DataFormats/JetReco/interface/JetExtendedAssociation.h" #include "FWCore/Framework/interface/Frameworkfwd.h" #include "FWCore/Framework/interface/one/EDAnalyzer.h" #include "FWCore/Framework/interface/Event.h" #include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/Framework/interface/MakerMacros.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" // TFile Service #include "FWCore/ServiceRegistry/interface/Service.h" #include "CommonTools/UtilAlgos/interface/TFileService.h" #include "Geometry/CaloTopology/interface/EcalTrigTowerConstituentsMap.h" #include "Geometry/CaloGeometry/interface/CaloGeometry.h" #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h" #include "Geometry/Records/interface/CaloTopologyRecord.h" #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h" #include "Geometry/CaloTopology/interface/CaloTopology.h" #include "MagneticField/Engine/interface/MagneticField.h" #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h" #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h" #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h" #include "RecoCaloTools/Navigation/interface/CaloNavigator.h" // SimHit #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h" #include "SimDataFormats/Track/interface/SimTrack.h" #include "SimDataFormats/Track/interface/SimTrackContainer.h" #include "SimDataFormats/Vertex/interface/SimVertex.h" #include "SimDataFormats/Vertex/interface/SimVertexContainer.h" // track associator #include "SimTracker/TrackerHitAssociation/interface/TrackerHitAssociator.h" #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h" #include "TrackingTools/TrackAssociator/interface/TrackDetectorAssociator.h" class IsolatedTracksHcalScale : public edm::one::EDAnalyzer<edm::one::SharedResources> { public: explicit IsolatedTracksHcalScale(const edm::ParameterSet &); ~IsolatedTracksHcalScale() override {} static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); private: void beginJob() override; void analyze(const edm::Event &, const edm::EventSetup &) override; void endJob() override; void clearTreeVectors(); private: bool doMC_; int myverbose_; std::string theTrackQuality_, minQuality_; spr::trackSelectionParameters selectionParameters_; double a_mipR_, a_coneR_, a_charIsoR_, a_neutIsoR_; double tMinE_, tMaxE_; TrackerHitAssociator::Config trackerHitAssociatorConfig_; edm::EDGetTokenT<reco::TrackCollection> tok_genTrack_; edm::EDGetTokenT<reco::VertexCollection> tok_recVtx_; edm::EDGetTokenT<reco::BeamSpot> tok_bs_; edm::EDGetTokenT<EcalRecHitCollection> tok_EB_; edm::EDGetTokenT<EcalRecHitCollection> tok_EE_; edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_; edm::EDGetTokenT<edm::SimTrackContainer> tok_simTk_; edm::EDGetTokenT<edm::SimVertexContainer> tok_simVtx_; edm::EDGetTokenT<edm::PCaloHitContainer> tok_caloEB_; edm::EDGetTokenT<edm::PCaloHitContainer> tok_caloEE_; edm::EDGetTokenT<edm::PCaloHitContainer> tok_caloHH_; edm::ESGetToken<CaloGeometry, CaloGeometryRecord> tok_geom_; edm::ESGetToken<CaloTopology, CaloTopologyRecord> tok_caloTopology_; edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> tok_magField_; edm::ESGetToken<EcalChannelStatus, EcalChannelStatusRcd> tok_ecalChStatus_; edm::ESGetToken<EcalSeverityLevelAlgo, EcalSeverityLevelAlgoRcd> tok_sevlv_; int nEventProc_; TTree *tree_; int t_nTracks, t_RunNo, t_EvtNo, t_Lumi, t_Bunch; std::vector<double> *t_trackP, *t_trackPt, *t_trackEta, *t_trackPhi; std::vector<double> *t_trackHcalEta, *t_trackHcalPhi, *t_eHCALDR; std::vector<double> *t_hCone, *t_conehmaxNearP, *t_eMipDR, *t_eECALDR; std::vector<double> *t_e11x11_20Sig, *t_e15x15_20Sig; std::vector<double> *t_eMipDR_1, *t_eECALDR_1, *t_eMipDR_2, *t_eECALDR_2; std::vector<double> *t_hConeHB, *t_eHCALDRHB; std::vector<double> *t_hsimInfoMatched, *t_hsimInfoRest, *t_hsimInfoPhoton; std::vector<double> *t_hsimInfoNeutHad, *t_hsimInfoCharHad, *t_hsimInfoPdgMatched; std::vector<double> *t_hsimInfoTotal, *t_hsim; std::vector<int> *t_hsimInfoNMatched, *t_hsimInfoNTotal, *t_hsimInfoNNeutHad; std::vector<int> *t_hsimInfoNCharHad, *t_hsimInfoNPhoton, *t_hsimInfoNRest; std::vector<int> *t_nSimHits; }; IsolatedTracksHcalScale::IsolatedTracksHcalScale(const edm::ParameterSet &iConfig) : doMC_(iConfig.getUntrackedParameter<bool>("DoMC", false)), myverbose_(iConfig.getUntrackedParameter<int>("Verbosity", 5)), theTrackQuality_(iConfig.getUntrackedParameter<std::string>("TrackQuality", "highPurity")), a_mipR_(iConfig.getUntrackedParameter<double>("ConeRadiusMIP", 14.0)), a_coneR_(iConfig.getUntrackedParameter<double>("ConeRadius", 34.98)), tMinE_(iConfig.getUntrackedParameter<double>("TimeMinCutECAL", -500.)), tMaxE_(iConfig.getUntrackedParameter<double>("TimeMaxCutECAL", 500.)), trackerHitAssociatorConfig_(consumesCollector()) { usesResource(TFileService::kSharedResource); //now do what ever initialization is needed reco::TrackBase::TrackQuality trackQuality = reco::TrackBase::qualityByName(theTrackQuality_); selectionParameters_.minPt = iConfig.getUntrackedParameter<double>("MinTrackPt", 10.0); selectionParameters_.minQuality = trackQuality; selectionParameters_.maxDxyPV = iConfig.getUntrackedParameter<double>("MaxDxyPV", 0.2); selectionParameters_.maxDzPV = iConfig.getUntrackedParameter<double>("MaxDzPV", 5.0); selectionParameters_.maxChi2 = iConfig.getUntrackedParameter<double>("MaxChi2", 5.0); selectionParameters_.maxDpOverP = iConfig.getUntrackedParameter<double>("MaxDpOverP", 0.1); selectionParameters_.minOuterHit = iConfig.getUntrackedParameter<int>("MinOuterHit", 4); selectionParameters_.minLayerCrossed = iConfig.getUntrackedParameter<int>("MinLayerCrossed", 8); selectionParameters_.maxInMiss = iConfig.getUntrackedParameter<int>("MaxInMiss", 0); selectionParameters_.maxOutMiss = iConfig.getUntrackedParameter<int>("MaxOutMiss", 0); a_charIsoR_ = a_coneR_ + 28.9; a_neutIsoR_ = a_charIsoR_ * 0.726; tok_genTrack_ = consumes<reco::TrackCollection>(edm::InputTag("generalTracks")); tok_recVtx_ = consumes<reco::VertexCollection>(edm::InputTag("offlinePrimaryVertices")); tok_bs_ = consumes<reco::BeamSpot>(edm::InputTag("offlineBeamSpot")); tok_EB_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit", "EcalRecHitsEB")); tok_EE_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit", "EcalRecHitsEE")); tok_hbhe_ = consumes<HBHERecHitCollection>(edm::InputTag("hbhereco")); tok_simTk_ = consumes<edm::SimTrackContainer>(edm::InputTag("g4SimHits")); tok_simVtx_ = consumes<edm::SimVertexContainer>(edm::InputTag("g4SimHits")); tok_caloEB_ = consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsEB")); tok_caloEE_ = consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsEE")); tok_caloHH_ = consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "HcalHits")); if (myverbose_ >= 0) { edm::LogVerbatim("IsoTrack") << "Parameters read from config file \n" << " doMC " << doMC_ << "\t myverbose " << myverbose_ << "\t minPt " << selectionParameters_.minPt << "\t theTrackQuality " << theTrackQuality_ << "\t minQuality " << selectionParameters_.minQuality << "\t maxDxyPV " << selectionParameters_.maxDxyPV << "\t maxDzPV " << selectionParameters_.maxDzPV << "\t maxChi2 " << selectionParameters_.maxChi2 << "\t maxDpOverP " << selectionParameters_.maxDpOverP << "\t minOuterHit " << selectionParameters_.minOuterHit << "\t minLayerCrossed " << selectionParameters_.minLayerCrossed << "\t maxInMiss " << selectionParameters_.maxInMiss << "\t maxOutMiss " << selectionParameters_.maxOutMiss << "\t a_coneR " << a_coneR_ << "\t a_charIsoR " << a_charIsoR_ << "\t a_neutIsoR " << a_neutIsoR_ << "\t a_mipR " << a_mipR_ << "\t time Range (" << tMinE_ << ":" << tMaxE_ << ")"; } tok_geom_ = esConsumes<CaloGeometry, CaloGeometryRecord>(); tok_caloTopology_ = esConsumes<CaloTopology, CaloTopologyRecord>(); tok_magField_ = esConsumes<MagneticField, IdealMagneticFieldRecord>(); tok_ecalChStatus_ = esConsumes<EcalChannelStatus, EcalChannelStatusRcd>(); tok_sevlv_ = esConsumes<EcalSeverityLevelAlgo, EcalSeverityLevelAlgoRcd>(); } void IsolatedTracksHcalScale::fillDescriptions(edm::ConfigurationDescriptions &descriptions) { edm::ParameterSetDescription desc; desc.addUntracked<bool>("doMC", false); desc.addUntracked<int>("Verbosity", 0); desc.addUntracked<std::string>("TrackQuality", "highPurity"); desc.addUntracked<double>("MinTrackPt", 10.0); desc.addUntracked<double>("MaxDxyPV", 0.02); desc.addUntracked<double>("MaxDzPV", 0.02); desc.addUntracked<double>("MaxChi2", 5.0); desc.addUntracked<double>("MaxDpOverP", 0.1); desc.addUntracked<int>("MinOuterHit", 4); desc.addUntracked<int>("MinLayerCrossed", 8); desc.addUntracked<int>("MaxInMiss", 0); desc.addUntracked<int>("MaxOutMiss", 0); desc.addUntracked<double>("ConeRadius", 34.98); desc.addUntracked<double>("ConeRadiusMIP", 14.0); desc.addUntracked<double>("TimeMinCutECAL", -500.0); desc.addUntracked<double>("TimeMaxCutECAL", 500.0); descriptions.add("isolatedTracksHcalScale", desc); } void IsolatedTracksHcalScale::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) { const CaloGeometry *geo = &iSetup.getData(tok_geom_); const MagneticField *bField = &iSetup.getData(tok_magField_); const EcalChannelStatus *theEcalChStatus = &iSetup.getData(tok_ecalChStatus_); const EcalSeverityLevelAlgo *sevlv = &iSetup.getData(tok_sevlv_); const CaloTopology *caloTopology = &iSetup.getData(tok_caloTopology_); clearTreeVectors(); ++nEventProc_; t_RunNo = iEvent.id().run(); t_EvtNo = iEvent.id().event(); t_Lumi = iEvent.luminosityBlock(); t_Bunch = iEvent.bunchCrossing(); if (myverbose_ > 0) edm::LogVerbatim("IsoTrack") << nEventProc_ << " Run " << t_RunNo << " Event " << t_EvtNo << " Lumi " << t_Lumi << " Bunch " << t_Bunch; edm::Handle<reco::TrackCollection> trkCollection; iEvent.getByToken(tok_genTrack_, trkCollection); edm::Handle<reco::VertexCollection> recVtxs; iEvent.getByToken(tok_recVtx_, recVtxs); // Get the beamspot edm::Handle<reco::BeamSpot> beamSpotH; iEvent.getByToken(tok_bs_, beamSpotH); math::XYZPoint leadPV(0, 0, 0); if (!recVtxs->empty() && !((*recVtxs)[0].isFake())) { leadPV = math::XYZPoint((*recVtxs)[0].x(), (*recVtxs)[0].y(), (*recVtxs)[0].z()); } else if (beamSpotH.isValid()) { leadPV = beamSpotH->position(); } if (myverbose_ > 0) { edm::LogVerbatim("IsoTrack") << "Primary Vertex " << leadPV; if (beamSpotH.isValid()) edm::LogVerbatim("IsoTrack") << "Beam Spot " << beamSpotH->position(); } std::vector<spr::propagatedTrackDirection> trkCaloDirections; spr::propagateCALO(trkCollection, geo, bField, theTrackQuality_, trkCaloDirections, (myverbose_ > 2)); std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr; edm::Handle<EcalRecHitCollection> barrelRecHitsHandle; edm::Handle<EcalRecHitCollection> endcapRecHitsHandle; iEvent.getByToken(tok_EB_, barrelRecHitsHandle); iEvent.getByToken(tok_EE_, endcapRecHitsHandle); edm::Handle<HBHERecHitCollection> hbhe; iEvent.getByToken(tok_hbhe_, hbhe); const HBHERecHitCollection Hithbhe = *(hbhe.product()); //get Handles to SimTracks and SimHits edm::Handle<edm::SimTrackContainer> SimTk; edm::Handle<edm::SimVertexContainer> SimVtx; //get Handles to PCaloHitContainers of eb/ee/hbhe edm::Handle<edm::PCaloHitContainer> pcaloeb; edm::Handle<edm::PCaloHitContainer> pcaloee; edm::Handle<edm::PCaloHitContainer> pcalohh; //associates tracker rechits/simhits to a track std::unique_ptr<TrackerHitAssociator> associate; if (doMC_) { iEvent.getByToken(tok_simTk_, SimTk); iEvent.getByToken(tok_simVtx_, SimVtx); iEvent.getByToken(tok_caloEB_, pcaloeb); iEvent.getByToken(tok_caloEE_, pcaloee); iEvent.getByToken(tok_caloHH_, pcalohh); associate = std::make_unique<TrackerHitAssociator>(iEvent, trackerHitAssociatorConfig_); } unsigned int nTracks = 0; for (trkDetItr = trkCaloDirections.begin(), nTracks = 0; trkDetItr != trkCaloDirections.end(); trkDetItr++, nTracks++) { const reco::Track *pTrack = &(*(trkDetItr->trkItr)); if (spr::goodTrack(pTrack, leadPV, selectionParameters_, (myverbose_ > 2)) && trkDetItr->okECAL && trkDetItr->okHCAL) { int nRH_eMipDR = 0, nRH_eDR = 0, nNearTRKs = 0, nRecHitsCone = -99; double distFromHotCell = -99.0, distFromHotCell2 = -99.0; int ietaHotCell = -99, iphiHotCell = -99; int ietaHotCell2 = -99, iphiHotCell2 = -99; GlobalPoint gposHotCell(0., 0., 0.), gposHotCell2(0., 0., 0.); std::vector<DetId> coneRecHitDetIds, coneRecHitDetIds2; std::pair<double, bool> e11x11_20SigP, e15x15_20SigP; double hCone = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_coneR_, trkDetItr->directionHCAL, nRecHitsCone, coneRecHitDetIds, distFromHotCell, ietaHotCell, iphiHotCell, gposHotCell, -1); double hConeHB = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_coneR_, trkDetItr->directionHCAL, nRecHitsCone, coneRecHitDetIds, distFromHotCell, ietaHotCell, iphiHotCell, gposHotCell, (int)(HcalBarrel)); double eHCALDR = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_charIsoR_, trkDetItr->directionHCAL, nRecHitsCone, coneRecHitDetIds2, distFromHotCell2, ietaHotCell2, iphiHotCell2, gposHotCell2, -1); double eHCALDRHB = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_charIsoR_, trkDetItr->directionHCAL, nRecHitsCone, coneRecHitDetIds2, distFromHotCell2, ietaHotCell2, iphiHotCell2, gposHotCell2, (int)(HcalBarrel)); double conehmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections, a_charIsoR_, nNearTRKs, (myverbose_ > 3)); double eMipDR = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_mipR_, trkDetItr->directionECAL, nRH_eMipDR); double eECALDR = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_neutIsoR_, trkDetItr->directionECAL, nRH_eDR); double eMipDR_1 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_mipR_, trkDetItr->directionECAL, nRH_eMipDR, 0.030, 0.150); double eECALDR_1 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_neutIsoR_, trkDetItr->directionECAL, nRH_eDR, 0.030, 0.150); double eMipDR_2 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_mipR_, trkDetItr->directionECAL, nRH_eMipDR, 0.060, 0.300); double eECALDR_2 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_neutIsoR_, trkDetItr->directionECAL, nRH_eDR, 0.060, 0.300); HcalDetId closestCell = (HcalDetId)(trkDetItr->detIdHCAL); e11x11_20SigP = spr::eECALmatrix(trkDetItr->detIdECAL, barrelRecHitsHandle, endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv, 5, 5, 0.060, 0.300, tMinE_, tMaxE_); e15x15_20SigP = spr::eECALmatrix(trkDetItr->detIdECAL, barrelRecHitsHandle, endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv, 7, 7, 0.060, 0.300, tMinE_, tMaxE_); // Fill the tree Branches here t_trackP->push_back(pTrack->p()); t_trackPt->push_back(pTrack->pt()); t_trackEta->push_back(pTrack->momentum().eta()); t_trackPhi->push_back(pTrack->momentum().phi()); t_trackHcalEta->push_back(closestCell.ieta()); t_trackHcalPhi->push_back(closestCell.iphi()); t_hCone->push_back(hCone); t_conehmaxNearP->push_back(conehmaxNearP); t_eMipDR->push_back(eMipDR); t_eECALDR->push_back(eECALDR); t_eHCALDR->push_back(eHCALDR); t_e11x11_20Sig->push_back(e11x11_20SigP.first); t_e15x15_20Sig->push_back(e15x15_20SigP.first); t_eMipDR_1->push_back(eMipDR_1); t_eECALDR_1->push_back(eECALDR_1); t_eMipDR_2->push_back(eMipDR_2); t_eECALDR_2->push_back(eECALDR_2); t_hConeHB->push_back(hConeHB); t_eHCALDRHB->push_back(eHCALDRHB); if (myverbose_ > 0) { edm::LogVerbatim("IsoTrack") << "Track p " << pTrack->p() << " pt " << pTrack->pt() << " eta " << pTrack->momentum().eta() << " phi " << pTrack->momentum().phi() << " ieta/iphi (" << closestCell.ieta() << ", " << closestCell.iphi() << ") Energy in cone " << hCone << " Charge Isolation " << conehmaxNearP << " eMIP (" << eMipDR << ", " << eMipDR_1 << ", " << eMipDR_2 << ")" << " Neutral isolation (ECAL) (" << eECALDR - eMipDR << ", " << eECALDR_1 - eMipDR_1 << ", " << eECALDR_2 - eMipDR_2 << ") (ECAL NxN) " << e15x15_20SigP.first - e11x11_20SigP.first << " (HCAL) " << eHCALDR - hCone; } if (doMC_) { int nSimHits = -999; double hsim; std::map<std::string, double> hsimInfo; std::vector<int> multiplicity; hsim = spr::eCone_hcal( geo, pcalohh, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_coneR_, trkDetItr->directionHCAL, nSimHits); hsimInfo = spr::eHCALSimInfoCone(iEvent, pcalohh, SimTk, SimVtx, pTrack, *associate, geo, trkDetItr->pointHCAL, trkDetItr->pointECAL, a_coneR_, trkDetItr->directionHCAL, multiplicity); t_hsimInfoMatched->push_back(hsimInfo["eMatched"]); t_hsimInfoRest->push_back(hsimInfo["eRest"]); t_hsimInfoPhoton->push_back(hsimInfo["eGamma"]); t_hsimInfoNeutHad->push_back(hsimInfo["eNeutralHad"]); t_hsimInfoCharHad->push_back(hsimInfo["eChargedHad"]); t_hsimInfoPdgMatched->push_back(hsimInfo["pdgMatched"]); t_hsimInfoTotal->push_back(hsimInfo["eTotal"]); t_hsimInfoNMatched->push_back(multiplicity.at(0)); t_hsimInfoNTotal->push_back(multiplicity.at(1)); t_hsimInfoNNeutHad->push_back(multiplicity.at(2)); t_hsimInfoNCharHad->push_back(multiplicity.at(3)); t_hsimInfoNPhoton->push_back(multiplicity.at(4)); t_hsimInfoNRest->push_back(multiplicity.at(5)); t_hsim->push_back(hsim); t_nSimHits->push_back(nSimHits); if (myverbose_ > 0) { edm::LogVerbatim("IsoTrack") << "Matched (E) " << hsimInfo["eMatched"] << " (N) " << multiplicity.at(0) << " Rest (E) " << hsimInfo["eRest"] << " (N) " << multiplicity.at(1) << " Gamma (E) " << hsimInfo["eGamma"] << " (N) " << multiplicity.at(2) << " Neutral Had (E) " << hsimInfo["eNeutralHad"] << " (N) " << multiplicity.at(3) << " Charged Had (E) " << hsimInfo["eChargedHad"] << " (N) " << multiplicity.at(4) << " Total (E) " << hsimInfo["eTotal"] << " (N) " << multiplicity.at(5) << " PDG " << hsimInfo["pdgMatched"] << " Total E " << hsim << " NHit " << nSimHits; } } } } tree_->Fill(); } void IsolatedTracksHcalScale::beginJob() { nEventProc_ = 0; edm::Service<TFileService> fs; //////Book Tree tree_ = fs->make<TTree>("tree", "tree"); tree_->SetAutoSave(10000); tree_->Branch("t_RunNo", &t_RunNo, "t_RunNo/I"); tree_->Branch("t_Lumi", &t_Lumi, "t_Lumi/I"); tree_->Branch("t_Bunch", &t_Bunch, "t_Bunch/I"); t_trackP = new std::vector<double>(); t_trackPt = new std::vector<double>(); t_trackEta = new std::vector<double>(); t_trackPhi = new std::vector<double>(); t_trackHcalEta = new std::vector<double>(); t_trackHcalPhi = new std::vector<double>(); t_hCone = new std::vector<double>(); t_conehmaxNearP = new std::vector<double>(); t_eMipDR = new std::vector<double>(); t_eECALDR = new std::vector<double>(); t_eHCALDR = new std::vector<double>(); t_e11x11_20Sig = new std::vector<double>(); t_e15x15_20Sig = new std::vector<double>(); t_eMipDR_1 = new std::vector<double>(); t_eECALDR_1 = new std::vector<double>(); t_eMipDR_2 = new std::vector<double>(); t_eECALDR_2 = new std::vector<double>(); t_hConeHB = new std::vector<double>(); t_eHCALDRHB = new std::vector<double>(); tree_->Branch("t_trackP", "std::vector<double>", &t_trackP); tree_->Branch("t_trackPt", "std::vector<double>", &t_trackPt); tree_->Branch("t_trackEta", "std::vector<double>", &t_trackEta); tree_->Branch("t_trackPhi", "std::vector<double>", &t_trackPhi); tree_->Branch("t_trackHcalEta", "std::vector<double>", &t_trackHcalEta); tree_->Branch("t_trackHcalPhi", "std::vector<double>", &t_trackHcalPhi); tree_->Branch("t_hCone", "std::vector<double>", &t_hCone); tree_->Branch("t_conehmaxNearP", "std::vector<double>", &t_conehmaxNearP); tree_->Branch("t_eMipDR", "std::vector<double>", &t_eMipDR); tree_->Branch("t_eECALDR", "std::vector<double>", &t_eECALDR); tree_->Branch("t_eHCALDR", "std::vector<double>", &t_eHCALDR); tree_->Branch("t_e11x11_20Sig", "std::vector<double>", &t_e11x11_20Sig); tree_->Branch("t_e15x15_20Sig", "std::vector<double>", &t_e15x15_20Sig); tree_->Branch("t_eMipDR_1", "std::vector<double>", &t_eMipDR_1); tree_->Branch("t_eECALDR_1", "std::vector<double>", &t_eECALDR_1); tree_->Branch("t_eMipDR_2", "std::vector<double>", &t_eMipDR_2); tree_->Branch("t_eECALDR_2", "std::vector<double>", &t_eECALDR_2); tree_->Branch("t_hConeHB", "std::vector<double>", &t_hConeHB); tree_->Branch("t_eHCALDRHB", "std::vector<double>", &t_eHCALDRHB); if (doMC_) { t_hsimInfoMatched = new std::vector<double>(); t_hsimInfoRest = new std::vector<double>(); t_hsimInfoPhoton = new std::vector<double>(); t_hsimInfoNeutHad = new std::vector<double>(); t_hsimInfoCharHad = new std::vector<double>(); t_hsimInfoPdgMatched = new std::vector<double>(); t_hsimInfoTotal = new std::vector<double>(); t_hsimInfoNMatched = new std::vector<int>(); t_hsimInfoNTotal = new std::vector<int>(); t_hsimInfoNNeutHad = new std::vector<int>(); t_hsimInfoNCharHad = new std::vector<int>(); t_hsimInfoNPhoton = new std::vector<int>(); t_hsimInfoNRest = new std::vector<int>(); t_hsim = new std::vector<double>(); t_nSimHits = new std::vector<int>(); tree_->Branch("t_hsimInfoMatched", "std::vector<double>", &t_hsimInfoMatched); tree_->Branch("t_hsimInfoRest", "std::vector<double>", &t_hsimInfoRest); tree_->Branch("t_hsimInfoPhoton", "std::vector<double>", &t_hsimInfoPhoton); tree_->Branch("t_hsimInfoNeutHad", "std::vector<double>", &t_hsimInfoNeutHad); tree_->Branch("t_hsimInfoCharHad", "std::vector<double>", &t_hsimInfoCharHad); tree_->Branch("t_hsimInfoPdgMatched", "std::vector<double>", &t_hsimInfoPdgMatched); tree_->Branch("t_hsimInfoTotal", "std::vector<double>", &t_hsimInfoTotal); tree_->Branch("t_hsimInfoNMatched", "std::vector<int>", &t_hsimInfoNMatched); tree_->Branch("t_hsimInfoNTotal", "std::vector<int>", &t_hsimInfoNTotal); tree_->Branch("t_hsimInfoNNeutHad", "std::vector<int>", &t_hsimInfoNNeutHad); tree_->Branch("t_hsimInfoNCharHad", "std::vector<int>", &t_hsimInfoNCharHad); tree_->Branch("t_hsimInfoNPhoton", "std::vector<int>", &t_hsimInfoNPhoton); tree_->Branch("t_hsimInfoNRest", "std::vector<int>", &t_hsimInfoNRest); tree_->Branch("t_hsim", "std::vector<double>", &t_hsim); tree_->Branch("t_nSimHits", "std::vector<int>", &t_nSimHits); } } void IsolatedTracksHcalScale::endJob() { edm::LogVerbatim("IsoTrack") << "Number of Events Processed " << nEventProc_; } void IsolatedTracksHcalScale::clearTreeVectors() { t_trackP->clear(); t_trackPt->clear(); t_trackEta->clear(); t_trackPhi->clear(); t_trackHcalEta->clear(); t_trackHcalPhi->clear(); t_hCone->clear(); t_conehmaxNearP->clear(); t_eMipDR->clear(); t_eECALDR->clear(); t_eHCALDR->clear(); t_e11x11_20Sig->clear(); t_e15x15_20Sig->clear(); t_eMipDR_1->clear(); t_eECALDR_1->clear(); t_eMipDR_2->clear(); t_eECALDR_2->clear(); t_hConeHB->clear(); t_eHCALDRHB->clear(); if (doMC_) { t_hsimInfoMatched->clear(); t_hsimInfoRest->clear(); t_hsimInfoPhoton->clear(); t_hsimInfoNeutHad->clear(); t_hsimInfoCharHad->clear(); t_hsimInfoPdgMatched->clear(); t_hsimInfoTotal->clear(); t_hsimInfoNMatched->clear(); t_hsimInfoNTotal->clear(); t_hsimInfoNNeutHad->clear(); t_hsimInfoNCharHad->clear(); t_hsimInfoNPhoton->clear(); t_hsimInfoNRest->clear(); t_hsim->clear(); t_nSimHits->clear(); } } //define this as a plug-in DEFINE_FWK_MODULE(IsolatedTracksHcalScale);