71 SiPixelPhase1TrackEfficiency::SiPixelPhase1TrackEfficiency(
const edm::ParameterSet& iConfig) :
75 tracksToken_ = consumes<reco::TrackCollection>(iConfig.getParameter<
edm::InputTag>(
"tracks"));
76 vtxToken_ = consumes<reco::VertexCollection>(iConfig.getParameter<
edm::InputTag>(
"primaryvertices"));
77 applyVertexCut_=iConfig.getUntrackedParameter<
bool>(
"VertexCut",
true);
78 trajTrackCollectionToken_ = consumes<TrajTrackAssociationCollection>(iConfig.getParameter<
edm::InputTag>(
"trajectoryInput")); clustersToken_=consumes<edmNew::DetSetVector<SiPixelCluster> >(iConfig.getParameter<
edm::InputTag>(
"clusters"));
82 if( !checktrigger(iEvent,iSetup,DCS) )
return;
96 trackerTopology_ = trackerTopologyHandle.
product();
101 std::unique_ptr<Propagator> propagatorUniquePtr(propagatorHandle.
product() ->
clone());
102 trackerPropagator_ = propagatorUniquePtr.get();
108 chi2MeasurementEstimator_ = chi2MeasurementEstimatorHandle.
product();
118 if (!vertices.isValid())
return;
120 if (applyVertexCut_ && vertices->empty())
return;
132 iEvent.
getByToken(trajTrackCollectionToken_, trajTrackCollectionHandle);
133 if (!trajTrackCollectionHandle.
isValid())
return;
137 iEvent.
getByToken(clustersToken_, siPixelClusters);
138 if(!siPixelClusters.
isValid())
return;
143 if(!cpEstimator.
isValid())
return;
151 int TRACK_QUALITY_HIGH_PURITY_BIT = 2;
152 int TRACK_QUALITY_HIGH_PURITY_MASK = 1 << TRACK_QUALITY_HIGH_PURITY_BIT;
155 float TRACK_PT_CUT_VAL = 1.0f;
158 int TRACK_NSTRIP_CUT_VAL = 10;
161 std::array<float, 4> TRACK_D0_CUT_BARREL_VAL = {{0.01f, 0.02f, 0.02f, 0.02f}};
162 float TRACK_D0_CUT_FORWARD_VAL = 0.05f;
165 float TRACK_DZ_CUT_BARREL_VAL = 0.01f;
166 float TRACK_DZ_CUT_FORWARD_VAL = 0.5f;
168 bool isBpixtrack =
false, isFpixtrack =
false;
177 bool passcuts =
true;
178 bool passcuts_hit =
true;
182 bool valid_layerFrom =
false;
189 std::vector<TrajectoryMeasurement> expTrajMeasurements;
190 std::vector<std::pair<int,bool[3]>> eff_pxb1_vector;
193 for(
const auto &pair: *trajTrackCollectionHandle) {
197 expTrajMeasurements.clear();
198 eff_pxb1_vector.clear();
200 if (applyVertexCut_ && (track->pt() < 0.75 ||
std::abs( track->dxy(vertices->at(0).position()) ) > 5*track->dxyError()))
continue;
202 isBpixtrack =
false, isFpixtrack =
false;
219 auto hb = track->recHitsBegin();
220 for(
unsigned int h=0;
h<track->recHitsSize();
h++){
223 if(!
hit->isValid())
continue;
226 uint32_t subdetid = (
id.subdetId());
232 if(trackerTopology_ -> pxbLayer(
id) == 1) nBpixL1Hits++;
233 if(trackerTopology_ -> pxbLayer(
id) == 2) nBpixL2Hits++;
234 if(trackerTopology_ -> pxbLayer(
id) == 3) nBpixL3Hits++;
235 if(trackerTopology_ -> pxbLayer(
id) == 4) nBpixL4Hits++;
240 if(trackerTopology_ -> pxfDisk(
id) == 1) nFpixD1Hits++;
241 if(trackerTopology_ -> pxfDisk(
id) == 2) nFpixD2Hits++;
242 if(trackerTopology_ -> pxfDisk(
id) == 3) nFpixD3Hits++;
256 if (!isBpixtrack && !isFpixtrack)
continue;
259 if(!((track->qualityMask() & TRACK_QUALITY_HIGH_PURITY_MASK) >> TRACK_QUALITY_HIGH_PURITY_BIT)) passcuts =
false;
262 if(!(TRACK_PT_CUT_VAL < track->
pt())) passcuts =
false;
265 if(!(TRACK_NSTRIP_CUT_VAL < nStripHits)) passcuts =
false;
269 for(
unsigned int h=0;
h<track->recHitsSize();
h++){
275 uint32_t subdetid = (
id.subdetId());
286 {
if(!((
std::abs( track->dxy(vertices->at(0).position()) ) * -1.0) < TRACK_D0_CUT_BARREL_VAL[trackerTopology_ -> pxbLayer(
id) -1])) passcuts_hit =
false;}
288 {
if(!((
std::abs( track->dxy(vertices->at(0).position()) ) * -1.0) < TRACK_D0_CUT_FORWARD_VAL)) passcuts_hit =
false;}
293 {
if(!(
std::abs( track->dz(vertices->at(0).position())) < TRACK_DZ_CUT_BARREL_VAL)) passcuts_hit =
false;}
295 {
if(!(
std::abs( track->dz(vertices->at(0).position())) < TRACK_DZ_CUT_FORWARD_VAL)) passcuts_hit =
false;}
302 if(trackerTopology_ -> pxbLayer(
id) == 1)
if(!(
303 (nBpixL2Hits > 0 && nBpixL3Hits > 0 && nBpixL4Hits > 0) ||
304 (nBpixL2Hits > 0 && nBpixL3Hits > 0 && nFpixD1Hits > 0) ||
305 (nBpixL2Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0) ||
306 (nFpixD1Hits > 0 && nFpixD2Hits > 0 && nFpixD3Hits > 0))) passcuts_hit =
false;
307 if(trackerTopology_ -> pxbLayer(
id) == 2)
if(!(
308 (nBpixL1Hits > 0 && nBpixL3Hits > 0 && nBpixL4Hits > 0) ||
309 (nBpixL1Hits > 0 && nBpixL3Hits > 0 && nFpixD1Hits > 0) ||
310 (nBpixL1Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0))) passcuts_hit =
false;
311 if(trackerTopology_ -> pxbLayer(
id) == 3)
if(!(
312 (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nBpixL4Hits > 0) ||
313 (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nFpixD1Hits > 0))) passcuts_hit =
false;
314 if(trackerTopology_ -> pxbLayer(
id) == 4)
if(!(
315 (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nBpixL3Hits > 0))) passcuts_hit =
false;
319 if(trackerTopology_ -> pxfDisk(
id) == 1)
if(!(
320 (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nBpixL3Hits > 0) ||
321 (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nFpixD2Hits > 0) ||
322 (nBpixL1Hits > 0 && nFpixD2Hits > 0 && nFpixD3Hits > 0))) passcuts_hit =
false;
323 if(trackerTopology_ -> pxfDisk(
id) == 2)
if(!(
324 (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nFpixD1Hits > 0) ||
325 (nBpixL1Hits > 0 && nFpixD1Hits > 0 && nFpixD3Hits > 0))) passcuts_hit =
false;
326 if(trackerTopology_ -> pxfDisk(
id) == 3)
if(!(
327 (nBpixL1Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0))) passcuts_hit =
false;
351 if (passcuts_hit ==
true && passcuts){
356 histo[VALID].fill(
id, &iEvent);
357 histo[EFFICIENCY].fill(1,
id, &iEvent);
360 histo[MISSING].fill(
id, &iEvent);
361 histo[EFFICIENCY].fill(0,
id, &iEvent);
364 histo[INACTIVE].fill(
id, &iEvent);
375 valid_layerFrom =
false;
378 for (
const auto &tm : traj->measurements()) {
379 if (tm.recHit().get() && tm.recHitR().isValid()) {
380 DetId where = tm.recHitR().geographicalId();
384 int source_layer = trackerTopology_ -> pxbLayer(where);
385 if (source_layer==2){
386 if (tm.updatedState().isValid()) {tsosPXB2 = tm.updatedState(); valid_layerFrom=
true;}
391 int source_layer = trackerTopology_ -> pxfDisk(where);
392 if (source_layer==1){
393 if (tm.updatedState().isValid()) {tsosPXB2 = tm.updatedState(); valid_layerFrom=
true;}
402 if (!valid_layerFrom)
continue;
403 if (!tsosPXB2.
isValid())
continue;
408 expTrajMeasurements = theLayerMeasurements_->
measurements(*pxbLayer1_, tsosPXB2, *trackerPropagator_, *chi2MeasurementEstimator_);
409 auto compDets = pxbLayer1_->compatibleDets(tsosPXB2, *trackerPropagator_, *chi2MeasurementEstimator_);
410 std::pair<int,bool[3]> eff_map;
412 bool missing =
false;
417 for(
uint p=0;
p<expTrajMeasurements.size();
p++){
420 const auto& pxb1Hit = pxb1TM.recHit();
421 int detidHit= pxb1Hit->geographicalId();
422 if (detidHit==0)
continue;
434 int row = (
int) mp.
x() % 80;
440 if (!((col < (centercol + 10)) && (col > (centercol - 10)) && (row < (centerrow + 10)) && (row > (centerrow - 10 )))) passcuts_hit =
false;
444 for (
const auto & detAndState : compDets) {
446 const auto & pXb1_lpos = detAndState.second.localPosition();
447 int detid = detAndState.first->geographicalId().rawId();
449 DetId detId(iter_cl->id());
450 float minD[2]; minD[0]=minD[1]=10000.;
451 if(detId.rawId()!=detAndState.first->geographicalId().rawId())
continue;
452 if(pxb1Hit->geographicalId().rawId()!=detAndState.first->geographicalId().rawId())
continue;
456 for( ; itCluster!=iter_cl->
end(); ++itCluster){
458 LocalPoint lp(itCluster->x(), itCluster->y(), 0.);
460 lp=std::get<0>(params);
462 float Xdist =
abs(lp.
x()-pXb1_lpos.x());
463 float Ydist =
abs(lp.
y()-pXb1_lpos.y());
472 if ((minD[0] < 0.02) && (minD[1] < 0.02)) {
487 if(!((
std::abs( track->dxy(vertices->at(0).position()) ) * -1.0) < TRACK_D0_CUT_BARREL_VAL[trackerTopology_ -> pxbLayer(detid) -1])) passcuts_hit =
false;
489 if(!(
std::abs( track->dz(vertices->at(0).position())) < TRACK_DZ_CUT_BARREL_VAL)) passcuts_hit =
false;
491 if(!((nBpixL2Hits > 0 && nBpixL3Hits > 0 && nBpixL4Hits > 0) || (nBpixL2Hits > 0 && nBpixL3Hits > 0 && nFpixD1Hits > 0) ||
492 (nBpixL2Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0) || (nFpixD1Hits > 0 && nFpixD2Hits > 0 && nFpixD3Hits > 0))) passcuts_hit =
false;
493 bool found_det =
false;
495 if (passcuts && passcuts_hit){
496 for (
unsigned int i_eff=0; i_eff<eff_pxb1_vector.size(); i_eff++){
498 if (eff_pxb1_vector[i_eff].
first==detid){
501 if (eff_pxb1_vector[i_eff].
second[0]==
false && valid==
true){
502 eff_pxb1_vector[i_eff].second[0]=valid;
503 eff_pxb1_vector[i_eff].second[1]=missing;
512 eff_map.second[0]=valid;
513 eff_map.second[1]=missing;
514 eff_pxb1_vector.push_back(eff_map);
524 for(
uint p=0;
p<expTrajMeasurements.size();
p++){
527 const auto& pxb1Hit = pxb1TM.recHit();
529 int detid= pxb1Hit->geographicalId();
530 bool found_det =
false;
532 if (passcuts && passcuts_hit){
533 for (
unsigned int i_eff=0; i_eff<eff_pxb1_vector.size(); i_eff++){
535 if (eff_pxb1_vector[i_eff].
first==detid){
538 if (eff_pxb1_vector[i_eff].
second[0]==
false && valid==
true){
539 eff_pxb1_vector[i_eff].second[2]=inactive;
549 eff_map.second[2]=inactive;
550 eff_pxb1_vector.push_back(eff_map);
560 if (eff_pxb1_vector.size() == 1) {
564 if (eff_pxb1_vector[0].
second[0]) {
565 histo[VALID].fill(eff_pxb1_vector[0].
first, &iEvent);
566 histo[EFFICIENCY].fill(1, eff_pxb1_vector[0].first, &iEvent);
569 if (eff_pxb1_vector[0].second[1]) {
570 histo[MISSING].fill(eff_pxb1_vector[0].first, &iEvent);
571 histo[EFFICIENCY].fill(0, eff_pxb1_vector[0].first, &iEvent);
575 if (eff_pxb1_vector[0].second[2]) {
576 histo[INACTIVE].fill(eff_pxb1_vector[0].first, &iEvent);
585 histo[VALID ].executePerEventHarvesting(&iEvent);
586 histo[MISSING ].executePerEventHarvesting(&iEvent);
587 histo[INACTIVE].executePerEventHarvesting(&iEvent);
T getParameter(std::string const &) const
std::vector< TrajectoryMeasurement > measurements(const DetLayer &layer, const TrajectoryStateOnSurface &startingState, const Propagator &prop, const MeasurementEstimator &est) const
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
static const char tracker_[]
FWCore Framework interface EventSetupRecordImplementation h
Helper function to determine trigger accepts.
virtual example_stream void analyze(const edm::Event &, const edm::EventSetup &) override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
#define DEFINE_FWK_MODULE(type)
data_type const * const_iterator
key_type key() const
Accessor for product key.
U second(std::pair< T, U > const &p)
const TrackerGeomDet * idToDetUnit(DetId) const override
Return the pointer to the GeomDetUnit corresponding to a given DetId.
virtual MeasurementPoint measurementPosition(const LocalPoint &) const =0
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector's numbering enum) ...
Abs< T >::type abs(const T &t)
std::vector< BarrelDetLayer const * > const & pixelBarrelLayers() const
void analyze(edm::Event const &e, edm::EventSetup const &eSetup) override=0
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
TEveGeoShape * clone(const TEveElement *element, TEveElement *parent)
const TrackerGeomDet * idToDet(DetId) const override
LocalPoint localPosition() const final
std::tuple< LocalPoint, LocalError, SiPixelRecHitQuality::QualWordType > ReturnType
T const * product() const
const GeometricSearchTracker * geometricSearchTracker() const