27 #include <Math/GenVector/VectorUtil.h>
112 #include "TProfile.h"
113 #include "TDirectory.h"
214 std::unique_ptr<std::vector<double> >
t_trkP;
215 std::unique_ptr<std::vector<double> >
t_simP;
216 std::unique_ptr<std::vector<double> >
t_trkPt;
219 std::unique_ptr<std::vector<double> >
t_e3x3;
221 std::unique_ptr<std::vector<std::vector<double> > >
t_v_eDR;
222 std::unique_ptr<std::vector<std::vector<double> > >
t_v_eMipDR;
224 std::unique_ptr<std::vector<double> >
t_h3x3;
225 std::unique_ptr<std::vector<double> >
t_h5x5;
287 std::unique_ptr<std::vector<std::vector<double> > >
t_v_hCone;
326 std::unique_ptr<std::vector<unsigned int> >
t_irun;
327 std::unique_ptr<std::vector<unsigned int> >
t_ievt;
328 std::unique_ptr<std::vector<unsigned int> >
t_ilum;
332 : doMC_(iConfig.getUntrackedParameter<
bool>(
"doMC",
false)),
333 myverbose_(iConfig.getUntrackedParameter<
int>(
"verbosity", 5)),
334 useJetTrigger_(iConfig.getUntrackedParameter<
bool>(
"useJetTrigger",
false)),
335 drLeadJetVeto_(iConfig.getUntrackedParameter<double>(
"drLeadJetVeto", 1.2)),
336 ptMinLeadJet_(iConfig.getUntrackedParameter<double>(
"ptMinLeadJet", 15.0)),
337 debugTrks_(iConfig.getUntrackedParameter<
int>(
"debugTracks")),
338 printTrkHitPattern_(iConfig.getUntrackedParameter<
bool>(
"printTrkHitPattern")),
339 trackerHitAssociatorConfig_(consumesCollector()),
361 minTrackP_(iConfig.getUntrackedParameter<double>(
"minTrackP", 10.0)),
362 maxTrackEta_(iConfig.getUntrackedParameter<double>(
"maxTrackEta", 5.0)),
363 maxNearTrackP_(iConfig.getUntrackedParameter<double>(
"maxNearTrackP", 1.0)),
364 debugEcalSimInfo_(iConfig.getUntrackedParameter<
int>(
"debugEcalSimInfo")),
365 applyEcalIsolation_(iConfig.getUntrackedParameter<
bool>(
"applyEcalIsolation")) {
389 desc.addUntracked<
bool>(
"doMC",
false);
390 desc.addUntracked<
int>(
"verbosity", 1);
391 desc.addUntracked<
bool>(
"useJetTrigger",
false);
392 desc.addUntracked<
double>(
"drLeadJetVeto", 1.2);
393 desc.addUntracked<
double>(
"ptMinLeadJet", 15.0);
394 desc.addUntracked<
int>(
"debugTracks", 0);
395 desc.addUntracked<
bool>(
"printTrkHitPattern",
true);
396 desc.addUntracked<
double>(
"minTrackP", 1.0);
397 desc.addUntracked<
double>(
"maxTrackEta", 2.6);
398 desc.addUntracked<
double>(
"maxNearTrackP", 1.0);
399 desc.addUntracked<
bool>(
"debugEcalSimInfo",
false);
400 desc.addUntracked<
bool>(
"applyEcalIsolation",
true);
401 desc.addUntracked<
bool>(
"debugL1Info",
false);
407 desc_TrackAssoc.
add<
double>(
"muonMaxDistanceSigmaX", 0.0);
408 desc_TrackAssoc.add<
double>(
"muonMaxDistanceSigmaY", 0.0);
410 desc_TrackAssoc.add<
double>(
"dRHcal", 9999.0);
411 desc_TrackAssoc.add<
double>(
"dREcal", 9999.0);
413 desc_TrackAssoc.add<
bool>(
"useEcal",
true);
414 desc_TrackAssoc.add<
double>(
"dREcalPreselection", 0.05);
416 desc_TrackAssoc.add<
double>(
"dRMuon", 9999.0);
417 desc_TrackAssoc.add<
std::string>(
"crossedEnergyType",
"SinglePointAlongTrajectory");
418 desc_TrackAssoc.add<
double>(
"muonMaxDistanceX", 5.0);
419 desc_TrackAssoc.add<
double>(
"muonMaxDistanceY", 5.0);
420 desc_TrackAssoc.add<
bool>(
"useHO",
false);
421 desc_TrackAssoc.add<
bool>(
"accountForTrajectoryChangeCalo",
false);
424 desc_TrackAssoc.add<
double>(
"dRHcalPreselection", 0.2);
425 desc_TrackAssoc.add<
bool>(
"useMuon",
false);
426 desc_TrackAssoc.add<
bool>(
"useCalo",
true);
428 desc_TrackAssoc.add<
double>(
"dRMuonPreselection", 0.2);
429 desc_TrackAssoc.add<
bool>(
"truthMatch",
false);
431 desc_TrackAssoc.add<
bool>(
"useHcal",
true);
432 desc_TrackAssoc.add<
bool>(
"usePreshower",
false);
433 desc_TrackAssoc.add<
double>(
"dRPreshowerPreselection", 0.2);
434 desc_TrackAssoc.add<
double>(
"trajectoryUncertaintyTolerance", 1.0);
436 descriptions.
add(
"isolatedTracksCone",
desc);
440 unsigned int irun = (
unsigned int)
iEvent.id().run();
441 unsigned int ilum = (
unsigned int)
iEvent.getLuminosityBlock().luminosityBlock();
442 unsigned int ievt = (
unsigned int)
iEvent.id().event();
454 l1extra::L1JetParticleCollection::const_iterator itr;
455 for (itr = l1TauHandle->begin(); itr != l1TauHandle->end(); ++itr) {
464 for (itr = l1CenJetHandle->begin(); itr != l1CenJetHandle->end(); ++itr) {
473 for (itr = l1FwdJetHandle->begin(); itr != l1FwdJetHandle->end(); ++itr) {
521 reco::TrackCollection::const_iterator trkItr;
524 edm::LogVerbatim(
"IsoTrack") <<
"Number of Tracks " << trkCollection->size();
558 std::vector<int> v_hlTriggers;
568 int hltL1SingleEG5(-99);
569 int hltZeroBias(-99);
570 int hltMinBiasHcal(-99);
571 int hltMinBiasEcal(-99);
572 int hltMinBiasPixel(-99);
573 int hltSingleIsoTau30_Trk5(-99);
574 int hltDoubleLooseIsoTau15_Trk5(-99);
610 if (
triggerNames.triggerName(
i) ==
"HLT_SingleIsoTau30_Trk5")
612 if (
triggerNames.triggerName(
i) ==
"HLT_DoubleLooseIsoTau15_Trk5")
620 std::unique_ptr<TrackerHitAssociator> associate;
633 for (trkItr = trkCollection->begin(); trkItr != trkCollection->end(); ++trkItr) {
642 bool trkQual = pTrack->
quality(trackQuality_);
647 double phi1 = pTrack->
momentum().phi();
648 double pt1 = pTrack->
pt();
649 double p1 = pTrack->
p();
666 if (!goodEta || !goodPt || !trkQual)
695 int nOuterHits = hitp.stripTOBLayersWithMeasurement() + hitp.stripTECLayersWithMeasurement();
699 edm::SimTrackContainer::const_iterator matchedSimTrk =
701 simP = matchedSimTrk->momentum().P();
726 double drFromLeadJet = 999.0;
730 drFromLeadJet =
sqrt(dphi * dphi + deta * deta);
738 const int a_size = 7;
739 double a_coneR[a_size];
740 double a_charIsoR[a_size];
741 double a_neutIsoR[a_size];
751 for (
int i = 0;
i < a_size;
i++) {
752 a_charIsoR[
i] = a_coneR[
i] + 28.9;
753 a_neutIsoR[
i] = a_charIsoR[
i] * 0.726;
761 double e3x3 = -999.0;
762 double trkEcalEne = -999.0;
764 if (
std::abs(point1.eta()) < 1.479) {
767 isoCell, barrelRecHitsHandle, endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv, 1, 1)
773 isoCell, barrelRecHitsHandle, endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv, 1, 1)
781 const int a_mip_size = 5;
782 double a_mipR[a_mip_size];
789 std::vector<double> v_eDR;
790 for (
int i = 0;
i < a_size;
i++) {
795 geo, barrelRecHitsHandle, endcapRecHitsHandle, hpoint1, point1, a_neutIsoR[
i], trackMomAtEcal, nRH_eDR);
796 v_eDR.push_back(eDR);
799 std::vector<double> v_eMipDR;
800 for (
int i = 0;
i < a_mip_size;
i++) {
803 geo, barrelRecHitsHandle, endcapRecHitsHandle, hpoint1, point1, a_mipR[
i], trackMomAtEcal, nRH_eMipDR);
805 v_eMipDR.push_back(eMipDR);
813 std::vector<double> v_hmaxNearP_goodTrk;
814 std::vector<double> v_hmaxNearP;
815 std::vector<int> v_hnNearTRKs;
816 std::vector<int> v_hnLayers_maxNearP;
817 std::vector<int> v_htrkQual_maxNearP;
819 std::vector<double> v_cone_hmaxNearP_goodTrk;
820 std::vector<double> v_cone_hmaxNearP;
821 std::vector<int> v_cone_hnNearTRKs;
822 std::vector<int> v_cone_hnLayers_maxNearP;
823 std::vector<int> v_cone_htrkQual_maxNearP;
825 for (
int i = 0;
i < a_size;
i++) {
826 double hmaxNearP = -999.0;
828 int hnLayers_maxNearP = 0;
829 int htrkQual_maxNearP = -1;
830 double hmaxNearP_goodTrk = -999.0;
832 double conehmaxNearP = -999.0;
833 int conehnNearTRKs = 0;
834 int conehnLayers_maxNearP = 0;
835 int conehtrkQual_maxNearP = -1;
836 double conehmaxNearP_goodTrk = -999.0;
846 conehnLayers_maxNearP,
847 conehtrkQual_maxNearP,
848 conehmaxNearP_goodTrk,
853 v_hmaxNearP_goodTrk.push_back(hmaxNearP_goodTrk);
854 v_hmaxNearP.push_back(hmaxNearP);
855 v_hnNearTRKs.push_back(hnNearTRKs);
856 v_hnLayers_maxNearP.push_back(hnLayers_maxNearP);
857 v_htrkQual_maxNearP.push_back(htrkQual_maxNearP);
859 v_cone_hmaxNearP_goodTrk.push_back(conehmaxNearP_goodTrk);
860 v_cone_hmaxNearP.push_back(conehmaxNearP);
861 v_cone_hnNearTRKs.push_back(conehnNearTRKs);
862 v_cone_hnLayers_maxNearP.push_back(conehnLayers_maxNearP);
863 v_cone_htrkQual_maxNearP.push_back(conehtrkQual_maxNearP);
866 double h3x3 = -999.0, h5x5 = -999.0;
867 double hsim3x3 = -999.0, hsim5x5 = -999.0, trkHcalEne = -999.0;
868 std::map<std::string, double> hsimInfo3x3, hsimInfo5x5;
869 double distFromHotCell_h3x3 = -99.;
870 int ietaFromHotCell_h3x3 = -99;
871 int iphiFromHotCell_h3x3 = -99;
872 double distFromHotCell_h5x5 = -99.;
873 int ietaFromHotCell_h5x5 = -99;
874 int iphiFromHotCell_h5x5 = -99;
879 int nRH_h3x3(0), nRH_h5x5(0);
894 std::vector<int> v_RH_h3x3_ieta;
895 std::vector<int> v_RH_h3x3_iphi;
896 std::vector<double> v_RH_h3x3_ene;
897 std::vector<int> v_RH_h5x5_ieta;
898 std::vector<int> v_RH_h5x5_iphi;
899 std::vector<double> v_RH_h5x5_ene;
929 std::vector<int> multiplicity3x3;
930 std::vector<int> multiplicity5x5;
936 iEvent, theHBHETopology, ClosestCell, pcalohh, SimTk, SimVtx, pTrack, *associate, 1, 1, multiplicity3x3);
938 iEvent, theHBHETopology, ClosestCell, pcalohh, SimTk, SimVtx, pTrack, *associate, 2, 2, multiplicity5x5);
945 std::vector<double> v_hsimInfoConeMatched;
946 std::vector<double> v_hsimInfoConeRest;
947 std::vector<double> v_hsimInfoConePhoton;
948 std::vector<double> v_hsimInfoConeNeutHad;
949 std::vector<double> v_hsimInfoConeCharHad;
950 std::vector<double> v_hsimInfoConePdgMatched;
951 std::vector<double> v_hsimInfoConeTotal;
953 std::vector<int> v_hsimInfoConeNMatched;
954 std::vector<int> v_hsimInfoConeNTotal;
955 std::vector<int> v_hsimInfoConeNNeutHad;
956 std::vector<int> v_hsimInfoConeNCharHad;
957 std::vector<int> v_hsimInfoConeNPhoton;
958 std::vector<int> v_hsimInfoConeNRest;
960 std::vector<double> v_hsimCone;
961 std::vector<double> v_hCone;
963 std::vector<int> v_nRecHitsCone;
964 std::vector<int> v_nSimHitsCone;
966 std::vector<double> v_distFromHotCell;
967 std::vector<int> v_ietaFromHotCell;
968 std::vector<int> v_iphiFromHotCell;
971 std::vector<int> v_RH_r26_ieta;
972 std::vector<int> v_RH_r26_iphi;
973 std::vector<double> v_RH_r26_ene;
974 std::vector<int> v_RH_r44_ieta;
975 std::vector<int> v_RH_r44_iphi;
976 std::vector<double> v_RH_r44_ene;
978 for (
int i = 0;
i < a_size;
i++) {
979 std::map<std::string, double> hsimInfoCone;
980 double hsimCone = -999.0, hCone = -999.0;
981 double distFromHotCell = -99.0;
982 int ietaFromHotCell = -99;
983 int iphiFromHotCell = -99;
984 int ietaHotCell = -99;
985 int iphiHotCell = -99;
986 int nRecHitsCone = -999;
987 int nSimHitsCone = -999;
989 std::vector<int> multiplicityCone;
990 std::vector<DetId> coneRecHitDetIds;
992 hsimCone =
spr::eCone_hcal(geo, pcalohh, hpoint1, point1, a_coneR[
i], trackMomAtHcal, nSimHitsCone);
996 bool makeHitmaps =
false;
997 if (a_coneR[
i] == 26.23 && makeHitmaps) {
1013 }
else if (a_coneR[
i] == 43.72 && makeHitmaps) {
1044 if (ietaHotCell != 99) {
1045 ietaFromHotCell = ietaHotCell - ClosestCell_HcalDetId.ieta();
1046 iphiFromHotCell = iphiHotCell - ClosestCell_HcalDetId.iphi();
1051 hsimInfoCone = spr::eHCALSimInfoCone(
iEvent,
1072 v_hsimInfoConeMatched.push_back(hsimInfoCone[
"eMatched"]);
1073 v_hsimInfoConeRest.push_back(hsimInfoCone[
"eRest"]);
1074 v_hsimInfoConePhoton.push_back(hsimInfoCone[
"eGamma"]);
1075 v_hsimInfoConeNeutHad.push_back(hsimInfoCone[
"eNeutralHad"]);
1076 v_hsimInfoConeCharHad.push_back(hsimInfoCone[
"eChargedHad"]);
1077 v_hsimInfoConePdgMatched.push_back(hsimInfoCone[
"pdgMatched"]);
1078 v_hsimInfoConeTotal.push_back(hsimInfoCone[
"eTotal"]);
1080 v_hsimInfoConeNMatched.push_back(multiplicityCone.at(0));
1082 v_hsimInfoConeNTotal.push_back(multiplicityCone.at(1));
1083 v_hsimInfoConeNNeutHad.push_back(multiplicityCone.at(2));
1084 v_hsimInfoConeNCharHad.push_back(multiplicityCone.at(3));
1085 v_hsimInfoConeNPhoton.push_back(multiplicityCone.at(4));
1086 v_hsimInfoConeNRest.push_back(multiplicityCone.at(5));
1088 v_hsimCone.push_back(hsimCone);
1089 v_nSimHitsCone.push_back(nSimHitsCone);
1091 v_hCone.push_back(hCone);
1092 v_nRecHitsCone.push_back(nRecHitsCone);
1094 v_distFromHotCell.push_back(distFromHotCell);
1095 v_ietaFromHotCell.push_back(ietaFromHotCell);
1096 v_iphiFromHotCell.push_back(iphiFromHotCell);
1257 genPartPBins_ = {{0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0,
1258 12.0, 15.0, 20.0, 25.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 100.0}};
1262 t_v_hnNearTRKs = std::make_unique<std::vector<std::vector<int> > >();
1266 t_v_hmaxNearP = std::make_unique<std::vector<std::vector<double> > >();
1278 t_trkP = std::make_unique<std::vector<double> >();
1279 t_trkPt = std::make_unique<std::vector<double> >();
1280 t_trkEta = std::make_unique<std::vector<double> >();
1281 t_trkPhi = std::make_unique<std::vector<double> >();
1283 t_e3x3 = std::make_unique<std::vector<double> >();
1284 t_v_eDR = std::make_unique<std::vector<std::vector<double> > >();
1285 t_v_eMipDR = std::make_unique<std::vector<std::vector<double> > >();
1287 t_h3x3 = std::make_unique<std::vector<double> >();
1288 t_h5x5 = std::make_unique<std::vector<double> >();
1290 t_nRH_h3x3 = std::make_unique<std::vector<double> >();
1291 t_nRH_h5x5 = std::make_unique<std::vector<double> >();
1294 t_simP = std::make_unique<std::vector<double> >();
1295 t_hsim3x3 = std::make_unique<std::vector<double> >();
1296 t_hsim5x5 = std::make_unique<std::vector<double> >();
1327 t_trkHcalEne = std::make_unique<std::vector<double> >();
1328 t_trkEcalEne = std::make_unique<std::vector<double> >();
1354 t_v_hsimCone = std::make_unique<std::vector<std::vector<double> > >();
1357 t_v_hCone = std::make_unique<std::vector<std::vector<double> > >();
1367 t_v_RH_h3x3_ene = std::make_unique<std::vector<std::vector<double> > >();
1370 t_v_RH_h5x5_ene = std::make_unique<std::vector<std::vector<double> > >();
1371 t_v_RH_r26_ieta = std::make_unique<std::vector<std::vector<int> > >();
1372 t_v_RH_r26_iphi = std::make_unique<std::vector<std::vector<int> > >();
1373 t_v_RH_r26_ene = std::make_unique<std::vector<std::vector<double> > >();
1374 t_v_RH_r44_ieta = std::make_unique<std::vector<std::vector<int> > >();
1375 t_v_RH_r44_iphi = std::make_unique<std::vector<std::vector<int> > >();
1376 t_v_RH_r44_ene = std::make_unique<std::vector<std::vector<double> > >();
1378 t_v_hlTriggers = std::make_unique<std::vector<std::vector<int> > >();
1380 t_hltHE = std::make_unique<std::vector<int> >();
1381 t_hltHB = std::make_unique<std::vector<int> >();
1383 t_hltJet30 = std::make_unique<std::vector<int> >();
1384 t_hltJet50 = std::make_unique<std::vector<int> >();
1385 t_hltJet80 = std::make_unique<std::vector<int> >();
1386 t_hltJet110 = std::make_unique<std::vector<int> >();
1387 t_hltJet140 = std::make_unique<std::vector<int> >();
1388 t_hltJet180 = std::make_unique<std::vector<int> >();
1397 t_irun = std::make_unique<std::vector<unsigned int> >();
1398 t_ievt = std::make_unique<std::vector<unsigned int> >();
1399 t_ilum = std::make_unique<std::vector<unsigned int> >();
1547 h_RawPt = fs->
make<TH1F>(
"hRawPt",
"hRawPt", 100, 0.0, 100.0);
1548 h_RawP = fs->
make<TH1F>(
"hRawP",
"hRawP", 100, 0.0, 100.0);
1549 h_RawEta = fs->
make<TH1F>(
"hRawEta",
"hRawEta", 15, 0.0, 3.0);
1550 h_RawPhi = fs->
make<TH1F>(
"hRawPhi",
"hRawPhi", 100, -3.2, 3.2);
1552 ntp_ = fs->
make<TTree>(
"ntp",
"ntp");
1555 ntp_->Branch(
"nEVT", &
nEVT,
"nEVT/I");
1559 ntp_->Branch(
"nTRK", &
nTRK,
"nTRK/I");
1583 ntp_->Branch(
"trkP",
"std::vector<double>", &
t_trkP);
1584 ntp_->Branch(
"trkPt",
"std::vector<double>", &
t_trkPt);
1585 ntp_->Branch(
"trkEta",
"std::vector<double>", &
t_trkEta);
1586 ntp_->Branch(
"trkPhi",
"std::vector<double>", &
t_trkPhi);
1587 ntp_->Branch(
"e3x3",
"std::vector<double>", &
t_e3x3);
1589 ntp_->Branch(
"e3x3",
"std::vector<double>", &
t_e3x3);
1590 ntp_->Branch(
"v_eDR",
"std::vector<std::vector<double> >", &
t_v_eDR);
1591 ntp_->Branch(
"v_eMipDR",
"std::vector<std::vector<double> >", &
t_v_eMipDR);
1593 ntp_->Branch(
"h3x3",
"std::vector<double>", &
t_h3x3);
1594 ntp_->Branch(
"h5x5",
"std::vector<double>", &
t_h5x5);
1599 ntp_->Branch(
"simP",
"std::vector<double>", &
t_simP);
1659 ntp_->Branch(
"v_hsimCone",
"std::vector<std::vector<double> >", &
t_v_hsimCone);
1662 ntp_->Branch(
"v_hCone",
"std::vector<std::vector<double> >", &
t_v_hCone);
1684 ntp_->Branch(
"v_hltHB",
"std::vector<int>", &
t_hltHB);
1685 ntp_->Branch(
"v_hltHE",
"std::vector<int>", &
t_hltHE);
1701 ntp_->Branch(
"irun",
"std::vector<unsigned int>", &
t_irun);
1702 ntp_->Branch(
"ievt",
"std::vector<unsigned int>", &
t_ievt);
1703 ntp_->Branch(
"ilum",
"std::vector<unsigned int>", &
t_ilum);
1711 <<
" TrackMmentum " << pTrack->
momentum() <<
" (pt,eta,phi)(" << pTrack->
pt() <<
","
1712 << pTrack->
eta() <<
"," << pTrack->
phi() <<
")"
1713 <<
" p " << pTrack->
p() <<
"\n"
1716 << pTrack->
d0() <<
"\n"
1719 <<
" " << pTrack->
quality(trackQuality_);