File indexing completed on 2024-04-06 11:57:21
0001 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0002 #include "DataFormats/GeometrySurface/interface/Surface.h"
0003 #include "DataFormats/CLHEP/interface/AlgebraicObjects.h"
0004 #include "DataFormats/Math/interface/Vector.h"
0005 #include "DataFormats/Math/interface/Error.h"
0006 #include "RecoVertex/VertexTools/interface/PerigeeLinearizedTrackState.h"
0007 #include "Alignment/ReferenceTrajectories/interface/TrajectoryFactoryPlugin.h"
0008
0009 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0010 #include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
0011 #include "Geometry/CommonDetUnit/interface/TrackerGeomDet.h"
0012
0013 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0014
0015 #include "FWCore/Framework/interface/ESHandle.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017
0018 #include "Alignment/ReferenceTrajectories/interface/TrajectoryFactoryBase.h"
0019 #include "Alignment/ReferenceTrajectories/interface/TwoBodyDecayTrajectory.h"
0020
0021 #include "Alignment/TwoBodyDecay/interface/TwoBodyDecayFitter.h"
0022
0023
0024
0025
0026
0027
0028 class TwoBodyDecayTrajectoryFactory : public TrajectoryFactoryBase {
0029 public:
0030 typedef TwoBodyDecayVirtualMeasurement VirtualMeasurement;
0031 typedef TwoBodyDecayTrajectoryState::TsosContainer TsosContainer;
0032 typedef TwoBodyDecayTrajectory::ConstRecHitCollection ConstRecHitCollection;
0033
0034 TwoBodyDecayTrajectoryFactory(const edm::ParameterSet &config, edm::ConsumesCollector &iC);
0035 ~TwoBodyDecayTrajectoryFactory() override;
0036 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> m_MagFieldToken;
0037 const edm::ESGetToken<GlobalTrackingGeometry, GlobalTrackingGeometryRecord> m_globTackingToken;
0038
0039
0040 const ReferenceTrajectoryCollection trajectories(const edm::EventSetup &setup,
0041 const ConstTrajTrackPairCollection &tracks,
0042 const reco::BeamSpot &beamSpot) const override;
0043
0044 const ReferenceTrajectoryCollection trajectories(const edm::EventSetup &setup,
0045 const ConstTrajTrackPairCollection &tracks,
0046 const ExternalPredictionCollection &external,
0047 const reco::BeamSpot &beamSpot) const override;
0048
0049 TwoBodyDecayTrajectoryFactory *clone() const override { return new TwoBodyDecayTrajectoryFactory(*this); }
0050
0051 protected:
0052 const ReferenceTrajectoryCollection constructTrajectories(const ConstTrajTrackPairCollection &tracks,
0053 const TwoBodyDecay &tbd,
0054 const MagneticField *magField,
0055 const reco::BeamSpot &beamSpot,
0056 bool setParameterErrors) const;
0057
0058 bool match(const TrajectoryStateOnSurface &state, const TransientTrackingRecHit::ConstRecHitPointer &recHit) const;
0059
0060 TwoBodyDecayFitter theFitter;
0061
0062 double thePrimaryMass;
0063 double thePrimaryWidth;
0064 double theSecondaryMass;
0065
0066 double theNSigmaCutValue;
0067 double theChi2CutValue;
0068 bool theUseRefittedStateFlag;
0069 bool theConstructTsosWithErrorsFlag;
0070 };
0071
0072
0073
0074
0075
0076 TwoBodyDecayTrajectoryFactory::TwoBodyDecayTrajectoryFactory(const edm::ParameterSet &config,
0077 edm::ConsumesCollector &iC)
0078 : TrajectoryFactoryBase(config, 2, iC),
0079 m_MagFieldToken(iC.esConsumes()),
0080 m_globTackingToken(iC.esConsumes()),
0081 theFitter(config) {
0082 const edm::ParameterSet ppc = config.getParameter<edm::ParameterSet>("ParticleProperties");
0083 thePrimaryMass = ppc.getParameter<double>("PrimaryMass");
0084 thePrimaryWidth = ppc.getParameter<double>("PrimaryWidth");
0085 theSecondaryMass = ppc.getParameter<double>("SecondaryMass");
0086
0087 theNSigmaCutValue = config.getParameter<double>("NSigmaCut");
0088 theChi2CutValue = config.getParameter<double>("Chi2Cut");
0089 theUseRefittedStateFlag = config.getParameter<bool>("UseRefittedState");
0090 theConstructTsosWithErrorsFlag = config.getParameter<bool>("ConstructTsosWithErrors");
0091 }
0092
0093 TwoBodyDecayTrajectoryFactory::~TwoBodyDecayTrajectoryFactory() {}
0094
0095 const TrajectoryFactoryBase::ReferenceTrajectoryCollection TwoBodyDecayTrajectoryFactory::trajectories(
0096 const edm::EventSetup &setup, const ConstTrajTrackPairCollection &tracks, const reco::BeamSpot &beamSpot) const {
0097 ReferenceTrajectoryCollection trajectories;
0098
0099 const MagneticField *magneticField = &setup.getData(m_MagFieldToken);
0100 const GlobalTrackingGeometry *trackingGeometry = &setup.getData(m_globTackingToken);
0101 if (tracks.size() == 2) {
0102
0103 std::vector<reco::TransientTrack> transientTracks(2);
0104
0105 transientTracks[0] = reco::TransientTrack(*tracks[0].second, magneticField);
0106 transientTracks[0].setTrackingGeometry(trackingGeometry);
0107
0108 transientTracks[1] = reco::TransientTrack(*tracks[1].second, magneticField);
0109 transientTracks[1].setTrackingGeometry(trackingGeometry);
0110
0111
0112 VirtualMeasurement vm(thePrimaryMass, thePrimaryWidth, theSecondaryMass, beamSpot);
0113 TwoBodyDecay tbd = theFitter.estimate(transientTracks, vm);
0114
0115 if (!tbd.isValid() || (tbd.chi2() > theChi2CutValue)) {
0116 trajectories.push_back(ReferenceTrajectoryPtr(new TwoBodyDecayTrajectory()));
0117 return trajectories;
0118 }
0119
0120 return constructTrajectories(tracks, tbd, magneticField, beamSpot, false);
0121 } else {
0122 edm::LogInfo("ReferenceTrajectories") << "@SUB=TwoBodyDecayTrajectoryFactory::trajectories"
0123 << "Need 2 tracks, got " << tracks.size() << ".\n";
0124 }
0125
0126 return trajectories;
0127 }
0128
0129 const TrajectoryFactoryBase::ReferenceTrajectoryCollection TwoBodyDecayTrajectoryFactory::trajectories(
0130 const edm::EventSetup &setup,
0131 const ConstTrajTrackPairCollection &tracks,
0132 const ExternalPredictionCollection &external,
0133 const reco::BeamSpot &beamSpot) const {
0134 ReferenceTrajectoryCollection trajectories;
0135
0136 const MagneticField *magneticField = &setup.getData(m_MagFieldToken);
0137
0138 const GlobalTrackingGeometry *trackingGeometry = &setup.getData(m_globTackingToken);
0139
0140 if (tracks.size() == 2 && external.size() == 2) {
0141 if (external[0].isValid() && external[1].isValid())
0142 {
0143
0144 std::vector<reco::TransientTrack> transientTracks(2);
0145
0146 transientTracks[0] = reco::TransientTrack(*tracks[0].second, magneticField);
0147 transientTracks[0].setTrackingGeometry(trackingGeometry);
0148
0149 transientTracks[1] = reco::TransientTrack(*tracks[1].second, magneticField);
0150 transientTracks[1].setTrackingGeometry(trackingGeometry);
0151
0152
0153
0154
0155 VirtualMeasurement vm(thePrimaryMass, thePrimaryWidth, theSecondaryMass, beamSpot);
0156 TwoBodyDecay tbd = theFitter.estimate(transientTracks, external, vm);
0157
0158 if (!tbd.isValid() || (tbd.chi2() > theChi2CutValue)) {
0159 trajectories.push_back(ReferenceTrajectoryPtr(new TwoBodyDecayTrajectory()));
0160 return trajectories;
0161 }
0162
0163 return constructTrajectories(tracks, tbd, magneticField, beamSpot, true);
0164 } else {
0165
0166 trajectories = this->trajectories(setup, tracks, beamSpot);
0167 }
0168 } else {
0169 edm::LogInfo("ReferenceTrajectories") << "@SUB=TwoBodyDecayTrajectoryFactory::trajectories"
0170 << "Need 2 tracks, got " << tracks.size() << ".\n";
0171 }
0172
0173 return trajectories;
0174 }
0175
0176 const TwoBodyDecayTrajectoryFactory::ReferenceTrajectoryCollection TwoBodyDecayTrajectoryFactory::constructTrajectories(
0177 const ConstTrajTrackPairCollection &tracks,
0178 const TwoBodyDecay &tbd,
0179 const MagneticField *magField,
0180 const reco::BeamSpot &beamSpot,
0181 bool setParameterErrors) const {
0182 ReferenceTrajectoryCollection trajectories;
0183
0184
0185 TrajectoryInput input1 = this->innermostStateAndRecHits(tracks[0]);
0186 TrajectoryInput input2 = this->innermostStateAndRecHits(tracks[1]);
0187
0188 if (!(input1.first.isValid() && input2.first.isValid()))
0189 return trajectories;
0190
0191
0192 TsosContainer tsos(input1.first, input2.first);
0193 ConstRecHitCollection recHits(input1.second, input2.second);
0194 TwoBodyDecayTrajectoryState trajectoryState(tsos, tbd, theSecondaryMass, magField);
0195
0196 if (!trajectoryState.isValid()) {
0197 trajectories.push_back(ReferenceTrajectoryPtr(new TwoBodyDecayTrajectory()));
0198 return trajectories;
0199 }
0200
0201
0202
0203
0204
0205
0206 TrackingRecHit::ConstRecHitPointer updatedRecHit1(recHits.first.front());
0207 TrackingRecHit::ConstRecHitPointer updatedRecHit2(recHits.second.front());
0208
0209 bool valid1 = match(trajectoryState.trajectoryStates(true).first, updatedRecHit1);
0210
0211 bool valid2 = match(trajectoryState.trajectoryStates(true).second, updatedRecHit2);
0212
0213 if (!valid1 || !valid2) {
0214 trajectories.push_back(ReferenceTrajectoryPtr(new TwoBodyDecayTrajectory()));
0215 return trajectories;
0216 }
0217
0218 ReferenceTrajectoryBase::Config config(materialEffects(), propagationDirection());
0219 config.useBeamSpot = useBeamSpot_;
0220 config.includeAPEs = includeAPEs_;
0221 config.allowZeroMaterial = allowZeroMaterial_;
0222 config.useRefittedState = theUseRefittedStateFlag;
0223 config.constructTsosWithErrors = theConstructTsosWithErrorsFlag;
0224
0225 config.hitsAreReverse = false;
0226 TwoBodyDecayTrajectory *result = new TwoBodyDecayTrajectory(trajectoryState, recHits, magField, beamSpot, config);
0227 if (setParameterErrors && tbd.hasError())
0228 result->setParameterErrors(tbd.covariance());
0229 trajectories.push_back(ReferenceTrajectoryPtr(result));
0230 return trajectories;
0231 }
0232
0233 bool TwoBodyDecayTrajectoryFactory::match(const TrajectoryStateOnSurface &state,
0234 const TransientTrackingRecHit::ConstRecHitPointer &recHit) const {
0235 LocalPoint lp1 = state.localPosition();
0236 LocalPoint lp2 = recHit->localPosition();
0237
0238 double deltaX = lp1.x() - lp2.x();
0239 double deltaY = lp1.y() - lp2.y();
0240
0241 LocalError le = recHit->localPositionError();
0242
0243 double varX = le.xx();
0244 double varY = le.yy();
0245
0246 return ((fabs(deltaX) / sqrt(varX) < theNSigmaCutValue) && (fabs(deltaY) / sqrt(varY) < theNSigmaCutValue));
0247 }
0248
0249 DEFINE_EDM_PLUGIN(TrajectoryFactoryPlugin, TwoBodyDecayTrajectoryFactory, "TwoBodyDecayTrajectoryFactory");