112 iEvent.
getByLabel(
"generalTracks", trkCollection);
115 iEvent.
getByLabel(
"offlinePrimaryVertices",recVtxs);
119 iEvent.
getByLabel(
"offlineBeamSpot", beamSpotH);
122 if (recVtxs->size()>0 && !((*recVtxs)[0].isFake())) {
123 leadPV =
math::XYZPoint( (*recVtxs)[0].
x(),(*recVtxs)[0].
y(), (*recVtxs)[0].
z() );
124 }
else if (beamSpotH.
isValid()) {
125 leadPV = beamSpotH->position();
129 std::cout <<
"Primary Vertex " << leadPV;
130 if (beamSpotH.
isValid())
std::cout <<
" Beam Spot " << beamSpotH->position();
134 std::vector<spr::propagatedTrackDirection> trkCaloDirections;
136 std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
140 iEvent.
getByLabel(
"ecalRecHit",
"EcalRecHitsEB",barrelRecHitsHandle);
141 iEvent.
getByLabel(
"ecalRecHit",
"EcalRecHitsEE",endcapRecHitsHandle);
149 edm::SimTrackContainer::const_iterator simTrkItr;
151 edm::SimVertexContainer::const_iterator vtxItr = SimVtx->begin();
164 iEvent.
getByLabel(
"g4SimHits",
"EcalHitsEB", pcaloeb);
165 iEvent.
getByLabel(
"g4SimHits",
"EcalHitsEE", pcaloee);
166 iEvent.
getByLabel(
"g4SimHits",
"HcalHits", pcalohh);
170 unsigned int nTracks=0;
171 for (trkDetItr = trkCaloDirections.begin(),nTracks=0; trkDetItr != trkCaloDirections.end(); trkDetItr++,nTracks++){
172 const reco::Track* pTrack = &(*(trkDetItr->trkItr));
174 int nRH_eMipDR=0, nRH_eDR=0, nNearTRKs=0, nRecHitsCone=-99;
175 double distFromHotCell=-99.0, distFromHotCell2=-99.0;
176 int ietaHotCell=-99, iphiHotCell=-99;
177 int ietaHotCell2=-99, iphiHotCell2=-99;
178 GlobalPoint gposHotCell(0.,0.,0.), gposHotCell2(0.,0.,0.);
179 std::vector<DetId> coneRecHitDetIds, coneRecHitDetIds2;
180 std::pair<double, bool> e11x11_20SigP, e15x15_20SigP;
181 double hCone =
spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL,
182 a_coneR, trkDetItr->directionHCAL, nRecHitsCone,
183 coneRecHitDetIds, distFromHotCell,
184 ietaHotCell, iphiHotCell, gposHotCell);
185 double eHCALDR =
spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL, trkDetItr->pointECAL,
186 a_charIsoR, trkDetItr->directionHCAL, nRecHitsCone,
187 coneRecHitDetIds2, distFromHotCell2,
188 ietaHotCell2, iphiHotCell2, gposHotCell2);
192 double eMipDR =
spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
193 trkDetItr->pointHCAL, trkDetItr->pointECAL,
194 a_mipR, trkDetItr->directionECAL, nRH_eMipDR);
195 double eECALDR =
spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
196 trkDetItr->pointHCAL, trkDetItr->pointECAL,
197 a_neutIsoR, trkDetItr->directionECAL, nRH_eDR);
204 e11x11_20SigP =
spr::eECALmatrix(trkDetItr->detIdECAL,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.
product(),5,5, 0.060, 0.300,
tMinE_,
tMaxE_);
205 e15x15_20SigP =
spr::eECALmatrix(trkDetItr->detIdECAL,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.
product(),7,7, 0.060, 0.300,
tMinE_,
tMaxE_);
223 std::cout <<
"Track p " << pTrack->
p() <<
" pt " << pTrack->
pt()
224 <<
" eta " << pTrack->
momentum().eta() <<
" phi "
225 << pTrack->
momentum().phi() <<
" ieta/iphi ("
226 << closestCell.
ieta() <<
", " << closestCell.
iphi()
227 <<
") Energy in cone " << hCone <<
" Charge Isolation "
228 << conehmaxNearP <<
" eMIP " << eMipDR
229 <<
" Neutral isolation (ECAL) " << eECALDR-eMipDR
230 <<
" (ECAL NxN) " << e15x15_20SigP.first-e11x11_20SigP.first
231 <<
" (HCAL) " << eHCALDR-hCone << std::endl;
237 std::map<std::string, double> hsimInfo;
238 std::vector<int> multiplicity;
241 trkDetItr->directionHCAL, nSimHits);
242 hsimInfo = spr::eHCALSimInfoCone(iEvent, pcalohh, SimTk, SimVtx,
243 pTrack, *associate, geo,
244 trkDetItr->pointHCAL,
246 trkDetItr->directionHCAL,
264 t_hsim ->push_back(hsim );
268 std::cout <<
"Matched (E) " << hsimInfo[
"eMatched"] <<
" (N) "
269 << multiplicity.at(0) <<
" Rest (E) " << hsimInfo[
"eRest"]
270 <<
" (N) " << multiplicity.at(1) <<
" Gamma (E) "
271 << hsimInfo[
"eGamma"] <<
" (N) " << multiplicity.at(2)
272 <<
" Neutral Had (E) " << hsimInfo[
"eNeutralHad"]
273 <<
" (N) " << multiplicity.at(3) <<
" Charged Had (E) "
274 << hsimInfo[
"eChargedHad"] <<
" (N) " << multiplicity.at(4)
275 <<
" Total (E) " << hsimInfo[
"eTotal"] <<
" (N) "
276 << multiplicity.at(5) <<
" PDG " << hsimInfo[
"pdgMatched"]
277 <<
" Total E " << hsim <<
" NHit " << nSimHits <<std::endl;
284 if (associate)
delete associate;
double p() const
momentum vector magnitude
std::vector< double > * t_hsimInfoMatched
EventNumber_t event() const
std::vector< double > * t_eHCALDR
std::vector< int > * t_hsimInfoNRest
std::string theTrackQuality
std::vector< double > * t_hCone
std::vector< spr::propagatedTrackID > propagateCALO(edm::Handle< reco::TrackCollection > &trkCollection, const CaloGeometry *geo, const MagneticField *bField, std::string &theTrackQuality, bool debug=false)
const Vector & momentum() const
track momentum vector
std::vector< int > * t_hsimInfoNTotal
std::vector< int > * t_hsimInfoNCharHad
std::vector< int > * t_hsimInfoNMatched
std::vector< double > * t_hsim
int bunchCrossing() const
std::vector< double > * t_hsimInfoNeutHad
edm::LuminosityBlockNumber_t luminosityBlock() const
std::vector< int > * t_hsimInfoNPhoton
double chargeIsolationCone(unsigned int trkIndex, std::vector< spr::propagatedTrackDirection > &trkDirs, double dR, int &nNearTRKs, bool debug=false)
double eCone_ecal(const CaloGeometry *geo, edm::Handle< T > &barrelhits, edm::Handle< T > &endcaphits, const GlobalPoint &hpoint1, const GlobalPoint &point1, double dR, const GlobalVector &trackMom, int &nRecHits)
std::vector< double > * t_e11x11_20Sig
std::vector< int > * t_hsimInfoNNeutHad
bool goodTrack(const reco::Track *pTrack, math::XYZPoint leadPV, trackSelectionParameters parameters, bool debug=false)
double pt() const
track transverse momentum
std::vector< double > * t_trackEta
int ieta() const
get the cell ieta
std::vector< double > * t_eECALDR
spr::trackSelectionParameters selectionParameters
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
int iphi() const
get the cell iphi
std::vector< double > * t_trackPt
const MagneticField * bField
std::vector< double > * t_trackHcalEta
std::vector< double > * t_trackP
XYZPointD XYZPoint
point in space with cartesian internal representation
std::vector< double > * t_trackHcalPhi
T const * product() const
T const * product() const
std::vector< double > * t_hsimInfoPhoton
std::vector< double > * t_hsimInfoTotal
std::vector< double > * t_hsimInfoCharHad
std::vector< double > * t_trackPhi
double eECALmatrix(const DetId &detId, edm::Handle< T > &hitsEB, edm::Handle< T > &hitsEE, const CaloGeometry *geo, const CaloTopology *caloTopology, const EcalSeverityLevelAlgo *sevlv, int ieta, int iphi, double ebThr=-100, double eeThr=-100, double tMin=-500, double tMax=500, bool debug=false)
double eCone_hcal(const CaloGeometry *geo, edm::Handle< T > &hits, const GlobalPoint &hpoint1, const GlobalPoint &point1, double dR, const GlobalVector &trackMom, int &nRecHits)
std::vector< double > * t_hsimInfoPdgMatched
std::vector< int > * t_nSimHits
std::vector< double > * t_hsimInfoRest
std::vector< double > * t_conehmaxNearP
std::vector< double > * t_e15x15_20Sig
std::vector< double > * t_eMipDR