File indexing completed on 2025-02-05 23:51:51
0001 #ifndef SimTracker_TrackHistory_TrackClassifier_h
0002 #define SimTracker_TrackHistory_TrackClassifier_h
0003
0004 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0005 #include "FWCore/Framework/interface/ConsumesCollector.h"
0006 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/Utilities/interface/InputTag.h"
0009 #include "MagneticField/Engine/interface/MagneticField.h"
0010 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0011 #include "SimGeneral/HepPDTRecord/interface/ParticleDataTable.h"
0012 #include "SimTracker/TrackHistory/interface/CMSProcessTypes.h"
0013 #include "SimTracker/TrackHistory/interface/TrackCategories.h"
0014 #include "SimTracker/TrackHistory/interface/TrackHistory.h"
0015 #include "SimTracker/TrackHistory/interface/TrackQuality.h"
0016 #include "TrackingTools/PatternTools/interface/TSCPBuilderNoMaterial.h"
0017 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0018 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0019 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0020
0021 class TrackTopology;
0022
0023
0024 class TrackClassifier : public TrackCategories {
0025 public:
0026
0027 typedef TrackCategories Categories;
0028
0029
0030 TrackClassifier(edm::ParameterSet const &, edm::ConsumesCollector &&);
0031
0032
0033 void newEvent(edm::Event const &, edm::EventSetup const &);
0034
0035
0036 TrackClassifier const &evaluate(reco::TrackBaseRef const &);
0037
0038
0039 TrackClassifier const &evaluate(TrackingParticleRef const &);
0040
0041
0042 TrackClassifier const &evaluate(reco::TrackRef const &track) { return evaluate(reco::TrackBaseRef(track)); }
0043
0044
0045 TrackHistory const &history() const { return tracer_; }
0046
0047
0048 TrackQuality const &quality() const { return quality_; }
0049
0050 static void fillPSetDescription(edm::ParameterSetDescription &desc);
0051
0052 private:
0053 const edm::InputTag hepMCLabel_;
0054 const edm::InputTag beamSpotLabel_;
0055
0056 double badPull_;
0057 double longLivedDecayLength_;
0058 double vertexClusteringSqDistance_;
0059 unsigned int numberOfInnerLayers_;
0060 unsigned int minTrackerSimHits_;
0061
0062 TrackHistory tracer_;
0063
0064 TrackQuality quality_;
0065
0066 const G4toCMSLegacyProcTypeMap g4toCMSProcMap_;
0067
0068 edm::ESHandle<MagneticField> magneticField_;
0069 edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magneticFieldToken_;
0070
0071 edm::Handle<edm::HepMCProduct> mcInformation_;
0072
0073 edm::ESHandle<ParticleDataTable> particleDataTable_;
0074 edm::ESGetToken<ParticleDataTable, PDTRecord> particleDataTableToken_;
0075
0076 edm::ESHandle<TransientTrackBuilder> transientTrackBuilder_;
0077 edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> transientTrackBuilderToken_;
0078
0079 edm::Handle<reco::BeamSpot> beamSpot_;
0080
0081 const TrackerTopology *tTopo_;
0082 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoHandToken_;
0083
0084
0085
0086 void reconstructionInformation(reco::TrackBaseRef const &);
0087
0088
0089 void simulationInformation();
0090
0091
0092 void qualityInformation(reco::TrackBaseRef const &);
0093
0094
0095 void hadronFlavor();
0096
0097
0098 void processesAtGenerator();
0099
0100
0101 void processesAtSimulation();
0102
0103
0104 void vertexInformation();
0105
0106
0107 struct GeneratedPrimaryVertex {
0108 GeneratedPrimaryVertex(double x1, double y1, double z1) : x(x1), y(y1), z(z1), ptsq(0), nGenTrk(0) {}
0109
0110 bool operator<(GeneratedPrimaryVertex const &reference) const { return ptsq < reference.ptsq; }
0111
0112 double x, y, z;
0113 double ptsq;
0114 int nGenTrk;
0115
0116 HepMC::FourVector ptot;
0117
0118 std::vector<int> finalstateParticles;
0119 std::vector<int> simTrackIndex;
0120 std::vector<int> genVertex;
0121 };
0122
0123 std::vector<GeneratedPrimaryVertex> genpvs_;
0124
0125
0126 bool isFinalstateParticle(const HepMC::GenParticle *);
0127 bool isCharged(const HepMC::GenParticle *);
0128 void genPrimaryVertices();
0129 };
0130
0131 #endif