File indexing completed on 2025-02-05 23:50:58
0001 #include "Alignment/CommonAlignmentProducer/interface/AlignmentTrackSelector.h"
0002
0003 #include "FWCore/Framework/interface/ESHandle.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006
0007 #include "FWCore/Framework/interface/Event.h"
0008
0009 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0010 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0011 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0012 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0013 #include "DataFormats/TrackerRecHit2D/interface/ProjectedSiStripRecHit2D.h"
0014 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0015 #include "DataFormats/TrackerRecHit2D/interface/SiTrackerMultiRecHit.h"
0016 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0017 #include "DataFormats/Alignment/interface/AlignmentClusterFlag.h"
0018 #include "DataFormats/Alignment/interface/AliClusterValueMap.h"
0019 #include "DataFormats/DetId/interface/DetId.h"
0020 #include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
0021 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0022
0023 #include <cmath>
0024
0025 const int kBPIX = PixelSubdetector::PixelBarrel;
0026 const int kFPIX = PixelSubdetector::PixelEndcap;
0027
0028
0029
0030 AlignmentTrackSelector::AlignmentTrackSelector(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC)
0031 : tTopoToken_(iC.esConsumes()),
0032 applyBasicCuts_(cfg.getParameter<bool>("applyBasicCuts")),
0033 applyNHighestPt_(cfg.getParameter<bool>("applyNHighestPt")),
0034 applyMultiplicityFilter_(cfg.getParameter<bool>("applyMultiplicityFilter")),
0035 seedOnlyFromAbove_(cfg.getParameter<int>("seedOnlyFrom")),
0036 applyIsolation_(cfg.getParameter<bool>("applyIsolationCut")),
0037 chargeCheck_(cfg.getParameter<bool>("applyChargeCheck")),
0038 nHighestPt_(cfg.getParameter<int>("nHighestPt")),
0039 minMultiplicity_(cfg.getParameter<int>("minMultiplicity")),
0040 maxMultiplicity_(cfg.getParameter<int>("maxMultiplicity")),
0041 multiplicityOnInput_(cfg.getParameter<bool>("multiplicityOnInput")),
0042 ptMin_(cfg.getParameter<double>("ptMin")),
0043 ptMax_(cfg.getParameter<double>("ptMax")),
0044 pMin_(cfg.getParameter<double>("pMin")),
0045 pMax_(cfg.getParameter<double>("pMax")),
0046 etaMin_(cfg.getParameter<double>("etaMin")),
0047 etaMax_(cfg.getParameter<double>("etaMax")),
0048 phiMin_(cfg.getParameter<double>("phiMin")),
0049 phiMax_(cfg.getParameter<double>("phiMax")),
0050 nHitMin_(cfg.getParameter<double>("nHitMin")),
0051 nHitMax_(cfg.getParameter<double>("nHitMax")),
0052 chi2nMax_(cfg.getParameter<double>("chi2nMax")),
0053 d0Min_(cfg.getParameter<double>("d0Min")),
0054 d0Max_(cfg.getParameter<double>("d0Max")),
0055 dzMin_(cfg.getParameter<double>("dzMin")),
0056 dzMax_(cfg.getParameter<double>("dzMax")),
0057 theCharge_(cfg.getParameter<int>("theCharge")),
0058 minHitChargeStrip_(cfg.getParameter<double>("minHitChargeStrip")),
0059 minHitIsolation_(cfg.getParameter<double>("minHitIsolation")),
0060 countStereoHitAs2D_(cfg.getParameter<bool>("countStereoHitAs2D")),
0061 nHitMin2D_(cfg.getParameter<unsigned int>("nHitMin2D")),
0062
0063 minHitsinTIB_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTIB")),
0064 minHitsinTOB_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTOB")),
0065 minHitsinTID_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTID")),
0066 minHitsinTEC_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTEC")),
0067 minHitsinBPIX_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inBPIX")),
0068 minHitsinFPIX_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inFPIX")),
0069 minHitsinPIX_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inPIXEL")),
0070 minHitsinTIDplus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTIDplus")),
0071 minHitsinTIDminus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTIDminus")),
0072 minHitsinTECplus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTECplus")),
0073 minHitsinTECminus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inTECminus")),
0074 minHitsinFPIXplus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inFPIXplus")),
0075 minHitsinFPIXminus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inFPIXminus")),
0076 minHitsinENDCAP_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inENDCAP")),
0077 minHitsinENDCAPplus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inENDCAPplus")),
0078 minHitsinENDCAPminus_(cfg.getParameter<edm::ParameterSet>("minHitsPerSubDet").getParameter<int>("inENDCAPminus")),
0079 maxHitDiffEndcaps_(cfg.getParameter<double>("maxHitDiffEndcaps")),
0080 nLostHitMax_(cfg.getParameter<double>("nLostHitMax")),
0081 RorZofFirstHitMin_(cfg.getParameter<std::vector<double>>("RorZofFirstHitMin")),
0082 RorZofFirstHitMax_(cfg.getParameter<std::vector<double>>("RorZofFirstHitMax")),
0083 RorZofLastHitMin_(cfg.getParameter<std::vector<double>>("RorZofLastHitMin")),
0084 RorZofLastHitMax_(cfg.getParameter<std::vector<double>>("RorZofLastHitMax")),
0085 clusterValueMapTag_(cfg.getParameter<edm::InputTag>("hitPrescaleMapTag")),
0086 minPrescaledHits_(cfg.getParameter<int>("minPrescaledHits")),
0087 applyPrescaledHitsFilter_(!clusterValueMapTag_.encode().empty() && minPrescaledHits_ > 0) {
0088 if (applyIsolation_) {
0089 rphirecHitsToken_ = iC.consumes<SiStripRecHit2DCollection>(cfg.getParameter<edm::InputTag>("rphirecHits"));
0090 matchedrecHitsToken_ =
0091 iC.consumes<SiStripMatchedRecHit2DCollection>(cfg.getParameter<edm::InputTag>("matchedrecHits"));
0092 }
0093
0094 if (applyPrescaledHitsFilter_) {
0095 clusterValueMapToken_ = iC.consumes<AliClusterValueMap>(clusterValueMapTag_);
0096 }
0097
0098
0099 std::vector<std::string> trkQualityStrings(cfg.getParameter<std::vector<std::string>>("trackQualities"));
0100 std::string qualities;
0101 if (!trkQualityStrings.empty()) {
0102 applyTrkQualityCheck_ = true;
0103 for (unsigned int i = 0; i < trkQualityStrings.size(); ++i) {
0104 (qualities += trkQualityStrings[i]) += ", ";
0105 trkQualities_.push_back(reco::TrackBase::qualityByName(trkQualityStrings[i]));
0106 }
0107 } else
0108 applyTrkQualityCheck_ = false;
0109
0110 std::vector<std::string> trkIterStrings(cfg.getParameter<std::vector<std::string>>("iterativeTrackingSteps"));
0111 if (!trkIterStrings.empty()) {
0112 applyIterStepCheck_ = true;
0113 std::string tracksteps;
0114 for (unsigned int i = 0; i < trkIterStrings.size(); ++i) {
0115 (tracksteps += trkIterStrings[i]) += ", ";
0116 trkSteps_.push_back(reco::TrackBase::algoByName(trkIterStrings[i]));
0117 }
0118 } else
0119 applyIterStepCheck_ = false;
0120
0121 if (applyBasicCuts_) {
0122 edm::LogInfo("AlignmentTrackSelector")
0123 << "applying basic track cuts ..."
0124 << "\nptmin,ptmax: " << ptMin_ << "," << ptMax_ << "\npmin,pmax: " << pMin_ << "," << pMax_
0125 << "\netamin,etamax: " << etaMin_ << "," << etaMax_ << "\nphimin,phimax: " << phiMin_ << "," << phiMax_
0126 << "\nnhitmin,nhitmax: " << nHitMin_ << "," << nHitMax_ << "\nnlosthitmax: " << nLostHitMax_
0127 << "\nnhitmin2D: " << nHitMin2D_ << (countStereoHitAs2D_ ? "," : ", not")
0128 << " counting hits on SiStrip stereo modules as 2D"
0129 << "\nchi2nmax: " << chi2nMax_;
0130
0131 if (applyIsolation_)
0132 edm::LogInfo("AlignmentTrackSelector")
0133 << "only retain tracks isolated at least by " << minHitIsolation_ << " cm from other rechits";
0134
0135 if (chargeCheck_)
0136 edm::LogInfo("AlignmentTrackSelector")
0137 << "only retain hits with at least " << minHitChargeStrip_ << " ADC counts of total cluster charge";
0138
0139 edm::LogInfo("AlignmentTrackSelector")
0140 << "Minimum number of hits in TIB/TID/TOB/TEC/BPIX/FPIX/PIXEL = " << minHitsinTIB_ << "/" << minHitsinTID_
0141 << "/" << minHitsinTOB_ << "/" << minHitsinTEC_ << "/" << minHitsinBPIX_ << "/" << minHitsinFPIX_ << "/"
0142 << minHitsinPIX_;
0143
0144 edm::LogInfo("AlignmentTrackSelector")
0145 << "Minimum number of hits in TID+/TID-/TEC+/TEC-/FPIX+/FPIX- = " << minHitsinTIDplus_ << "/"
0146 << minHitsinTIDminus_ << "/" << minHitsinTECplus_ << "/" << minHitsinTECminus_ << "/" << minHitsinFPIXplus_
0147 << "/" << minHitsinFPIXminus_;
0148
0149 edm::LogInfo("AlignmentTrackSelector")
0150 << "Minimum number of hits in EndCap (TID+TEC)/EndCap+/EndCap- = " << minHitsinENDCAP_ << "/"
0151 << minHitsinENDCAPplus_ << "/" << minHitsinENDCAPminus_;
0152
0153 edm::LogInfo("AlignmentTrackSelector")
0154 << "Max value of |nHitsinENDCAPplus - nHitsinENDCAPminus| = " << maxHitDiffEndcaps_;
0155
0156 if (!trkQualityStrings.empty()) {
0157 edm::LogInfo("AlignmentTrackSelector") << "Select tracks with these qualities: " << qualities;
0158 }
0159 }
0160
0161 if (applyNHighestPt_)
0162 edm::LogInfo("AlignmentTrackSelector") << "filter N tracks with highest Pt N=" << nHighestPt_;
0163
0164 if (applyMultiplicityFilter_)
0165 edm::LogInfo("AlignmentTrackSelector")
0166 << "apply multiplicity filter N>= " << minMultiplicity_ << "and N<= " << maxMultiplicity_ << " on "
0167 << (multiplicityOnInput_ ? "input" : "output");
0168
0169 if (applyPrescaledHitsFilter_) {
0170 edm::LogInfo("AlignmentTrackSelector") << "apply cut on number of prescaled hits N>= " << minPrescaledHits_
0171 << " (prescale info from " << clusterValueMapTag_ << ")";
0172 }
0173
0174
0175 if (RorZofFirstHitMin_.size() != 2) {
0176 throw cms::Exception("BadConfig") << "@SUB=AlignmentTrackSelector::AlignmentTrackSelector"
0177 << "Wrong configuration of 'RorZofFirstHitMin'."
0178 << " Must have exactly 2 values instead of configured "
0179 << RorZofFirstHitMin_.size() << ")";
0180 } else {
0181 RorZofFirstHitMin_.at(0) = std::fabs(RorZofFirstHitMin_.at(0));
0182 RorZofFirstHitMin_.at(1) = std::fabs(RorZofFirstHitMin_.at(1));
0183 }
0184 if (RorZofFirstHitMax_.size() != 2) {
0185 throw cms::Exception("BadConfig") << "@SUB=AlignmentTrackSelector::AlignmentTrackSelector"
0186 << "Wrong configuration of 'RorZofFirstHitMax'."
0187 << " Must have exactly 2 values instead of configured "
0188 << RorZofFirstHitMax_.size() << ")";
0189 } else {
0190 RorZofFirstHitMax_.at(0) = std::fabs(RorZofFirstHitMax_.at(0));
0191 RorZofFirstHitMax_.at(1) = std::fabs(RorZofFirstHitMax_.at(1));
0192 }
0193 if (RorZofLastHitMin_.size() != 2) {
0194 throw cms::Exception("BadConfig") << "@SUB=AlignmentTrackSelector::AlignmentTrackSelector"
0195 << "Wrong configuration of 'RorZofLastHitMin'."
0196 << " Must have exactly 2 values instead of configured "
0197 << RorZofLastHitMin_.size() << ")";
0198 } else {
0199 RorZofLastHitMin_.at(0) = std::fabs(RorZofLastHitMin_.at(0));
0200 RorZofLastHitMin_.at(1) = std::fabs(RorZofLastHitMin_.at(1));
0201 }
0202 if (RorZofLastHitMax_.size() != 2) {
0203 throw cms::Exception("BadConfig") << "@SUB=AlignmentTrackSelector::AlignmentTrackSelector"
0204 << "Wrong configuration of 'RorZofLastHitMax'."
0205 << " Must have exactly 2 values instead of configured "
0206 << RorZofLastHitMax_.size() << ")";
0207 } else {
0208 RorZofLastHitMax_.at(0) = std::fabs(RorZofLastHitMax_.at(0));
0209 RorZofLastHitMax_.at(1) = std::fabs(RorZofLastHitMax_.at(1));
0210 }
0211
0212 if (RorZofFirstHitMin_.at(0) > RorZofLastHitMax_.at(0) && RorZofFirstHitMin_.at(1) > RorZofLastHitMax_.at(1)) {
0213 throw cms::Exception("BadConfig") << "@SUB=AlignmentTrackSelector::AlignmentTrackSelector"
0214 << "Position of the first hit is set to larger distance than the last hit:."
0215 << " First hit(min): [" << RorZofFirstHitMin_.at(0) << ", "
0216 << RorZofFirstHitMin_.at(1) << "]; Last hit(max): [" << RorZofLastHitMax_.at(0)
0217 << ", " << RorZofLastHitMax_.at(1) << "];";
0218 }
0219 }
0220
0221 void AlignmentTrackSelector::fillPSetDescription(edm::ParameterSetDescription& desc) {
0222
0223 desc.add<bool>("applyBasicCuts", true);
0224 desc.add<double>("ptMin", 0.0);
0225 desc.add<double>("ptMax", 999.0);
0226 desc.add<double>("pMin", 0.0);
0227 desc.add<double>("pMax", 9999.0);
0228 desc.add<double>("etaMin", -2.6);
0229 desc.add<double>("etaMax", 2.6);
0230 desc.add<double>("phiMax", 3.1416);
0231 desc.add<double>("phiMin", -3.1416);
0232 desc.add<double>("chi2nMax", 999999.0);
0233 desc.add<int>("theCharge", 0);
0234 desc.add<double>("d0Min", -999999.0);
0235 desc.add<double>("d0Max", 999999.0);
0236 desc.add<double>("dzMin", -999999.0);
0237 desc.add<double>("dzMax", 999999.0);
0238 desc.add<double>("nHitMin", 0.0);
0239 desc.add<double>("nHitMax", 999.0);
0240 desc.add<double>("nLostHitMax", 999.0);
0241 desc.add<unsigned int>("nHitMin2D", 0);
0242 desc.add<std::vector<double>>("RorZofFirstHitMin", {0.0, 0.0});
0243 desc.add<std::vector<double>>("RorZofFirstHitMax", {999.0, 999.0});
0244 desc.add<std::vector<double>>("RorZofLastHitMin", {0.0, 0.0});
0245 desc.add<std::vector<double>>("RorZofLastHitMax", {999.0, 999.0});
0246 desc.add<bool>("countStereoHitAs2D", true);
0247
0248
0249 edm::ParameterSetDescription minHitsPerSubDetDesc;
0250 minHitsPerSubDetDesc.add<int>("inTEC", 0);
0251 minHitsPerSubDetDesc.add<int>("inTOB", 0);
0252 minHitsPerSubDetDesc.add<int>("inFPIX", 0);
0253 minHitsPerSubDetDesc.add<int>("inTID", 0);
0254 minHitsPerSubDetDesc.add<int>("inBPIX", 0);
0255 minHitsPerSubDetDesc.add<int>("inTIB", 0);
0256 minHitsPerSubDetDesc.add<int>("inPIXEL", 0);
0257 minHitsPerSubDetDesc.add<int>("inTIDplus", 0);
0258 minHitsPerSubDetDesc.add<int>("inTIDminus", 0);
0259 minHitsPerSubDetDesc.add<int>("inTECplus", 0);
0260 minHitsPerSubDetDesc.add<int>("inTECminus", 0);
0261 minHitsPerSubDetDesc.add<int>("inFPIXplus", 0);
0262 minHitsPerSubDetDesc.add<int>("inFPIXminus", 0);
0263 minHitsPerSubDetDesc.add<int>("inENDCAP", 0);
0264 minHitsPerSubDetDesc.add<int>("inENDCAPplus", 0);
0265 minHitsPerSubDetDesc.add<int>("inENDCAPminus", 0);
0266 desc.add<edm::ParameterSetDescription>("minHitsPerSubDet", minHitsPerSubDetDesc);
0267
0268 desc.add<double>("maxHitDiffEndcaps", 999.0);
0269 desc.add<int>("seedOnlyFrom", 0);
0270
0271
0272 desc.add<bool>("applyMultiplicityFilter", false);
0273 desc.add<int>("minMultiplicity", 1);
0274 desc.add<int>("maxMultiplicity", 999999);
0275 desc.add<bool>("multiplicityOnInput", false);
0276
0277
0278 desc.add<bool>("applyNHighestPt", false);
0279 desc.add<int>("nHighestPt", 2);
0280
0281
0282 desc.add<edm::InputTag>("rphirecHits", edm::InputTag("siStripMatchedRecHits", "rphiRecHit"));
0283 desc.add<edm::InputTag>("matchedrecHits", edm::InputTag("siStripMatchedRecHits", "matchedRecHit"));
0284 desc.add<bool>("applyIsolationCut", false);
0285 desc.add<double>("minHitIsolation", 0.01);
0286 desc.add<bool>("applyChargeCheck", false);
0287 desc.add<double>("minHitChargeStrip", 20.0);
0288
0289
0290 desc.add<std::vector<std::string>>("trackQualities", {});
0291 desc.add<std::vector<std::string>>("iterativeTrackingSteps", {});
0292
0293
0294 desc.add<edm::InputTag>("hitPrescaleMapTag", edm::InputTag(""));
0295 desc.add<int>("minPrescaledHits", -1);
0296 }
0297
0298
0299
0300 AlignmentTrackSelector::~AlignmentTrackSelector() = default;
0301
0302
0303
0304 AlignmentTrackSelector::Tracks AlignmentTrackSelector::select(const Tracks& tracks,
0305 const edm::Event& evt,
0306 const edm::EventSetup& eSetup) const {
0307 if (applyMultiplicityFilter_ && multiplicityOnInput_ &&
0308 (tracks.size() < static_cast<unsigned int>(minMultiplicity_) ||
0309 tracks.size() > static_cast<unsigned int>(maxMultiplicity_))) {
0310 return Tracks();
0311 }
0312
0313 Tracks result = tracks;
0314
0315 if (applyBasicCuts_)
0316 result = this->basicCuts(result, evt, eSetup);
0317
0318
0319 if (applyNHighestPt_)
0320 result = this->theNHighestPtTracks(result);
0321
0322
0323 if (applyMultiplicityFilter_ && !multiplicityOnInput_) {
0324 if (result.size() < static_cast<unsigned int>(minMultiplicity_) ||
0325 result.size() > static_cast<unsigned int>(maxMultiplicity_)) {
0326 result.clear();
0327 }
0328 }
0329
0330 if (applyPrescaledHitsFilter_) {
0331 result = this->checkPrescaledHits(result, evt);
0332 }
0333
0334 return result;
0335 }
0336
0337
0338 bool AlignmentTrackSelector::useThisFilter() {
0339 return applyMultiplicityFilter_ || applyBasicCuts_ || applyNHighestPt_ || applyPrescaledHitsFilter_;
0340 }
0341
0342
0343
0344 AlignmentTrackSelector::Tracks AlignmentTrackSelector::basicCuts(const Tracks& tracks,
0345 const edm::Event& evt,
0346 const edm::EventSetup& eSetup) const {
0347 Tracks result;
0348
0349 for (Tracks::const_iterator it = tracks.begin(); it != tracks.end(); ++it) {
0350 const reco::Track* trackp = *it;
0351 float pt = trackp->pt();
0352 float p = trackp->p();
0353 float eta = trackp->eta();
0354 float phi = trackp->phi();
0355 int nhit = trackp->numberOfValidHits();
0356 int nlosthit = trackp->numberOfLostHits();
0357 float chi2n = trackp->normalizedChi2();
0358
0359 int q = trackp->charge();
0360 bool isChargeOk = false;
0361 if (theCharge_ == -1 && q < 0)
0362 isChargeOk = true;
0363 else if (theCharge_ == 1 && q > 0)
0364 isChargeOk = true;
0365 else if (theCharge_ == 0)
0366 isChargeOk = true;
0367
0368 float d0 = trackp->d0();
0369 float dz = trackp->dz();
0370
0371
0372
0373
0374 if (pt > ptMin_ && pt < ptMax_ && p > pMin_ && p < pMax_ && eta > etaMin_ && eta < etaMax_ && phi > phiMin_ &&
0375 phi < phiMax_ && nhit >= nHitMin_ && nhit <= nHitMax_ && nlosthit <= nLostHitMax_ && chi2n < chi2nMax_ &&
0376 isChargeOk && d0 >= d0Min_ && d0 <= d0Max_ && dz >= dzMin_ && dz <= dzMax_) {
0377 bool trkQualityOk = false;
0378 if (!applyTrkQualityCheck_ && !applyIterStepCheck_)
0379 trkQualityOk = true;
0380 else
0381 trkQualityOk = this->isOkTrkQuality(trackp);
0382
0383 bool hitsCheckOk = this->detailedHitsCheck(trackp, evt, eSetup);
0384
0385 if (trkQualityOk && hitsCheckOk)
0386 result.push_back(trackp);
0387 }
0388 }
0389
0390 return result;
0391 }
0392
0393
0394
0395 bool AlignmentTrackSelector::detailedHitsCheck(const reco::Track* trackp,
0396 const edm::Event& evt,
0397 const edm::EventSetup& eSetup) const {
0398
0399 const TrackerTopology* const tTopo = &eSetup.getData(tTopoToken_);
0400
0401
0402
0403 if (minHitsinTIB_ || minHitsinTOB_ || minHitsinTID_ || minHitsinTEC_ || minHitsinENDCAP_ || minHitsinENDCAPplus_ ||
0404 minHitsinENDCAPminus_ || (maxHitDiffEndcaps_ < 999) || minHitsinTIDplus_ || minHitsinTIDminus_ ||
0405 minHitsinFPIXplus_ || minHitsinFPIXminus_ || minHitsinTECplus_ || minHitsinTECminus_ || minHitsinFPIX_ ||
0406 minHitsinBPIX_ || minHitsinPIX_ || nHitMin2D_ || chargeCheck_ || applyIsolation_ ||
0407 (seedOnlyFromAbove_ == 1 || seedOnlyFromAbove_ == 2) || !RorZofFirstHitMin_.empty() ||
0408 !RorZofFirstHitMax_.empty() || !RorZofLastHitMin_.empty() || !RorZofLastHitMax_.empty()) {
0409
0410
0411 int nhitinTIB = 0, nhitinTOB = 0, nhitinTID = 0;
0412 int nhitinTEC = 0, nhitinBPIX = 0, nhitinFPIX = 0, nhitinPIXEL = 0;
0413 int nhitinENDCAP = 0, nhitinENDCAPplus = 0, nhitinENDCAPminus = 0;
0414 int nhitinTIDplus = 0, nhitinTIDminus = 0;
0415 int nhitinFPIXplus = 0, nhitinFPIXminus = 0;
0416 int nhitinTECplus = 0, nhitinTECminus = 0;
0417 unsigned int nHit2D = 0;
0418 unsigned int thishit = 0;
0419
0420 for (auto const& hit : trackp->recHits()) {
0421 thishit++;
0422 const DetId detId(hit->geographicalId());
0423 const int subdetId = detId.subdetId();
0424
0425
0426
0427
0428
0429
0430 if (seedOnlyFromAbove_ == 1 && thishit == 1 &&
0431 (subdetId == int(SiStripDetId::TOB) || subdetId == int(SiStripDetId::TEC))) {
0432 return false;
0433 }
0434 if (seedOnlyFromAbove_ == 2 && thishit == 1 && subdetId == int(SiStripDetId::TIB)) {
0435 return false;
0436 }
0437
0438 if (!hit->isValid())
0439 continue;
0440 if (detId.det() != DetId::Tracker) {
0441 edm::LogError("DetectorMismatch")
0442 << "@SUB=AlignmentTrackSelector::detailedHitsCheck"
0443 << "DetId.det() != DetId::Tracker (=" << DetId::Tracker << "), but " << detId.det() << ".";
0444 }
0445 if (chargeCheck_ && !(this->isOkCharge(hit)))
0446 return false;
0447 if (applyIsolation_ && (!this->isIsolated(hit, evt)))
0448 return false;
0449 if (SiStripDetId::TIB == subdetId)
0450 ++nhitinTIB;
0451 else if (SiStripDetId::TOB == subdetId)
0452 ++nhitinTOB;
0453 else if (SiStripDetId::TID == subdetId) {
0454 ++nhitinTID;
0455 ++nhitinENDCAP;
0456
0457 if (tTopo->tidIsZMinusSide(detId)) {
0458 ++nhitinTIDminus;
0459 ++nhitinENDCAPminus;
0460 } else if (tTopo->tidIsZPlusSide(detId)) {
0461 ++nhitinTIDplus;
0462 ++nhitinENDCAPplus;
0463 }
0464 } else if (SiStripDetId::TEC == subdetId) {
0465 ++nhitinTEC;
0466 ++nhitinENDCAP;
0467
0468 if (tTopo->tecIsZMinusSide(detId)) {
0469 ++nhitinTECminus;
0470 ++nhitinENDCAPminus;
0471 } else if (tTopo->tecIsZPlusSide(detId)) {
0472 ++nhitinTECplus;
0473 ++nhitinENDCAPplus;
0474 }
0475 } else if (kBPIX == subdetId) {
0476 ++nhitinBPIX;
0477 ++nhitinPIXEL;
0478 } else if (kFPIX == subdetId) {
0479 ++nhitinFPIX;
0480 ++nhitinPIXEL;
0481
0482 if (tTopo->pxfSide(detId) == 1)
0483 ++nhitinFPIXminus;
0484 else if (tTopo->pxfSide(detId) == 2)
0485 ++nhitinFPIXplus;
0486 }
0487
0488 if (nHit2D < nHitMin2D_ && this->isHit2D(*hit))
0489 ++nHit2D;
0490 }
0491
0492
0493 bool passedLastHitPositionR = true;
0494 bool passedLastHitPositionZ = true;
0495 bool passedFirstHitPositionR = true;
0496 bool passedFirstHitPositionZ = true;
0497
0498 if (RorZofFirstHitMin_.at(0) != 0.0 || RorZofFirstHitMin_.at(1) != 0.0 || RorZofFirstHitMax_.at(0) != 999.0 ||
0499 RorZofFirstHitMax_.at(1) != 999.0) {
0500 const reco::TrackBase::Point& firstPoint(trackp->innerPosition());
0501
0502 if ((std::fabs(firstPoint.R()) < RorZofFirstHitMin_.at(0)))
0503 passedFirstHitPositionR = false;
0504 if ((std::fabs(firstPoint.R()) > RorZofFirstHitMax_.at(0)))
0505 passedFirstHitPositionR = false;
0506 if ((std::fabs(firstPoint.Z()) < RorZofFirstHitMin_.at(1)))
0507 passedFirstHitPositionZ = false;
0508 if ((std::fabs(firstPoint.Z()) > RorZofFirstHitMax_.at(1)))
0509 passedFirstHitPositionZ = false;
0510 }
0511
0512 if (RorZofLastHitMin_.at(0) != 0.0 || RorZofLastHitMin_.at(1) != 0.0 || RorZofLastHitMax_.at(0) != 999.0 ||
0513 RorZofLastHitMax_.at(1) != 999.0) {
0514 const reco::TrackBase::Point& lastPoint(trackp->outerPosition());
0515
0516 if ((std::fabs(lastPoint.R()) < RorZofLastHitMin_.at(0)))
0517 passedLastHitPositionR = false;
0518 if ((std::fabs(lastPoint.R()) > RorZofLastHitMax_.at(0)))
0519 passedLastHitPositionR = false;
0520 if ((std::fabs(lastPoint.Z()) < RorZofLastHitMin_.at(1)))
0521 passedLastHitPositionZ = false;
0522 if ((std::fabs(lastPoint.Z()) > RorZofLastHitMax_.at(1)))
0523 passedLastHitPositionZ = false;
0524 }
0525
0526 bool passedFirstHitPosition = passedFirstHitPositionR || passedFirstHitPositionZ;
0527 bool passedLastHitPosition = passedLastHitPositionR || passedLastHitPositionZ;
0528
0529 return (nhitinTIB >= minHitsinTIB_ && nhitinTOB >= minHitsinTOB_ && nhitinTID >= minHitsinTID_ &&
0530 nhitinTEC >= minHitsinTEC_ && nhitinENDCAP >= minHitsinENDCAP_ &&
0531 nhitinENDCAPplus >= minHitsinENDCAPplus_ && nhitinENDCAPminus >= minHitsinENDCAPminus_ &&
0532 std::abs(nhitinENDCAPplus - nhitinENDCAPminus) <= maxHitDiffEndcaps_ &&
0533 nhitinTIDplus >= minHitsinTIDplus_ && nhitinTIDminus >= minHitsinTIDminus_ &&
0534 nhitinFPIXplus >= minHitsinFPIXplus_ && nhitinFPIXminus >= minHitsinFPIXminus_ &&
0535 nhitinTECplus >= minHitsinTECplus_ && nhitinTECminus >= minHitsinTECminus_ &&
0536 nhitinBPIX >= minHitsinBPIX_ && nhitinFPIX >= minHitsinFPIX_ && nhitinPIXEL >= minHitsinPIX_ &&
0537 nHit2D >= nHitMin2D_ && passedFirstHitPosition && passedLastHitPosition);
0538 } else {
0539 return true;
0540 }
0541 }
0542
0543
0544
0545 bool AlignmentTrackSelector::isHit2D(const TrackingRecHit& hit) const {
0546
0547
0548 if (!hit.isValid() || (hit.dimension() < 2 && !countStereoHitAs2D_ && !dynamic_cast<const SiStripRecHit1D*>(&hit))) {
0549 return false;
0550 } else {
0551 const DetId detId(hit.geographicalId());
0552 if (detId.det() == DetId::Tracker) {
0553 if (detId.subdetId() == kBPIX || detId.subdetId() == kFPIX) {
0554 return true;
0555 } else {
0556 const SiStripDetId stripId(detId);
0557 if (stripId.stereo())
0558 return countStereoHitAs2D_;
0559 else if (dynamic_cast<const SiStripRecHit1D*>(&hit) || dynamic_cast<const SiStripRecHit2D*>(&hit))
0560 return false;
0561
0562 else if (dynamic_cast<const SiStripMatchedRecHit2D*>(&hit))
0563 return true;
0564 else if (dynamic_cast<const ProjectedSiStripRecHit2D*>(&hit)) {
0565 const ProjectedSiStripRecHit2D* pH = static_cast<const ProjectedSiStripRecHit2D*>(&hit);
0566 return (countStereoHitAs2D_ && this->isHit2D(pH->originalHit()));
0567 } else {
0568 edm::LogError("UnknownType") << "@SUB=AlignmentTrackSelector::isHit2D"
0569 << "Tracker hit not in pixel, neither SiStripRecHit[12]D nor "
0570 << "SiStripMatchedRecHit2D nor ProjectedSiStripRecHit2D.";
0571 return false;
0572 }
0573 }
0574 } else {
0575 edm::LogWarning("DetectorMismatch") << "@SUB=AlignmentTrackSelector::isHit2D"
0576 << "Hit not in tracker with 'official' dimension >=2.";
0577 return true;
0578 }
0579 }
0580
0581 }
0582
0583
0584
0585 bool AlignmentTrackSelector::isOkCharge(const TrackingRecHit* hit) const {
0586 if (!hit || !hit->isValid())
0587 return true;
0588
0589
0590 const DetId id(hit->geographicalId());
0591 if (id.det() != DetId::Tracker) {
0592 edm::LogWarning("DetectorMismatch") << "@SUB=isOkCharge"
0593 << "Hit not in tracker!";
0594 return true;
0595 }
0596 if (id.subdetId() == kFPIX || id.subdetId() == kBPIX) {
0597 return true;
0598 }
0599
0600
0601 const std::type_info& type = typeid(*hit);
0602
0603 if (type == typeid(SiStripRecHit2D)) {
0604 const SiStripRecHit2D* stripHit2D = dynamic_cast<const SiStripRecHit2D*>(hit);
0605 if (stripHit2D) {
0606 return this->isOkChargeStripHit(*stripHit2D);
0607 }
0608 } else if (type == typeid(SiStripRecHit1D)) {
0609 const SiStripRecHit1D* stripHit1D = dynamic_cast<const SiStripRecHit1D*>(hit);
0610 if (stripHit1D) {
0611 return this->isOkChargeStripHit(*stripHit1D);
0612 }
0613 }
0614
0615 else if (type ==
0616 typeid(SiStripMatchedRecHit2D)) {
0617 const SiStripMatchedRecHit2D* matchedHit = dynamic_cast<const SiStripMatchedRecHit2D*>(hit);
0618 if (matchedHit) {
0619 return (this->isOkChargeStripHit(matchedHit->monoHit()) && this->isOkChargeStripHit(matchedHit->stereoHit()));
0620 }
0621 } else if (type == typeid(ProjectedSiStripRecHit2D)) {
0622
0623 const ProjectedSiStripRecHit2D* projHit = dynamic_cast<const ProjectedSiStripRecHit2D*>(hit);
0624 if (projHit) {
0625 return this->isOkChargeStripHit(projHit->originalHit());
0626 }
0627 } else {
0628 edm::LogError("AlignmentTrackSelector") << "@SUB=isOkCharge"
0629 << "Unknown type of a valid tracker hit in Strips "
0630 << " SubDet = " << id.subdetId();
0631 return false;
0632 }
0633
0634
0635
0636 edm::LogError("AlignmentTrackSelector") << "@SUB=isOkCharge"
0637 << "Unknown valid tracker hit not in pixel, subdet " << id.subdetId()
0638 << ", SiTrackerMultiRecHit " << dynamic_cast<const SiTrackerMultiRecHit*>(hit)
0639 << ", BaseTrackerRecHit " << dynamic_cast<const BaseTrackerRecHit*>(hit);
0640
0641 return true;
0642 }
0643
0644
0645
0646 bool AlignmentTrackSelector::isOkChargeStripHit(const SiStripRecHit2D& siStripRecHit2D) const {
0647 double charge = 0.;
0648
0649 SiStripRecHit2D::ClusterRef cluster(siStripRecHit2D.cluster());
0650 const auto& amplitudes = cluster->amplitudes();
0651
0652 for (size_t ia = 0; ia < amplitudes.size(); ++ia) {
0653 charge += amplitudes[ia];
0654 }
0655
0656 return (charge >= minHitChargeStrip_);
0657 }
0658
0659
0660
0661 bool AlignmentTrackSelector::isOkChargeStripHit(const SiStripRecHit1D& siStripRecHit1D) const {
0662 double charge = 0.;
0663
0664 SiStripRecHit1D::ClusterRef cluster(siStripRecHit1D.cluster());
0665 const auto& amplitudes = cluster->amplitudes();
0666
0667 for (size_t ia = 0; ia < amplitudes.size(); ++ia) {
0668 charge += amplitudes[ia];
0669 }
0670
0671 return (charge >= minHitChargeStrip_);
0672 }
0673
0674
0675
0676 bool AlignmentTrackSelector::isIsolated(const TrackingRecHit* hit, const edm::Event& evt) const {
0677
0678
0679
0680
0681 edm::Handle<SiStripRecHit2DCollection> rphirecHits;
0682 edm::Handle<SiStripMatchedRecHit2DCollection> matchedrecHits;
0683
0684 evt.getByToken(rphirecHitsToken_, rphirecHits);
0685 evt.getByToken(matchedrecHitsToken_, matchedrecHits);
0686
0687 SiStripRecHit2DCollection::DataContainer::const_iterator istripSt;
0688 SiStripMatchedRecHit2DCollection::DataContainer::const_iterator istripStm;
0689 const SiStripRecHit2DCollection::DataContainer& stripcollSt = rphirecHits->data();
0690 const SiStripMatchedRecHit2DCollection::DataContainer& stripcollStm = matchedrecHits->data();
0691
0692 DetId idet = hit->geographicalId();
0693
0694
0695
0696
0697 for (istripSt = stripcollSt.begin(); istripSt != stripcollSt.end(); ++istripSt) {
0698 const SiStripRecHit2D* aHit = &*(istripSt);
0699 DetId mydet1 = aHit->geographicalId();
0700 if (idet.rawId() != mydet1.rawId())
0701 continue;
0702 float theDistance = (hit->localPosition() - aHit->localPosition()).mag();
0703
0704 if (theDistance > 0.001 && theDistance < minHitIsolation_)
0705 return false;
0706 }
0707
0708
0709 for (istripStm = stripcollStm.begin(); istripStm != stripcollStm.end(); ++istripStm) {
0710 const SiStripMatchedRecHit2D* aHit = &*(istripStm);
0711 DetId mydet2 = aHit->geographicalId();
0712 if (idet.rawId() != mydet2.rawId())
0713 continue;
0714 float theDistance = (hit->localPosition() - aHit->localPosition()).mag();
0715
0716 if (theDistance > 0.001 && theDistance < minHitIsolation_)
0717 return false;
0718 }
0719 return true;
0720 }
0721
0722
0723
0724 AlignmentTrackSelector::Tracks AlignmentTrackSelector::theNHighestPtTracks(const Tracks& tracks) const {
0725 Tracks sortedTracks = tracks;
0726 Tracks result;
0727
0728
0729 std::sort(sortedTracks.begin(), sortedTracks.end(), ptComparator);
0730
0731
0732 int n = 0;
0733 for (Tracks::const_iterator it = sortedTracks.begin(); it != sortedTracks.end(); ++it) {
0734 if (n < nHighestPt_) {
0735 result.push_back(*it);
0736 n++;
0737 }
0738 }
0739
0740 return result;
0741 }
0742
0743
0744
0745 AlignmentTrackSelector::Tracks AlignmentTrackSelector::checkPrescaledHits(const Tracks& tracks,
0746 const edm::Event& evt) const {
0747 Tracks result;
0748
0749
0750 edm::Handle<AliClusterValueMap> fMap;
0751 evt.getByToken(clusterValueMapToken_, fMap);
0752 const AliClusterValueMap& flagMap = *fMap;
0753
0754
0755 for (auto const& trackp : tracks) {
0756 int ntakenhits = 0;
0757
0758
0759 for (auto const& hit : trackp->recHits()) {
0760 if (!hit->isValid())
0761 continue;
0762 DetId detid = hit->geographicalId();
0763 int subDet = detid.subdetId();
0764 AlignmentClusterFlag flag;
0765
0766 bool isPixelHit = (subDet == kFPIX || subDet == kBPIX);
0767
0768 if (!isPixelHit) {
0769 const std::type_info& type = typeid(*hit);
0770
0771 if (type == typeid(SiStripRecHit2D)) {
0772 const SiStripRecHit2D* striphit = dynamic_cast<const SiStripRecHit2D*>(hit);
0773 if (striphit != nullptr) {
0774 SiStripRecHit2D::ClusterRef stripclust(striphit->cluster());
0775 flag = flagMap[stripclust];
0776 }
0777 } else if (type == typeid(SiStripRecHit1D)) {
0778 const SiStripRecHit1D* striphit = dynamic_cast<const SiStripRecHit1D*>(hit);
0779 if (striphit != nullptr) {
0780 SiStripRecHit1D::ClusterRef stripclust(striphit->cluster());
0781 flag = flagMap[stripclust];
0782 }
0783 } else {
0784 edm::LogError("AlignmentTrackSelector")
0785 << "ERROR in <AlignmentTrackSelector::checkPrescaledHits>: Dynamic cast of Strip RecHit failed!"
0786 << " Skipping this hit.";
0787 continue;
0788 }
0789
0790 }
0791 else {
0792 const SiPixelRecHit* pixelhit = dynamic_cast<const SiPixelRecHit*>(hit);
0793 if (pixelhit != nullptr) {
0794 SiPixelRecHit::ClusterRef pixclust(pixelhit->cluster());
0795 flag = flagMap[pixclust];
0796 } else {
0797 edm::LogError("AlignmentTrackSelector")
0798 << "ERROR in <AlignmentTrackSelector::checkPrescaledHits>: Dynamic cast of Pixel RecHit failed! ";
0799 }
0800 }
0801
0802 if (flag.isTaken())
0803 ntakenhits++;
0804
0805 }
0806 if (ntakenhits >= minPrescaledHits_)
0807 result.push_back(trackp);
0808 }
0809
0810 return result;
0811 }
0812
0813
0814
0815 bool AlignmentTrackSelector::isOkTrkQuality(const reco::Track* track) const {
0816 bool qualityOk = false;
0817 bool iterStepOk = false;
0818
0819
0820 if (applyIterStepCheck_) {
0821 for (unsigned int i = 0; i < trkSteps_.size(); ++i) {
0822 if (track->algo() == (trkSteps_[i])) {
0823 iterStepOk = true;
0824 }
0825 }
0826 } else
0827 iterStepOk = true;
0828
0829
0830 if (applyTrkQualityCheck_) {
0831 for (unsigned int i = 0; i < trkQualities_.size(); ++i) {
0832 if (track->quality(trkQualities_[i])) {
0833 qualityOk = true;
0834 }
0835 }
0836 } else
0837 qualityOk = true;
0838
0839 return qualityOk && iterStepOk;
0840 }