File indexing completed on 2025-01-09 23:33:54
0001 #include <sstream>
0002
0003 #include "FWCore/Framework/interface/ESWatcher.h"
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0007 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0008 #include "RecoTracker/CkfPattern/interface/BaseCkfTrajectoryBuilderFactory.h"
0009 #include "RecoTracker/CkfPattern/interface/CkfTrajectoryBuilder.h"
0010 #include "RecoTracker/MeasurementDet/interface/MeasurementTracker.h"
0011 #include "RecoTracker/MeasurementDet/interface/MeasurementTrackerEvent.h"
0012 #include "RecoTracker/TkDetLayers/interface/GeometricSearchTracker.h"
0013 #include "TrackingTools/DetLayers/interface/NavigationSchool.h"
0014 #include "TrackingTools/KalmanUpdators/interface/Chi2MeasurementEstimatorBase.h"
0015 #include "TrackingTools/MeasurementDet/interface/LayerMeasurements.h"
0016 #include "TrackingTools/PatternTools/interface/TrajMeasLessEstim.h"
0017 #include "TrackingTools/PatternTools/interface/TransverseImpactPointExtrapolator.h"
0018 #include "TrackingTools/TrajectoryFiltering/interface/TrajectoryFilter.h"
0019 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0020
0021 #include "EtaPhiEstimator.h"
0022
0023 class TrackingComponentsRecord;
0024
0025 class MuonCkfTrajectoryBuilder : public CkfTrajectoryBuilder {
0026 public:
0027 MuonCkfTrajectoryBuilder(const edm::ParameterSet& conf, edm::ConsumesCollector& iC);
0028 ~MuonCkfTrajectoryBuilder() override;
0029
0030 static void fillPSetDescription(edm::ParameterSetDescription& iDesc);
0031
0032 protected:
0033 void setEvent_(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0034
0035 void collectMeasurement(const DetLayer* layer,
0036 const std::vector<const DetLayer*>& nl,
0037 const TrajectoryStateOnSurface& currentState,
0038 std::vector<TM>& result,
0039 int& invalidHits,
0040 const Propagator*) const;
0041
0042 void findCompatibleMeasurements(const TrajectorySeed& seed,
0043 const TempTrajectory& traj,
0044 std::vector<TrajectoryMeasurement>& result) const override;
0045
0046
0047 bool theUseSeedLayer;
0048 double theRescaleErrorIfFail;
0049 const double theDeltaEta;
0050 const double theDeltaPhi;
0051 const std::string theProximityPropagatorName;
0052 const Propagator* theProximityPropagator;
0053 const edm::ESGetToken<Propagator, TrackingComponentsRecord> thePropagatorToken;
0054 edm::ESWatcher<BaseCkfTrajectoryBuilder::Chi2MeasurementEstimatorRecord> theEstimatorWatcher;
0055 std::unique_ptr<Chi2MeasurementEstimatorBase> theEtaPhiEstimator;
0056 };
0057
0058 MuonCkfTrajectoryBuilder::MuonCkfTrajectoryBuilder(const edm::ParameterSet& conf, edm::ConsumesCollector& iC)
0059 : CkfTrajectoryBuilder(conf, iC),
0060 theDeltaEta(conf.getParameter<double>("deltaEta")),
0061 theDeltaPhi(conf.getParameter<double>("deltaPhi")),
0062 theProximityPropagatorName(conf.getParameter<std::string>("propagatorProximity")),
0063 theProximityPropagator(nullptr),
0064 thePropagatorToken(iC.esConsumes(edm::ESInputTag("", theProximityPropagatorName))),
0065 theEtaPhiEstimator(nullptr) {
0066
0067 theUseSeedLayer = conf.getParameter<bool>("useSeedLayer");
0068 theRescaleErrorIfFail = conf.getParameter<double>("rescaleErrorIfFail");
0069 }
0070
0071 MuonCkfTrajectoryBuilder::~MuonCkfTrajectoryBuilder() {}
0072
0073 void MuonCkfTrajectoryBuilder::fillPSetDescription(edm::ParameterSetDescription& iDesc) {
0074 CkfTrajectoryBuilder::fillPSetDescription(iDesc);
0075 iDesc.add<double>("deltaEta", .1);
0076 iDesc.add<double>("deltaPhi", .1);
0077 iDesc.add<std::string>("propagatorProximity", "SteppingHelixPropagatorAny");
0078 iDesc.add<bool>("useSeedLayer", false);
0079 iDesc.add<double>("rescaleErrorIfFail", 1.);
0080 }
0081
0082 void MuonCkfTrajectoryBuilder::setEvent_(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0083 CkfTrajectoryBuilder::setEvent_(iEvent, iSetup);
0084
0085 theProximityPropagator = &iSetup.getData(thePropagatorToken);
0086
0087
0088 if (theEstimatorWatcher.check(iSetup) && theDeltaEta > 0 && theDeltaPhi > 0)
0089 theEtaPhiEstimator = std::make_unique<EtaPhiEstimator>(theDeltaEta, theDeltaPhi, theEstimator);
0090 }
0091
0092
0093
0094
0095
0096
0097
0098
0099
0100
0101
0102
0103
0104
0105
0106
0107
0108
0109
0110
0111
0112
0113
0114
0115
0116
0117
0118
0119
0120
0121
0122 void MuonCkfTrajectoryBuilder::collectMeasurement(const DetLayer* layer,
0123 const std::vector<const DetLayer*>& nl,
0124 const TrajectoryStateOnSurface& currentState,
0125 std::vector<TM>& result,
0126 int& invalidHits,
0127 const Propagator* prop) const {
0128 for (std::vector<const DetLayer*>::const_iterator il = nl.begin(); il != nl.end(); il++) {
0129 TSOS stateToUse = currentState;
0130
0131 if (layer == (*il)) {
0132 LogDebug("CkfPattern") << " self propagating in findCompatibleMeasurements.\n from: \n" << stateToUse;
0133
0134
0135 TransverseImpactPointExtrapolator middle;
0136 GlobalPoint center(0, 0, 0);
0137 stateToUse = middle.extrapolate(stateToUse, center, *prop);
0138
0139 if (!stateToUse.isValid())
0140 continue;
0141 LogDebug("CkfPattern") << "to: " << stateToUse;
0142 }
0143
0144 LayerMeasurements layerMeasurements(theMeasurementTracker->measurementTracker(), *theMeasurementTracker);
0145 std::vector<TM> tmp = layerMeasurements.measurements((**il), stateToUse, *prop, *theEstimator);
0146
0147 if (tmp.size() == 1 && theEtaPhiEstimator) {
0148 LogDebug("CkfPattern") << "only an invalid hit is found. trying differently";
0149 tmp = layerMeasurements.measurements((**il), stateToUse, *prop, *theEtaPhiEstimator);
0150 }
0151 LogDebug("CkfPattern") << tmp.size() << " measurements returned by LayerMeasurements";
0152
0153 if (!tmp.empty()) {
0154
0155
0156
0157
0158 if (result.empty())
0159 result = tmp;
0160 else {
0161
0162 result.insert(result.end() - invalidHits, tmp.begin(), tmp.end());
0163 }
0164 invalidHits++;
0165 }
0166 }
0167
0168 LogDebug("CkfPattern") << "starting from:\n"
0169 << "x: " << currentState.globalPosition() << "\n"
0170 << "p: " << currentState.globalMomentum() << "\n"
0171 << PrintoutHelper::dumpMeasurements(result);
0172 }
0173
0174 void MuonCkfTrajectoryBuilder::findCompatibleMeasurements(const TrajectorySeed& seed,
0175 const TempTrajectory& traj,
0176 std::vector<TrajectoryMeasurement>& result) const {
0177 int invalidHits = 0;
0178
0179 std::vector<const DetLayer*> nl;
0180
0181 if (traj.empty()) {
0182 LogDebug("CkfPattern") << "using JR patch for no measurement case";
0183
0184
0185
0186 PTrajectoryStateOnDet ptod = seed.startingState();
0187 DetId id(ptod.detId());
0188 const GeomDet* g = theMeasurementTracker->geomTracker()->idToDet(id);
0189 const Surface* surface = &g->surface();
0190
0191 TrajectoryStateOnSurface currentState(
0192 trajectoryStateTransform::transientState(ptod, surface, forwardPropagator(seed)->magneticField()));
0193
0194
0195 const DetLayer* l = theMeasurementTracker->geometricSearchTracker()->detLayer(id);
0196
0197 if (theUseSeedLayer) {
0198 {
0199
0200 LogDebug("CkfPattern") << "using the layer of the seed first.";
0201 nl.push_back(l);
0202 collectMeasurement(l, nl, currentState, result, invalidHits, theProximityPropagator);
0203 }
0204
0205
0206 if ((unsigned int)invalidHits == result.size() && theRescaleErrorIfFail != 1.0 && !result.empty()) {
0207 result.clear();
0208 LogDebug("CkfPattern") << "using a rescale by " << theRescaleErrorIfFail << " to find measurements.";
0209 TrajectoryStateOnSurface rescaledCurrentState = currentState;
0210 rescaledCurrentState.rescaleError(theRescaleErrorIfFail);
0211 invalidHits = 0;
0212 collectMeasurement(l, nl, rescaledCurrentState, result, invalidHits, theProximityPropagator);
0213 }
0214 }
0215
0216
0217 if (result.empty() || (unsigned int)invalidHits == result.size()) {
0218 result.clear();
0219 LogDebug("CkfPattern") << "Need to go to next layer to get measurements";
0220
0221 nl = theNavigationSchool->nextLayers(*l, *currentState.freeState(), traj.direction());
0222 if (nl.empty()) {
0223 LogDebug("CkfPattern") << " there was no next layer with wellInside. Use the next with no check.";
0224
0225
0226 nl = theNavigationSchool->nextLayers(*l, ((traj.direction() == alongMomentum) ? insideOut : outsideIn));
0227 }
0228 invalidHits = 0;
0229 collectMeasurement(l, nl, currentState, result, invalidHits, forwardPropagator(seed));
0230 }
0231
0232
0233 if (!result.empty() && (unsigned int)invalidHits == result.size() && theRescaleErrorIfFail != 1.0) {
0234 result.clear();
0235 LogDebug("CkfPattern") << "using a rescale by " << theRescaleErrorIfFail
0236 << " to find measurements on next layers.";
0237 TrajectoryStateOnSurface rescaledCurrentState = currentState;
0238 rescaledCurrentState.rescaleError(theRescaleErrorIfFail);
0239 invalidHits = 0;
0240 collectMeasurement(l, nl, rescaledCurrentState, result, invalidHits, forwardPropagator(seed));
0241 }
0242
0243 } else
0244 {
0245 TSOS currentState(traj.lastMeasurement().updatedState());
0246
0247 nl = theNavigationSchool->nextLayers(*traj.lastLayer(), *currentState.freeState(), traj.direction());
0248 if (nl.empty()) {
0249 LogDebug("CkfPattern") << " no next layers... going " << traj.direction() << "\n from: \n"
0250 << currentState
0251 << "\n from detId: " << traj.lastMeasurement().recHit()->geographicalId().rawId();
0252 return;
0253 }
0254
0255 collectMeasurement(traj.lastLayer(), nl, currentState, result, invalidHits, forwardPropagator(seed));
0256 }
0257
0258
0259 if (result.size() > 1) {
0260 sort(result.begin(), result.end() - invalidHits, TrajMeasLessEstim());
0261 }
0262
0263 #ifdef DEBUG_INVALID
0264 bool afterInvalid = false;
0265 for (std::vector<TM>::const_iterator i = result.begin(); i != result.end(); i++) {
0266 if (!i->recHit().isValid())
0267 afterInvalid = true;
0268 if (afterInvalid && i->recHit().isValid()) {
0269 edm::LogError("CkfPattern") << "CkfTrajectoryBuilder error: valid hit after invalid!";
0270 }
0271 }
0272 #endif
0273
0274
0275 }
0276
0277 #include "FWCore/ParameterSet/interface/ValidatedPluginMacros.h"
0278 DEFINE_EDM_VALIDATED_PLUGIN(BaseCkfTrajectoryBuilderFactory, MuonCkfTrajectoryBuilder, "MuonCkfTrajectoryBuilder");