File indexing completed on 2025-06-03 00:12:26
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "FWCore/PluginManager/interface/ModuleDef.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0012 #include "FWCore/Framework/interface/Run.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/ESHandle.h"
0015 #include "FWCore/Framework/interface/EventSetup.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/Utilities/interface/InputTag.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0019 #include "FWCore/ServiceRegistry/interface/Service.h"
0020
0021
0022
0023 #include "DataFormats/Common/interface/Handle.h"
0024 #include "DataFormats/Common/interface/Ref.h"
0025
0026 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0027 #include "DataFormats/L1TrackTrigger/interface/TTCluster.h"
0028 #include "DataFormats/L1TrackTrigger/interface/TTStub.h"
0029 #include "DataFormats/L1TrackTrigger/interface/TTTrack.h"
0030 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0031 #include "SimDataFormats/TrackingAnalysis/interface/TrackingVertex.h"
0032 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0033 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0034 #include "SimDataFormats/Associations/interface/TTClusterAssociationMap.h"
0035 #include "SimDataFormats/Associations/interface/TTStubAssociationMap.h"
0036 #include "SimDataFormats/Associations/interface/TTTrackAssociationMap.h"
0037 #include "Geometry/Records/interface/StackedTrackerGeometryRecord.h"
0038
0039 #include "DataFormats/JetReco/interface/GenJetCollection.h"
0040 #include "DataFormats/JetReco/interface/GenJet.h"
0041
0042
0043
0044 #include "MagneticField/Engine/interface/MagneticField.h"
0045 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0046 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0047 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0048 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0049 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0050
0051 #include "Geometry/CommonTopologies/interface/PixelGeomDetUnit.h"
0052 #include "Geometry/CommonTopologies/interface/PixelGeomDetType.h"
0053 #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyBuilder.h"
0054 #include "Geometry/Records/interface/StackedTrackerGeometryRecord.h"
0055
0056
0057
0058 #include "L1Trigger/TrackTrigger/interface/Setup.h"
0059 #include "L1Trigger/TrackerTFP/interface/LayerEncoding.h"
0060 #include "L1Trigger/TrackFindingTracklet/interface/HitPatternHelper.h"
0061 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0062 #include "CLHEP/Units/PhysicalConstants.h"
0063
0064
0065
0066 #include <TROOT.h>
0067 #include <TCanvas.h>
0068 #include <TTree.h>
0069 #include <TFile.h>
0070 #include <TMath.h>
0071 #include <TF1.h>
0072 #include <TH2F.h>
0073 #include <TH1F.h>
0074
0075
0076
0077 #include <memory>
0078 #include <string>
0079 #include <iostream>
0080
0081
0082
0083 using namespace std;
0084 using namespace edm;
0085
0086
0087
0088
0089
0090
0091
0092 class L1TrackNtupleMaker : public one::EDAnalyzer<one::WatchRuns, one::SharedResources> {
0093 public:
0094
0095 explicit L1TrackNtupleMaker(const edm::ParameterSet& iConfig);
0096 ~L1TrackNtupleMaker() override;
0097
0098
0099 void beginJob() override;
0100 void endJob() override;
0101 void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0102 void beginRun(const Run& iEvent, const EventSetup& iSetup) override {}
0103 void endRun(const Run& iEvent, const EventSetup& iSetup) override {}
0104
0105 protected:
0106 private:
0107
0108
0109 edm::ParameterSet config;
0110
0111 int MyProcess;
0112 bool DebugMode;
0113 bool SaveAllTracks;
0114 bool SaveStubs;
0115 int L1Tk_nPar;
0116 int TP_minNStub;
0117 int TP_minNStubLayer;
0118 double TP_minPt;
0119 double TP_maxEta;
0120 double TP_maxZ0;
0121 int L1Tk_minNStub;
0122
0123 bool TrackingInJets;
0124
0125 edm::InputTag L1TrackInputTag;
0126 edm::InputTag MCTruthTrackInputTag;
0127 edm::InputTag MCTruthClusterInputTag;
0128 edm::InputTag L1StubInputTag;
0129 edm::InputTag MCTruthStubInputTag;
0130 edm::InputTag TrackingParticleInputTag;
0131 edm::InputTag TrackingVertexInputTag;
0132 edm::InputTag GenJetInputTag;
0133
0134 edm::EDGetTokenT<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>> ttClusterToken_;
0135 edm::EDGetTokenT<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>> ttStubToken_;
0136 edm::EDGetTokenT<TTClusterAssociationMap<Ref_Phase2TrackerDigi_>> ttClusterMCTruthToken_;
0137 edm::EDGetTokenT<TTStubAssociationMap<Ref_Phase2TrackerDigi_>> ttStubMCTruthToken_;
0138
0139 edm::EDGetTokenT<std::vector<TTTrack<Ref_Phase2TrackerDigi_>>> ttTrackToken_;
0140 edm::EDGetTokenT<TTTrackAssociationMap<Ref_Phase2TrackerDigi_>> ttTrackMCTruthToken_;
0141
0142 edm::EDGetTokenT<std::vector<TrackingParticle>> TrackingParticleToken_;
0143 edm::EDGetTokenT<std::vector<TrackingVertex>> TrackingVertexToken_;
0144
0145 edm::EDGetTokenT<std::vector<reco::GenJet>> GenJetToken_;
0146
0147 edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> getTokenTrackerGeom_;
0148 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> getTokenTrackerTopo_;
0149 edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> getTokenBField_;
0150 edm::ESGetToken<hph::Setup, hph::SetupRcd> getTokenHPHSetup_;
0151 edm::ESGetToken<tt::Setup, tt::SetupRcd> getTokenSetup_;
0152 edm::ESGetToken<trackerTFP::LayerEncoding, trackerTFP::DataFormatsRcd> getTokenLayerEncoding_;
0153
0154
0155
0156 bool available_;
0157
0158 TTree* eventTree;
0159
0160
0161 std::vector<float>* m_trk_pt;
0162 std::vector<float>* m_trk_eta;
0163 std::vector<float>* m_trk_phi;
0164 std::vector<float>* m_trk_d0;
0165 std::vector<float>* m_trk_z0;
0166 std::vector<float>* m_trk_chi2;
0167 std::vector<float>* m_trk_chi2_dof;
0168 std::vector<float>* m_trk_chi2rphi;
0169 std::vector<float>* m_trk_chi2rphi_dof;
0170 std::vector<float>* m_trk_chi2rz;
0171 std::vector<float>* m_trk_chi2rz_dof;
0172 std::vector<float>* m_trk_bendchi2;
0173 std::vector<int>* m_trk_nstub;
0174 std::vector<int>* m_trk_lhits;
0175 std::vector<int>* m_trk_dhits;
0176 std::vector<int>* m_trk_seed;
0177 std::vector<int>* m_trk_hitpattern;
0178 std::vector<int>* m_trk_lhits_hitpattern;
0179 std::vector<int>* m_trk_dhits_hitpattern;
0180 std::vector<int>* m_trk_nPSstub_hitpattern;
0181 std::vector<int>* m_trk_n2Sstub_hitpattern;
0182 std::vector<int>* m_trk_nLostPSstub_hitpattern;
0183 std::vector<int>* m_trk_nLost2Sstub_hitpattern;
0184 std::vector<int>* m_trk_nLoststub_V1_hitpattern;
0185 std::vector<int>* m_trk_nLoststub_V2_hitpattern;
0186 std::vector<int>* m_trk_charge;
0187 std::vector<unsigned int>* m_trk_phiSector;
0188 std::vector<int>* m_trk_etaSector;
0189 std::vector<int>* m_trk_genuine;
0190 std::vector<int>* m_trk_loose;
0191 std::vector<int>* m_trk_unknown;
0192 std::vector<int>* m_trk_combinatoric;
0193 std::vector<int>* m_trk_fake;
0194 std::vector<float>* m_trk_MVA1;
0195 std::vector<int>* m_trk_matchtp_pdgid;
0196 std::vector<float>* m_trk_matchtp_pt;
0197 std::vector<float>* m_trk_matchtp_eta;
0198 std::vector<float>* m_trk_matchtp_phi;
0199 std::vector<float>* m_trk_matchtp_z0;
0200 std::vector<float>* m_trk_matchtp_lxy;
0201 std::vector<float>* m_trk_matchtp_d0;
0202 std::vector<int>* m_trk_injet;
0203 std::vector<int>* m_trk_injet_highpt;
0204 std::vector<int>* m_trk_injet_vhighpt;
0205 std::vector<std::vector<int>>* m_trk_layers;
0206
0207
0208 std::vector<float>* m_tp_pt;
0209 std::vector<float>* m_tp_eta;
0210 std::vector<float>* m_tp_phi;
0211 std::vector<float>* m_tp_lxy;
0212 std::vector<float>* m_tp_d0;
0213 std::vector<float>* m_tp_z0;
0214 std::vector<float>* m_tp_d0_prod;
0215 std::vector<float>* m_tp_z0_prod;
0216 std::vector<int>* m_tp_pdgid;
0217 std::vector<int>* m_tp_nmatch;
0218 std::vector<int>* m_tp_nstub;
0219 std::vector<int>* m_tp_eventid;
0220 std::vector<int>* m_tp_charge;
0221 std::vector<int>* m_tp_injet;
0222 std::vector<int>* m_tp_injet_highpt;
0223 std::vector<int>* m_tp_injet_vhighpt;
0224
0225
0226 std::vector<float>* m_matchtrk_pt;
0227 std::vector<float>* m_matchtrk_eta;
0228 std::vector<float>* m_matchtrk_phi;
0229 std::vector<float>* m_matchtrk_d0;
0230 std::vector<float>* m_matchtrk_z0;
0231 std::vector<float>* m_matchtrk_chi2;
0232 std::vector<float>* m_matchtrk_chi2_dof;
0233 std::vector<float>* m_matchtrk_chi2rphi;
0234 std::vector<float>* m_matchtrk_chi2rphi_dof;
0235 std::vector<float>* m_matchtrk_chi2rz;
0236 std::vector<float>* m_matchtrk_chi2rz_dof;
0237 std::vector<float>* m_matchtrk_bendchi2;
0238 std::vector<float>* m_matchtrk_MVA1;
0239 std::vector<int>* m_matchtrk_nstub;
0240 std::vector<int>* m_matchtrk_lhits;
0241 std::vector<int>* m_matchtrk_dhits;
0242 std::vector<int>* m_matchtrk_seed;
0243 std::vector<int>* m_matchtrk_hitpattern;
0244 std::vector<int>* m_matchtrk_charge;
0245 std::vector<int>* m_matchtrk_injet;
0246 std::vector<int>* m_matchtrk_injet_highpt;
0247 std::vector<int>* m_matchtrk_injet_vhighpt;
0248
0249
0250 std::vector<float>* m_allstub_x;
0251 std::vector<float>* m_allstub_y;
0252 std::vector<float>* m_allstub_z;
0253
0254 std::vector<int>* m_allstub_isBarrel;
0255 std::vector<int>* m_allstub_layer;
0256 std::vector<int>* m_allstub_isPSmodule;
0257 std::vector<int>* m_allstub_isTiltedBarrel;
0258
0259 std::vector<float>* m_allstub_trigDisplace;
0260 std::vector<float>* m_allstub_trigOffset;
0261 std::vector<float>* m_allstub_trigPos;
0262 std::vector<float>* m_allstub_trigBend;
0263
0264
0265 std::vector<int>* m_allstub_matchTP_pdgid;
0266 std::vector<float>* m_allstub_matchTP_pt;
0267 std::vector<float>* m_allstub_matchTP_eta;
0268 std::vector<float>* m_allstub_matchTP_phi;
0269
0270 std::vector<int>* m_allstub_genuine;
0271
0272
0273 std::vector<float>* m_jet_eta;
0274 std::vector<float>* m_jet_phi;
0275 std::vector<float>* m_jet_pt;
0276 std::vector<float>* m_jet_tp_sumpt;
0277 std::vector<float>* m_jet_trk_sumpt;
0278 std::vector<float>* m_jet_matchtrk_sumpt;
0279 };
0280
0281
0282
0283
0284
0285
0286
0287
0288
0289 L1TrackNtupleMaker::L1TrackNtupleMaker(edm::ParameterSet const& iConfig) : config(iConfig) {
0290 usesResource("TFileService");
0291 MyProcess = iConfig.getParameter<int>("MyProcess");
0292 DebugMode = iConfig.getParameter<bool>("DebugMode");
0293 SaveAllTracks = iConfig.getParameter<bool>("SaveAllTracks");
0294 SaveStubs = iConfig.getParameter<bool>("SaveStubs");
0295 L1Tk_nPar = iConfig.getParameter<int>("L1Tk_nPar");
0296 TP_minNStub = iConfig.getParameter<int>("TP_minNStub");
0297 TP_minNStubLayer = iConfig.getParameter<int>("TP_minNStubLayer");
0298 TP_minPt = iConfig.getParameter<double>("TP_minPt");
0299 TP_maxEta = iConfig.getParameter<double>("TP_maxEta");
0300 TP_maxZ0 = iConfig.getParameter<double>("TP_maxZ0");
0301 L1TrackInputTag = iConfig.getParameter<edm::InputTag>("L1TrackInputTag");
0302 MCTruthTrackInputTag = iConfig.getParameter<edm::InputTag>("MCTruthTrackInputTag");
0303 L1Tk_minNStub = iConfig.getParameter<int>("L1Tk_minNStub");
0304
0305 TrackingInJets = iConfig.getParameter<bool>("TrackingInJets");
0306
0307 L1StubInputTag = iConfig.getParameter<edm::InputTag>("L1StubInputTag");
0308 MCTruthClusterInputTag = iConfig.getParameter<edm::InputTag>("MCTruthClusterInputTag");
0309 MCTruthStubInputTag = iConfig.getParameter<edm::InputTag>("MCTruthStubInputTag");
0310 TrackingParticleInputTag = iConfig.getParameter<edm::InputTag>("TrackingParticleInputTag");
0311 TrackingVertexInputTag = iConfig.getParameter<edm::InputTag>("TrackingVertexInputTag");
0312 GenJetInputTag = iConfig.getParameter<edm::InputTag>("GenJetInputTag");
0313
0314 ttTrackToken_ = consumes<std::vector<TTTrack<Ref_Phase2TrackerDigi_>>>(L1TrackInputTag);
0315 ttTrackMCTruthToken_ = consumes<TTTrackAssociationMap<Ref_Phase2TrackerDigi_>>(MCTruthTrackInputTag);
0316 ttStubToken_ = consumes<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>>(L1StubInputTag);
0317 ttClusterMCTruthToken_ = consumes<TTClusterAssociationMap<Ref_Phase2TrackerDigi_>>(MCTruthClusterInputTag);
0318 ttStubMCTruthToken_ = consumes<TTStubAssociationMap<Ref_Phase2TrackerDigi_>>(MCTruthStubInputTag);
0319
0320 TrackingParticleToken_ = consumes<std::vector<TrackingParticle>>(TrackingParticleInputTag);
0321 TrackingVertexToken_ = consumes<std::vector<TrackingVertex>>(TrackingVertexInputTag);
0322 GenJetToken_ = consumes<std::vector<reco::GenJet>>(GenJetInputTag);
0323
0324 getTokenTrackerGeom_ = esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>();
0325 getTokenTrackerTopo_ = esConsumes<TrackerTopology, TrackerTopologyRcd>();
0326 getTokenBField_ = esConsumes<MagneticField, IdealMagneticFieldRecord>();
0327 getTokenHPHSetup_ = esConsumes<hph::Setup, hph::SetupRcd>();
0328 getTokenSetup_ = esConsumes<tt::Setup, tt::SetupRcd>();
0329 getTokenLayerEncoding_ = esConsumes<trackerTFP::LayerEncoding, trackerTFP::DataFormatsRcd>();
0330 }
0331
0332
0333
0334 L1TrackNtupleMaker::~L1TrackNtupleMaker() {}
0335
0336
0337
0338 void L1TrackNtupleMaker::endJob() {
0339
0340 edm::LogVerbatim("Tracklet") << "L1TrackNtupleMaker::endJob";
0341
0342
0343 delete m_trk_pt;
0344 delete m_trk_eta;
0345 delete m_trk_phi;
0346 delete m_trk_z0;
0347 delete m_trk_d0;
0348 delete m_trk_chi2;
0349 delete m_trk_chi2_dof;
0350 delete m_trk_chi2rphi;
0351 delete m_trk_chi2rphi_dof;
0352 delete m_trk_chi2rz;
0353 delete m_trk_chi2rz_dof;
0354 delete m_trk_bendchi2;
0355 delete m_trk_nstub;
0356 delete m_trk_lhits;
0357 delete m_trk_dhits;
0358 delete m_trk_seed;
0359 delete m_trk_hitpattern;
0360 delete m_trk_lhits_hitpattern;
0361 delete m_trk_dhits_hitpattern;
0362 delete m_trk_nPSstub_hitpattern;
0363 delete m_trk_n2Sstub_hitpattern;
0364 delete m_trk_nLostPSstub_hitpattern;
0365 delete m_trk_nLost2Sstub_hitpattern;
0366 delete m_trk_nLoststub_V1_hitpattern;
0367 delete m_trk_nLoststub_V2_hitpattern;
0368 delete m_trk_charge;
0369 delete m_trk_phiSector;
0370 delete m_trk_etaSector;
0371 delete m_trk_genuine;
0372 delete m_trk_loose;
0373 delete m_trk_unknown;
0374 delete m_trk_combinatoric;
0375 delete m_trk_fake;
0376 delete m_trk_MVA1;
0377 delete m_trk_matchtp_pdgid;
0378 delete m_trk_matchtp_pt;
0379 delete m_trk_matchtp_eta;
0380 delete m_trk_matchtp_phi;
0381 delete m_trk_matchtp_z0;
0382 delete m_trk_matchtp_lxy;
0383 delete m_trk_matchtp_d0;
0384 delete m_trk_injet;
0385 delete m_trk_injet_highpt;
0386 delete m_trk_injet_vhighpt;
0387
0388 delete m_tp_pt;
0389 delete m_tp_eta;
0390 delete m_tp_phi;
0391 delete m_tp_lxy;
0392 delete m_tp_d0;
0393 delete m_tp_z0;
0394 delete m_tp_d0_prod;
0395 delete m_tp_z0_prod;
0396 delete m_tp_pdgid;
0397 delete m_tp_nmatch;
0398 delete m_tp_nstub;
0399 delete m_tp_eventid;
0400 delete m_tp_charge;
0401 delete m_tp_injet;
0402 delete m_tp_injet_highpt;
0403 delete m_tp_injet_vhighpt;
0404
0405 delete m_matchtrk_pt;
0406 delete m_matchtrk_eta;
0407 delete m_matchtrk_phi;
0408 delete m_matchtrk_z0;
0409 delete m_matchtrk_d0;
0410 delete m_matchtrk_chi2;
0411 delete m_matchtrk_chi2_dof;
0412 delete m_matchtrk_chi2rphi;
0413 delete m_matchtrk_chi2rphi_dof;
0414 delete m_matchtrk_chi2rz;
0415 delete m_matchtrk_chi2rz_dof;
0416 delete m_matchtrk_bendchi2;
0417 delete m_matchtrk_MVA1;
0418 delete m_matchtrk_nstub;
0419 delete m_matchtrk_dhits;
0420 delete m_matchtrk_lhits;
0421 delete m_matchtrk_seed;
0422 delete m_matchtrk_hitpattern;
0423 delete m_matchtrk_charge;
0424 delete m_matchtrk_injet;
0425 delete m_matchtrk_injet_highpt;
0426 delete m_matchtrk_injet_vhighpt;
0427
0428 delete m_allstub_x;
0429 delete m_allstub_y;
0430 delete m_allstub_z;
0431 delete m_allstub_isBarrel;
0432 delete m_allstub_layer;
0433 delete m_allstub_isPSmodule;
0434 delete m_allstub_trigDisplace;
0435 delete m_allstub_trigOffset;
0436 delete m_allstub_trigPos;
0437 delete m_allstub_trigBend;
0438 delete m_allstub_matchTP_pdgid;
0439 delete m_allstub_matchTP_pt;
0440 delete m_allstub_matchTP_eta;
0441 delete m_allstub_matchTP_phi;
0442 delete m_allstub_genuine;
0443
0444 delete m_jet_eta;
0445 delete m_jet_phi;
0446 delete m_jet_pt;
0447 delete m_jet_tp_sumpt;
0448 delete m_jet_trk_sumpt;
0449 delete m_jet_matchtrk_sumpt;
0450 }
0451
0452
0453
0454 void L1TrackNtupleMaker::beginJob() {
0455
0456 edm::LogVerbatim("Tracklet") << "L1TrackNtupleMaker::beginJob";
0457
0458
0459
0460 edm::Service<TFileService> fs;
0461 available_ = fs.isAvailable();
0462 if (not available_)
0463 return;
0464
0465
0466 m_trk_pt = new std::vector<float>;
0467 m_trk_eta = new std::vector<float>;
0468 m_trk_phi = new std::vector<float>;
0469 m_trk_z0 = new std::vector<float>;
0470 m_trk_d0 = new std::vector<float>;
0471 m_trk_chi2 = new std::vector<float>;
0472 m_trk_chi2_dof = new std::vector<float>;
0473 m_trk_chi2rphi = new std::vector<float>;
0474 m_trk_chi2rphi_dof = new std::vector<float>;
0475 m_trk_chi2rz = new std::vector<float>;
0476 m_trk_chi2rz_dof = new std::vector<float>;
0477 m_trk_bendchi2 = new std::vector<float>;
0478 m_trk_nstub = new std::vector<int>;
0479 m_trk_lhits = new std::vector<int>;
0480 m_trk_dhits = new std::vector<int>;
0481 m_trk_seed = new std::vector<int>;
0482 m_trk_hitpattern = new std::vector<int>;
0483 m_trk_lhits_hitpattern = new std::vector<int>;
0484 m_trk_dhits_hitpattern = new std::vector<int>;
0485 m_trk_nPSstub_hitpattern = new std::vector<int>;
0486 m_trk_n2Sstub_hitpattern = new std::vector<int>;
0487 m_trk_nLostPSstub_hitpattern = new std::vector<int>;
0488 m_trk_nLost2Sstub_hitpattern = new std::vector<int>;
0489 m_trk_nLoststub_V1_hitpattern = new std::vector<int>;
0490 m_trk_nLoststub_V2_hitpattern = new std::vector<int>;
0491 m_trk_charge = new std::vector<int>;
0492 m_trk_phiSector = new std::vector<unsigned int>;
0493 m_trk_etaSector = new std::vector<int>;
0494 m_trk_genuine = new std::vector<int>;
0495 m_trk_loose = new std::vector<int>;
0496 m_trk_unknown = new std::vector<int>;
0497 m_trk_combinatoric = new std::vector<int>;
0498 m_trk_fake = new std::vector<int>;
0499 m_trk_MVA1 = new std::vector<float>;
0500 m_trk_matchtp_pdgid = new std::vector<int>;
0501 m_trk_matchtp_pt = new std::vector<float>;
0502 m_trk_matchtp_eta = new std::vector<float>;
0503 m_trk_matchtp_phi = new std::vector<float>;
0504 m_trk_matchtp_z0 = new std::vector<float>;
0505 m_trk_matchtp_lxy = new std::vector<float>;
0506 m_trk_matchtp_d0 = new std::vector<float>;
0507 m_trk_injet = new std::vector<int>;
0508 m_trk_injet_highpt = new std::vector<int>;
0509 m_trk_injet_vhighpt = new std::vector<int>;
0510 m_trk_layers = new std::vector<std::vector<int>>;
0511
0512 m_tp_pt = new std::vector<float>;
0513 m_tp_eta = new std::vector<float>;
0514 m_tp_phi = new std::vector<float>;
0515 m_tp_lxy = new std::vector<float>;
0516 m_tp_d0 = new std::vector<float>;
0517 m_tp_z0 = new std::vector<float>;
0518 m_tp_d0_prod = new std::vector<float>;
0519 m_tp_z0_prod = new std::vector<float>;
0520 m_tp_pdgid = new std::vector<int>;
0521 m_tp_nmatch = new std::vector<int>;
0522 m_tp_nstub = new std::vector<int>;
0523 m_tp_eventid = new std::vector<int>;
0524 m_tp_charge = new std::vector<int>;
0525 m_tp_injet = new std::vector<int>;
0526 m_tp_injet_highpt = new std::vector<int>;
0527 m_tp_injet_vhighpt = new std::vector<int>;
0528
0529 m_matchtrk_pt = new std::vector<float>;
0530 m_matchtrk_eta = new std::vector<float>;
0531 m_matchtrk_phi = new std::vector<float>;
0532 m_matchtrk_z0 = new std::vector<float>;
0533 m_matchtrk_d0 = new std::vector<float>;
0534 m_matchtrk_chi2 = new std::vector<float>;
0535 m_matchtrk_chi2_dof = new std::vector<float>;
0536 m_matchtrk_chi2rphi = new std::vector<float>;
0537 m_matchtrk_chi2rphi_dof = new std::vector<float>;
0538 m_matchtrk_chi2rz = new std::vector<float>;
0539 m_matchtrk_chi2rz_dof = new std::vector<float>;
0540 m_matchtrk_bendchi2 = new std::vector<float>;
0541 m_matchtrk_MVA1 = new std::vector<float>;
0542 m_matchtrk_nstub = new std::vector<int>;
0543 m_matchtrk_dhits = new std::vector<int>;
0544 m_matchtrk_lhits = new std::vector<int>;
0545 m_matchtrk_seed = new std::vector<int>;
0546 m_matchtrk_hitpattern = new std::vector<int>;
0547 m_matchtrk_charge = new std::vector<int>;
0548 m_matchtrk_injet = new std::vector<int>;
0549 m_matchtrk_injet_highpt = new std::vector<int>;
0550 m_matchtrk_injet_vhighpt = new std::vector<int>;
0551
0552 m_allstub_x = new std::vector<float>;
0553 m_allstub_y = new std::vector<float>;
0554 m_allstub_z = new std::vector<float>;
0555
0556 m_allstub_isBarrel = new std::vector<int>;
0557 m_allstub_layer = new std::vector<int>;
0558 m_allstub_isPSmodule = new std::vector<int>;
0559 m_allstub_isTiltedBarrel = new std::vector<int>;
0560 m_allstub_trigDisplace = new std::vector<float>;
0561 m_allstub_trigOffset = new std::vector<float>;
0562 m_allstub_trigPos = new std::vector<float>;
0563 m_allstub_trigBend = new std::vector<float>;
0564
0565 m_allstub_matchTP_pdgid = new std::vector<int>;
0566 m_allstub_matchTP_pt = new std::vector<float>;
0567 m_allstub_matchTP_eta = new std::vector<float>;
0568 m_allstub_matchTP_phi = new std::vector<float>;
0569
0570 m_allstub_genuine = new std::vector<int>;
0571
0572 m_jet_eta = new std::vector<float>;
0573 m_jet_phi = new std::vector<float>;
0574 m_jet_pt = new std::vector<float>;
0575 m_jet_tp_sumpt = new std::vector<float>;
0576 m_jet_trk_sumpt = new std::vector<float>;
0577 m_jet_matchtrk_sumpt = new std::vector<float>;
0578
0579
0580 eventTree = fs->make<TTree>("eventTree", "Event tree");
0581
0582 if (SaveAllTracks) {
0583 eventTree->Branch("trk_pt", &m_trk_pt);
0584 eventTree->Branch("trk_eta", &m_trk_eta);
0585 eventTree->Branch("trk_phi", &m_trk_phi);
0586 eventTree->Branch("trk_d0", &m_trk_d0);
0587 eventTree->Branch("trk_z0", &m_trk_z0);
0588 eventTree->Branch("trk_chi2", &m_trk_chi2);
0589 eventTree->Branch("trk_chi2_dof", &m_trk_chi2_dof);
0590 eventTree->Branch("trk_chi2rphi", &m_trk_chi2rphi);
0591 eventTree->Branch("trk_chi2rphi_dof", &m_trk_chi2rphi_dof);
0592 eventTree->Branch("trk_chi2rz", &m_trk_chi2rz);
0593 eventTree->Branch("trk_chi2rz_dof", &m_trk_chi2rz_dof);
0594 eventTree->Branch("trk_bendchi2", &m_trk_bendchi2);
0595 eventTree->Branch("trk_nstub", &m_trk_nstub);
0596 eventTree->Branch("trk_lhits", &m_trk_lhits);
0597 eventTree->Branch("trk_dhits", &m_trk_dhits);
0598 eventTree->Branch("trk_seed", &m_trk_seed);
0599 eventTree->Branch("trk_hitpattern", &m_trk_hitpattern);
0600 eventTree->Branch("trk_lhits_hitpattern", &m_trk_lhits_hitpattern);
0601 eventTree->Branch("trk_dhits_hitpattern", &m_trk_dhits_hitpattern);
0602 eventTree->Branch("trk_nPSstub_hitpattern", &m_trk_nPSstub_hitpattern);
0603 eventTree->Branch("trk_n2Sstub_hitpattern", &m_trk_n2Sstub_hitpattern);
0604 eventTree->Branch("trk_nLostPSstub_hitpattern", &m_trk_nLostPSstub_hitpattern);
0605 eventTree->Branch("trk_nLost2Sstub_hitpattern", &m_trk_nLost2Sstub_hitpattern);
0606 eventTree->Branch("trk_nLoststub_V1_hitpattern", &m_trk_nLoststub_V1_hitpattern);
0607 eventTree->Branch("trk_nLoststub_V2_hitpattern", &m_trk_nLoststub_V2_hitpattern);
0608 eventTree->Branch("trk_charge", &m_trk_charge);
0609 eventTree->Branch("trk_phiSector", &m_trk_phiSector);
0610 eventTree->Branch("trk_etaSector", &m_trk_etaSector);
0611 eventTree->Branch("trk_genuine", &m_trk_genuine);
0612 eventTree->Branch("trk_loose", &m_trk_loose);
0613 eventTree->Branch("trk_unknown", &m_trk_unknown);
0614 eventTree->Branch("trk_combinatoric", &m_trk_combinatoric);
0615 eventTree->Branch("trk_fake", &m_trk_fake);
0616 eventTree->Branch("trk_MVA1", &m_trk_MVA1);
0617 eventTree->Branch("trk_matchtp_pdgid", &m_trk_matchtp_pdgid);
0618 eventTree->Branch("trk_matchtp_pt", &m_trk_matchtp_pt);
0619 eventTree->Branch("trk_matchtp_eta", &m_trk_matchtp_eta);
0620 eventTree->Branch("trk_matchtp_phi", &m_trk_matchtp_phi);
0621 eventTree->Branch("trk_matchtp_z0", &m_trk_matchtp_z0);
0622 eventTree->Branch("trk_matchtp_lxy", &m_trk_matchtp_lxy);
0623 eventTree->Branch("trk_matchtp_d0", &m_trk_matchtp_d0);
0624 if (TrackingInJets) {
0625 eventTree->Branch("trk_injet", &m_trk_injet);
0626 eventTree->Branch("trk_injet_highpt", &m_trk_injet_highpt);
0627 eventTree->Branch("trk_injet_vhighpt", &m_trk_injet_vhighpt);
0628 }
0629 eventTree->Branch("m_trk_layers", &m_trk_layers);
0630 }
0631
0632 eventTree->Branch("tp_pt", &m_tp_pt);
0633 eventTree->Branch("tp_eta", &m_tp_eta);
0634 eventTree->Branch("tp_phi", &m_tp_phi);
0635 eventTree->Branch("tp_lxy", &m_tp_lxy);
0636 eventTree->Branch("tp_d0", &m_tp_d0);
0637 eventTree->Branch("tp_z0", &m_tp_z0);
0638 eventTree->Branch("tp_d0_prod", &m_tp_d0_prod);
0639 eventTree->Branch("tp_z0_prod", &m_tp_z0_prod);
0640 eventTree->Branch("tp_pdgid", &m_tp_pdgid);
0641 eventTree->Branch("tp_nmatch", &m_tp_nmatch);
0642 eventTree->Branch("tp_nstub", &m_tp_nstub);
0643 eventTree->Branch("tp_eventid", &m_tp_eventid);
0644 eventTree->Branch("tp_charge", &m_tp_charge);
0645 if (TrackingInJets) {
0646 eventTree->Branch("tp_injet", &m_tp_injet);
0647 eventTree->Branch("tp_injet_highpt", &m_tp_injet_highpt);
0648 eventTree->Branch("tp_injet_vhighpt", &m_tp_injet_vhighpt);
0649 }
0650
0651 eventTree->Branch("matchtrk_pt", &m_matchtrk_pt);
0652 eventTree->Branch("matchtrk_eta", &m_matchtrk_eta);
0653 eventTree->Branch("matchtrk_phi", &m_matchtrk_phi);
0654 eventTree->Branch("matchtrk_z0", &m_matchtrk_z0);
0655 eventTree->Branch("matchtrk_d0", &m_matchtrk_d0);
0656 eventTree->Branch("matchtrk_chi2", &m_matchtrk_chi2);
0657 eventTree->Branch("matchtrk_chi2_dof", &m_matchtrk_chi2_dof);
0658 eventTree->Branch("matchtrk_chi2rphi", &m_matchtrk_chi2rphi);
0659 eventTree->Branch("matchtrk_chi2rphi_dof", &m_matchtrk_chi2rphi_dof);
0660 eventTree->Branch("matchtrk_chi2rz", &m_matchtrk_chi2rz);
0661 eventTree->Branch("matchtrk_chi2rz_dof", &m_matchtrk_chi2rz_dof);
0662 eventTree->Branch("matchtrk_bendchi2", &m_matchtrk_bendchi2);
0663 eventTree->Branch("matchtrk_MVA1", &m_matchtrk_MVA1);
0664 eventTree->Branch("matchtrk_nstub", &m_matchtrk_nstub);
0665 eventTree->Branch("matchtrk_lhits", &m_matchtrk_lhits);
0666 eventTree->Branch("matchtrk_dhits", &m_matchtrk_dhits);
0667 eventTree->Branch("matchtrk_seed", &m_matchtrk_seed);
0668 eventTree->Branch("matchtrk_hitpattern", &m_matchtrk_hitpattern);
0669 eventTree->Branch("matchtrk_charge", &m_matchtrk_charge);
0670 if (TrackingInJets) {
0671 eventTree->Branch("matchtrk_injet", &m_matchtrk_injet);
0672 eventTree->Branch("matchtrk_injet_highpt", &m_matchtrk_injet_highpt);
0673 eventTree->Branch("matchtrk_injet_vhighpt", &m_matchtrk_injet_vhighpt);
0674 }
0675
0676 if (SaveStubs) {
0677 eventTree->Branch("allstub_x", &m_allstub_x);
0678 eventTree->Branch("allstub_y", &m_allstub_y);
0679 eventTree->Branch("allstub_z", &m_allstub_z);
0680
0681 eventTree->Branch("allstub_isBarrel", &m_allstub_isBarrel);
0682 eventTree->Branch("allstub_layer", &m_allstub_layer);
0683 eventTree->Branch("allstub_isPSmodule", &m_allstub_isPSmodule);
0684 eventTree->Branch("allstub_isTiltedBarrel", &m_allstub_isTiltedBarrel);
0685
0686 eventTree->Branch("allstub_trigDisplace", &m_allstub_trigDisplace);
0687 eventTree->Branch("allstub_trigOffset", &m_allstub_trigOffset);
0688 eventTree->Branch("allstub_trigPos", &m_allstub_trigPos);
0689 eventTree->Branch("allstub_trigBend", &m_allstub_trigBend);
0690
0691 eventTree->Branch("allstub_matchTP_pdgid", &m_allstub_matchTP_pdgid);
0692 eventTree->Branch("allstub_matchTP_pt", &m_allstub_matchTP_pt);
0693 eventTree->Branch("allstub_matchTP_eta", &m_allstub_matchTP_eta);
0694 eventTree->Branch("allstub_matchTP_phi", &m_allstub_matchTP_phi);
0695
0696 eventTree->Branch("allstub_genuine", &m_allstub_genuine);
0697 }
0698
0699 if (TrackingInJets) {
0700 eventTree->Branch("jet_eta", &m_jet_eta);
0701 eventTree->Branch("jet_phi", &m_jet_phi);
0702 eventTree->Branch("jet_pt", &m_jet_pt);
0703 eventTree->Branch("jet_tp_sumpt", &m_jet_tp_sumpt);
0704 eventTree->Branch("jet_trk_sumpt", &m_jet_trk_sumpt);
0705 eventTree->Branch("jet_matchtrk_sumpt", &m_jet_matchtrk_sumpt);
0706 }
0707 }
0708
0709
0710
0711 void L1TrackNtupleMaker::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0712 if (not available_)
0713 return;
0714
0715 if (!(MyProcess == 13 || MyProcess == 11 || MyProcess == 211 || MyProcess == 6 || MyProcess == 15 ||
0716 MyProcess == 1)) {
0717 edm::LogVerbatim("Tracklet") << "The specified MyProcess is invalid! Exiting...";
0718 return;
0719 }
0720
0721 if (!(L1Tk_nPar == 4 || L1Tk_nPar == 5)) {
0722 edm::LogVerbatim("Tracklet") << "Invalid number of track parameters, specified L1Tk_nPar == " << L1Tk_nPar
0723 << " but only 4/5 are valid options! Exiting...";
0724 return;
0725 }
0726
0727
0728 if (SaveAllTracks) {
0729 m_trk_pt->clear();
0730 m_trk_eta->clear();
0731 m_trk_phi->clear();
0732 m_trk_d0->clear();
0733 m_trk_z0->clear();
0734 m_trk_chi2->clear();
0735 m_trk_chi2_dof->clear();
0736 m_trk_chi2rphi->clear();
0737 m_trk_chi2rphi_dof->clear();
0738 m_trk_chi2rz->clear();
0739 m_trk_chi2rz_dof->clear();
0740 m_trk_bendchi2->clear();
0741 m_trk_nstub->clear();
0742 m_trk_lhits->clear();
0743 m_trk_dhits->clear();
0744 m_trk_seed->clear();
0745 m_trk_hitpattern->clear();
0746 m_trk_lhits_hitpattern->clear();
0747 m_trk_dhits_hitpattern->clear();
0748 m_trk_nPSstub_hitpattern->clear();
0749 m_trk_n2Sstub_hitpattern->clear();
0750 m_trk_nLostPSstub_hitpattern->clear();
0751 m_trk_nLost2Sstub_hitpattern->clear();
0752 m_trk_nLoststub_V1_hitpattern->clear();
0753 m_trk_nLoststub_V2_hitpattern->clear();
0754 m_trk_charge->clear();
0755 m_trk_phiSector->clear();
0756 m_trk_etaSector->clear();
0757 m_trk_genuine->clear();
0758 m_trk_loose->clear();
0759 m_trk_unknown->clear();
0760 m_trk_combinatoric->clear();
0761 m_trk_fake->clear();
0762 m_trk_MVA1->clear();
0763 m_trk_matchtp_pdgid->clear();
0764 m_trk_matchtp_pt->clear();
0765 m_trk_matchtp_eta->clear();
0766 m_trk_matchtp_phi->clear();
0767 m_trk_matchtp_z0->clear();
0768 m_trk_matchtp_lxy->clear();
0769 m_trk_matchtp_d0->clear();
0770 m_trk_injet->clear();
0771 m_trk_injet_highpt->clear();
0772 m_trk_injet_vhighpt->clear();
0773 m_trk_layers->clear();
0774 }
0775
0776 m_tp_pt->clear();
0777 m_tp_eta->clear();
0778 m_tp_phi->clear();
0779 m_tp_lxy->clear();
0780 m_tp_d0->clear();
0781 m_tp_z0->clear();
0782 m_tp_d0_prod->clear();
0783 m_tp_z0_prod->clear();
0784 m_tp_pdgid->clear();
0785 m_tp_nmatch->clear();
0786 m_tp_nstub->clear();
0787 m_tp_eventid->clear();
0788 m_tp_charge->clear();
0789 m_tp_injet->clear();
0790 m_tp_injet_highpt->clear();
0791 m_tp_injet_vhighpt->clear();
0792
0793 m_matchtrk_pt->clear();
0794 m_matchtrk_eta->clear();
0795 m_matchtrk_phi->clear();
0796 m_matchtrk_z0->clear();
0797 m_matchtrk_d0->clear();
0798 m_matchtrk_chi2->clear();
0799 m_matchtrk_chi2_dof->clear();
0800 m_matchtrk_chi2rphi->clear();
0801 m_matchtrk_chi2rphi_dof->clear();
0802 m_matchtrk_chi2rz->clear();
0803 m_matchtrk_chi2rz_dof->clear();
0804 m_matchtrk_bendchi2->clear();
0805 m_matchtrk_MVA1->clear();
0806 m_matchtrk_nstub->clear();
0807 m_matchtrk_lhits->clear();
0808 m_matchtrk_dhits->clear();
0809 m_matchtrk_seed->clear();
0810 m_matchtrk_hitpattern->clear();
0811 m_matchtrk_charge->clear();
0812 m_matchtrk_injet->clear();
0813 m_matchtrk_injet_highpt->clear();
0814 m_matchtrk_injet_vhighpt->clear();
0815
0816 if (SaveStubs) {
0817 m_allstub_x->clear();
0818 m_allstub_y->clear();
0819 m_allstub_z->clear();
0820
0821 m_allstub_isBarrel->clear();
0822 m_allstub_layer->clear();
0823 m_allstub_isPSmodule->clear();
0824 m_allstub_isTiltedBarrel->clear();
0825
0826 m_allstub_trigDisplace->clear();
0827 m_allstub_trigOffset->clear();
0828 m_allstub_trigPos->clear();
0829 m_allstub_trigBend->clear();
0830
0831 m_allstub_matchTP_pdgid->clear();
0832 m_allstub_matchTP_pt->clear();
0833 m_allstub_matchTP_eta->clear();
0834 m_allstub_matchTP_phi->clear();
0835
0836 m_allstub_genuine->clear();
0837 }
0838
0839 m_jet_eta->clear();
0840 m_jet_phi->clear();
0841 m_jet_pt->clear();
0842 m_jet_tp_sumpt->clear();
0843 m_jet_trk_sumpt->clear();
0844 m_jet_matchtrk_sumpt->clear();
0845
0846
0847
0848
0849
0850
0851 edm::Handle<std::vector<TTTrack<Ref_Phase2TrackerDigi_>>> TTTrackHandle;
0852 iEvent.getByToken(ttTrackToken_, TTTrackHandle);
0853
0854
0855 edm::Handle<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>> TTStubHandle;
0856 if (SaveStubs)
0857 iEvent.getByToken(ttStubToken_, TTStubHandle);
0858
0859
0860 edm::Handle<TTClusterAssociationMap<Ref_Phase2TrackerDigi_>> MCTruthTTClusterHandle;
0861 iEvent.getByToken(ttClusterMCTruthToken_, MCTruthTTClusterHandle);
0862 edm::Handle<TTStubAssociationMap<Ref_Phase2TrackerDigi_>> MCTruthTTStubHandle;
0863 iEvent.getByToken(ttStubMCTruthToken_, MCTruthTTStubHandle);
0864 edm::Handle<TTTrackAssociationMap<Ref_Phase2TrackerDigi_>> MCTruthTTTrackHandle;
0865 iEvent.getByToken(ttTrackMCTruthToken_, MCTruthTTTrackHandle);
0866
0867
0868 edm::Handle<std::vector<TrackingParticle>> TrackingParticleHandle;
0869 edm::Handle<std::vector<TrackingVertex>> TrackingVertexHandle;
0870 iEvent.getByToken(TrackingParticleToken_, TrackingParticleHandle);
0871
0872
0873
0874
0875 edm::ESHandle<TrackerGeometry> tGeomHandle = iSetup.getHandle(getTokenTrackerGeom_);
0876
0877 edm::ESHandle<TrackerTopology> tTopoHandle = iSetup.getHandle(getTokenTrackerTopo_);
0878
0879 edm::ESHandle<MagneticField> bFieldHandle = iSetup.getHandle(getTokenBField_);
0880
0881 edm::ESHandle<hph::Setup> hphHandle = iSetup.getHandle(getTokenHPHSetup_);
0882 edm::ESHandle<tt::Setup> handleSetup = iSetup.getHandle(getTokenSetup_);
0883 edm::ESHandle<trackerTFP::LayerEncoding> handleLayerEncoding = iSetup.getHandle(getTokenLayerEncoding_);
0884
0885 const TrackerTopology* const tTopo = tTopoHandle.product();
0886 const TrackerGeometry* const theTrackerGeom = tGeomHandle.product();
0887 const hph::Setup* hphSetup = hphHandle.product();
0888 const tt::Setup* setup = handleSetup.product();
0889 const trackerTFP::LayerEncoding* layerEncoding = handleLayerEncoding.product();
0890
0891
0892
0893
0894
0895 if (SaveStubs) {
0896 for (auto gd = theTrackerGeom->dets().begin(); gd != theTrackerGeom->dets().end(); gd++) {
0897 DetId detid = (*gd)->geographicalId();
0898 if (detid.subdetId() != StripSubdetector::TOB && detid.subdetId() != StripSubdetector::TID)
0899 continue;
0900 if (!tTopo->isLower(detid))
0901 continue;
0902 DetId stackDetid = tTopo->stack(detid);
0903
0904 if (TTStubHandle->find(stackDetid) == TTStubHandle->end())
0905 continue;
0906
0907
0908 edmNew::DetSet<TTStub<Ref_Phase2TrackerDigi_>> stubs = (*TTStubHandle)[stackDetid];
0909 const GeomDetUnit* det0 = theTrackerGeom->idToDetUnit(detid);
0910 const auto* theGeomDet = dynamic_cast<const PixelGeomDetUnit*>(det0);
0911 const PixelTopology* topol = dynamic_cast<const PixelTopology*>(&(theGeomDet->specificTopology()));
0912
0913
0914 for (auto stubIter = stubs.begin(); stubIter != stubs.end(); ++stubIter) {
0915 edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>> tempStubPtr =
0916 edmNew::makeRefTo(TTStubHandle, stubIter);
0917
0918 int isBarrel = 0;
0919 int layer = -999999;
0920 if (detid.subdetId() == StripSubdetector::TOB) {
0921 isBarrel = 1;
0922 layer = static_cast<int>(tTopo->layer(detid));
0923 } else if (detid.subdetId() == StripSubdetector::TID) {
0924 isBarrel = 0;
0925 layer = static_cast<int>(tTopo->layer(detid));
0926 } else {
0927 edm::LogVerbatim("Tracklet") << "WARNING -- neither TOB or TID stub, shouldn't happen...";
0928 layer = -1;
0929 }
0930
0931 int isPSmodule = 0;
0932 if (topol->nrows() == 960)
0933 isPSmodule = 1;
0934
0935 const unsigned int tobSide = tTopo->tobSide(detid);
0936 int isTiltedBarrel = 0;
0937 if (isBarrel == 1 && (tobSide == 1 || tobSide == 2))
0938 isTiltedBarrel = 1;
0939
0940 MeasurementPoint coords = tempStubPtr->clusterRef(0)->findAverageLocalCoordinatesCentered();
0941 LocalPoint clustlp = topol->localPosition(coords);
0942 GlobalPoint posStub = theGeomDet->surface().toGlobal(clustlp);
0943
0944 double tmp_stub_x = posStub.x();
0945 double tmp_stub_y = posStub.y();
0946 double tmp_stub_z = posStub.z();
0947
0948 float trigDisplace = tempStubPtr->rawBend();
0949 float trigOffset = tempStubPtr->bendOffset();
0950 float trigPos = tempStubPtr->innerClusterPosition();
0951 float trigBend = tempStubPtr->bendFE();
0952
0953 m_allstub_x->push_back(tmp_stub_x);
0954 m_allstub_y->push_back(tmp_stub_y);
0955 m_allstub_z->push_back(tmp_stub_z);
0956
0957 m_allstub_isBarrel->push_back(isBarrel);
0958 m_allstub_layer->push_back(layer);
0959 m_allstub_isPSmodule->push_back(isPSmodule);
0960 m_allstub_isTiltedBarrel->push_back(isTiltedBarrel);
0961
0962 m_allstub_trigDisplace->push_back(trigDisplace);
0963 m_allstub_trigOffset->push_back(trigOffset);
0964 m_allstub_trigPos->push_back(trigPos);
0965 m_allstub_trigBend->push_back(trigBend);
0966
0967
0968 edm::Ptr<TrackingParticle> my_tp = MCTruthTTStubHandle->findTrackingParticlePtr(tempStubPtr);
0969
0970 int myTP_pdgid = -999;
0971 float myTP_pt = -999;
0972 float myTP_eta = -999;
0973 float myTP_phi = -999;
0974
0975 if (my_tp.isNull() == false) {
0976 int tmp_eventid = my_tp->eventId().event();
0977
0978 if (tmp_eventid > 0)
0979 continue;
0980
0981 myTP_pdgid = my_tp->pdgId();
0982 myTP_pt = my_tp->p4().pt();
0983 myTP_eta = my_tp->p4().eta();
0984 myTP_phi = my_tp->p4().phi();
0985 }
0986
0987 m_allstub_matchTP_pdgid->push_back(myTP_pdgid);
0988 m_allstub_matchTP_pt->push_back(myTP_pt);
0989 m_allstub_matchTP_eta->push_back(myTP_eta);
0990 m_allstub_matchTP_phi->push_back(myTP_phi);
0991
0992 int tmp_stub_genuine = 0;
0993 if (MCTruthTTStubHandle->isGenuine(tempStubPtr))
0994 tmp_stub_genuine = 1;
0995
0996 m_allstub_genuine->push_back(tmp_stub_genuine);
0997 }
0998 }
0999 }
1000
1001
1002
1003
1004
1005 std::vector<math::XYZTLorentzVector> v_jets;
1006 std::vector<int> v_jets_highpt;
1007 std::vector<int> v_jets_vhighpt;
1008
1009 if (TrackingInJets) {
1010
1011 if (DebugMode)
1012 edm::LogVerbatim("Tracklet") << "get genjets";
1013 edm::Handle<std::vector<reco::GenJet>> GenJetHandle;
1014 iEvent.getByToken(GenJetToken_, GenJetHandle);
1015
1016 if (GenJetHandle.isValid()) {
1017 if (DebugMode)
1018 edm::LogVerbatim("Tracklet") << "loop over genjets";
1019 std::vector<reco::GenJet>::const_iterator iterGenJet;
1020 for (iterGenJet = GenJetHandle->begin(); iterGenJet != GenJetHandle->end(); ++iterGenJet) {
1021 reco::GenJet myJet = reco::GenJet(*iterGenJet);
1022
1023 if (myJet.pt() < 30.0)
1024 continue;
1025 if (std::abs(myJet.eta()) > 2.5)
1026 continue;
1027
1028 if (DebugMode)
1029 edm::LogVerbatim("Tracklet") << "genjet pt = " << myJet.pt() << ", eta = " << myJet.eta();
1030
1031 bool ishighpt = false;
1032 bool isveryhighpt = false;
1033 if (myJet.pt() > 100.0)
1034 ishighpt = true;
1035 if (myJet.pt() > 200.0)
1036 isveryhighpt = true;
1037
1038 const math::XYZTLorentzVector& jetP4 = myJet.p4();
1039 v_jets.push_back(jetP4);
1040 if (ishighpt)
1041 v_jets_highpt.push_back(1);
1042 else
1043 v_jets_highpt.push_back(0);
1044 if (isveryhighpt)
1045 v_jets_vhighpt.push_back(1);
1046 else
1047 v_jets_vhighpt.push_back(0);
1048
1049 }
1050 }
1051
1052 }
1053
1054 const int NJETS = 10;
1055 float jets_tp_sumpt[NJETS] = {0};
1056 float jets_matchtrk_sumpt[NJETS] = {0};
1057 float jets_trk_sumpt[NJETS] = {0};
1058
1059
1060
1061
1062
1063 if (SaveAllTracks) {
1064 if (DebugMode) {
1065 edm::LogVerbatim("Tracklet") << "\n Loop over L1 tracks!";
1066 edm::LogVerbatim("Tracklet") << "\n Looking at " << L1Tk_nPar << "-parameter tracks!";
1067 }
1068
1069 int this_l1track = 0;
1070 std::vector<TTTrack<Ref_Phase2TrackerDigi_>>::const_iterator iterL1Track;
1071 for (iterL1Track = TTTrackHandle->begin(); iterL1Track != TTTrackHandle->end(); iterL1Track++) {
1072 edm::Ptr<TTTrack<Ref_Phase2TrackerDigi_>> l1track_ptr(TTTrackHandle, this_l1track);
1073 this_l1track++;
1074
1075 float tmp_trk_pt = iterL1Track->momentum().perp();
1076 float tmp_trk_eta = iterL1Track->momentum().eta();
1077 float tmp_trk_phi = iterL1Track->momentum().phi();
1078 float tmp_trk_z0 = iterL1Track->z0();
1079 float tmp_trk_tanL = iterL1Track->tanL();
1080 int tmp_trk_charge = (int)TMath::Sign(1, iterL1Track->rInv());
1081 bool usingNewKF = hphSetup->useNewKF();
1082 if (usingNewKF) {
1083
1084 constexpr float crazy_z0_cut = 30.;
1085 if (fabs(tmp_trk_z0) > crazy_z0_cut)
1086 continue;
1087 }
1088
1089 int tmp_trk_hitpattern = 0;
1090 tmp_trk_hitpattern = (int)iterL1Track->hitPattern();
1091 hph::HitPatternHelper hph(hphSetup, tmp_trk_hitpattern, tmp_trk_tanL, tmp_trk_z0);
1092 std::vector<int> hitpattern_expanded_binary = hph.binary();
1093 int tmp_trk_lhits_hitpattern = 0;
1094 int tmp_trk_dhits_hitpattern = 0;
1095 for (int i = 0; i < (int)hitpattern_expanded_binary.size(); i++) {
1096 if (hitpattern_expanded_binary[i]) {
1097 if (i < 6) {
1098 tmp_trk_lhits_hitpattern += pow(10, i);
1099 } else {
1100 tmp_trk_dhits_hitpattern += pow(10, i - 6);
1101 }
1102 }
1103 }
1104 int tmp_trk_nPSstub_hitpattern = hph.numPS();
1105 int tmp_trk_n2Sstub_hitpattern = hph.num2S();
1106 int tmp_trk_nLostPSstub_hitpattern = hph.numMissingPS();
1107 int tmp_trk_nLost2Sstub_hitpattern = hph.numMissing2S();
1108 int tmp_trk_nLoststub_V1_hitpattern = hph.numMissingInterior1();
1109 int tmp_trk_nLoststub_V2_hitpattern = hph.numMissingInterior2();
1110
1111 float tmp_trk_d0 = -999;
1112 if (L1Tk_nPar == 5) {
1113 float tmp_trk_x0 = iterL1Track->POCA().x();
1114 float tmp_trk_y0 = iterL1Track->POCA().y();
1115 tmp_trk_d0 = tmp_trk_x0 * sin(tmp_trk_phi) - tmp_trk_y0 * cos(tmp_trk_phi);
1116 }
1117
1118 float tmp_trk_chi2 = iterL1Track->chi2();
1119 float tmp_trk_chi2rphi = iterL1Track->chi2XY();
1120 float tmp_trk_chi2rz = iterL1Track->chi2Z();
1121 float tmp_trk_bendchi2 = iterL1Track->stubPtConsistency();
1122 float tmp_trk_MVA1 = iterL1Track->trkMVA1();
1123
1124 std::vector<edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>>>
1125 stubRefs = iterL1Track->getStubRefs();
1126 int tmp_trk_nstub = (int)stubRefs.size();
1127 int ndof = 2 * tmp_trk_nstub - L1Tk_nPar;
1128 int ndofrphi = tmp_trk_nstub - L1Tk_nPar + 2;
1129 int ndofrz = tmp_trk_nstub - 2;
1130 float tmp_trk_chi2_dof = (float)tmp_trk_chi2 / ndof;
1131 float tmp_trk_chi2rphi_dof = (float)tmp_trk_chi2rphi / ndofrphi;
1132 float tmp_trk_chi2rz_dof = (float)tmp_trk_chi2rz / ndofrz;
1133
1134 int tmp_trk_seed = 0;
1135 tmp_trk_seed = (int)iterL1Track->trackSeedType();
1136
1137 unsigned int tmp_trk_phiSector = iterL1Track->phiSector();
1138 int tmp_trk_etaSector = hph.etaSector();
1139
1140
1141
1142
1143
1144 int tmp_trk_dhits = 0;
1145 int tmp_trk_lhits = 0;
1146
1147 if (true) {
1148
1149 for (int is = 0; is < tmp_trk_nstub; is++) {
1150
1151 DetId detIdStub = theTrackerGeom->idToDet((stubRefs.at(is)->clusterRef(0))->getDetId())->geographicalId();
1152
1153 MeasurementPoint coords = stubRefs.at(is)->clusterRef(0)->findAverageLocalCoordinatesCentered();
1154 const GeomDet* theGeomDet = theTrackerGeom->idToDet(detIdStub);
1155 Global3DPoint posStub = theGeomDet->surface().toGlobal(theGeomDet->topology().localPosition(coords));
1156
1157 double x = posStub.x();
1158 double y = posStub.y();
1159 double z = posStub.z();
1160
1161 int layer = -999999;
1162 if (detIdStub.subdetId() == StripSubdetector::TOB) {
1163 layer = static_cast<int>(tTopo->layer(detIdStub));
1164 if (DebugMode)
1165 edm::LogVerbatim("Tracklet")
1166 << " stub in layer " << layer << " at position x y z = " << x << " " << y << " " << z;
1167 tmp_trk_lhits += pow(10, layer - 1);
1168 } else if (detIdStub.subdetId() == StripSubdetector::TID) {
1169 layer = static_cast<int>(tTopo->layer(detIdStub));
1170 if (DebugMode)
1171 edm::LogVerbatim("Tracklet")
1172 << " stub in disk " << layer << " at position x y z = " << x << " " << y << " " << z;
1173 tmp_trk_dhits += pow(10, layer - 1);
1174 }
1175
1176 }
1177 }
1178
1179
1180 int tmp_trk_genuine = 0;
1181 int tmp_trk_loose = 0;
1182 int tmp_trk_unknown = 0;
1183 int tmp_trk_combinatoric = 0;
1184 if (MCTruthTTTrackHandle->isLooselyGenuine(l1track_ptr))
1185 tmp_trk_loose = 1;
1186 if (MCTruthTTTrackHandle->isGenuine(l1track_ptr))
1187 tmp_trk_genuine = 1;
1188 if (MCTruthTTTrackHandle->isUnknown(l1track_ptr))
1189 tmp_trk_unknown = 1;
1190 if (MCTruthTTTrackHandle->isCombinatoric(l1track_ptr))
1191 tmp_trk_combinatoric = 1;
1192
1193 if (DebugMode) {
1194 edm::LogVerbatim("Tracklet") << "L1 track,"
1195 << " pt: " << tmp_trk_pt << " eta: " << tmp_trk_eta << " phi: " << tmp_trk_phi
1196 << " z0: " << tmp_trk_z0 << " chi2: " << tmp_trk_chi2
1197 << " chi2rphi: " << tmp_trk_chi2rphi << " chi2rz: " << tmp_trk_chi2rz
1198 << " nstub: " << tmp_trk_nstub;
1199 if (tmp_trk_genuine)
1200 edm::LogVerbatim("Tracklet") << " (is genuine)";
1201 if (tmp_trk_unknown)
1202 edm::LogVerbatim("Tracklet") << " (is unknown)";
1203 if (tmp_trk_combinatoric)
1204 edm::LogVerbatim("Tracklet") << " (is combinatoric)";
1205 }
1206
1207 m_trk_pt->push_back(tmp_trk_pt);
1208 m_trk_eta->push_back(tmp_trk_eta);
1209 m_trk_phi->push_back(tmp_trk_phi);
1210 m_trk_z0->push_back(tmp_trk_z0);
1211 if (L1Tk_nPar == 5)
1212 m_trk_d0->push_back(tmp_trk_d0);
1213 else
1214 m_trk_d0->push_back(999.);
1215 m_trk_chi2->push_back(tmp_trk_chi2);
1216 m_trk_chi2_dof->push_back(tmp_trk_chi2_dof);
1217 m_trk_chi2rphi->push_back(tmp_trk_chi2rphi);
1218 m_trk_chi2rphi_dof->push_back(tmp_trk_chi2rphi_dof);
1219 m_trk_chi2rz->push_back(tmp_trk_chi2rz);
1220 m_trk_chi2rz_dof->push_back(tmp_trk_chi2rz_dof);
1221 m_trk_bendchi2->push_back(tmp_trk_bendchi2);
1222 m_trk_MVA1->push_back(tmp_trk_MVA1);
1223 m_trk_nstub->push_back(tmp_trk_nstub);
1224 m_trk_dhits->push_back(tmp_trk_dhits);
1225 m_trk_lhits->push_back(tmp_trk_lhits);
1226 m_trk_seed->push_back(tmp_trk_seed);
1227 m_trk_hitpattern->push_back(tmp_trk_hitpattern);
1228 m_trk_lhits_hitpattern->push_back(tmp_trk_lhits_hitpattern);
1229 m_trk_dhits_hitpattern->push_back(tmp_trk_dhits_hitpattern);
1230 m_trk_nPSstub_hitpattern->push_back(tmp_trk_nPSstub_hitpattern);
1231 m_trk_n2Sstub_hitpattern->push_back(tmp_trk_n2Sstub_hitpattern);
1232 m_trk_nLostPSstub_hitpattern->push_back(tmp_trk_nLostPSstub_hitpattern);
1233 m_trk_nLost2Sstub_hitpattern->push_back(tmp_trk_nLost2Sstub_hitpattern);
1234 m_trk_nLoststub_V1_hitpattern->push_back(tmp_trk_nLoststub_V1_hitpattern);
1235 m_trk_nLoststub_V2_hitpattern->push_back(tmp_trk_nLoststub_V2_hitpattern);
1236 m_trk_charge->push_back(tmp_trk_charge);
1237 m_trk_phiSector->push_back(tmp_trk_phiSector);
1238 m_trk_etaSector->push_back(tmp_trk_etaSector);
1239 m_trk_genuine->push_back(tmp_trk_genuine);
1240 m_trk_loose->push_back(tmp_trk_loose);
1241 m_trk_unknown->push_back(tmp_trk_unknown);
1242 m_trk_combinatoric->push_back(tmp_trk_combinatoric);
1243
1244
1245
1246
1247
1248 edm::Ptr<TrackingParticle> my_tp = MCTruthTTTrackHandle->findTrackingParticlePtr(l1track_ptr);
1249
1250 int myFake = 0;
1251
1252 int tmp_matchtp_pdgid = -999;
1253 float tmp_matchtp_pt = -999;
1254 float tmp_matchtp_eta = -999;
1255 float tmp_matchtp_phi = -999;
1256 float tmp_matchtp_z0 = -999;
1257 float tmp_matchtp_lxy = -999;
1258 float tmp_matchtp_d0 = -999;
1259
1260 if (my_tp.isNull())
1261 myFake = 0;
1262 else {
1263 int tmp_eventid = my_tp->eventId().event();
1264
1265 if (tmp_eventid > 0)
1266 myFake = 2;
1267 else
1268 myFake = 1;
1269
1270 tmp_matchtp_pdgid = my_tp->pdgId();
1271 tmp_matchtp_pt = my_tp->pt();
1272 tmp_matchtp_eta = my_tp->eta();
1273 tmp_matchtp_phi = my_tp->phi();
1274
1275 float tmp_matchtp_vz = my_tp->vz();
1276 float tmp_matchtp_vx = my_tp->vx();
1277 float tmp_matchtp_vy = my_tp->vy();
1278 tmp_matchtp_lxy = sqrt(tmp_matchtp_vx * tmp_matchtp_vx + tmp_matchtp_vy * tmp_matchtp_vy);
1279
1280
1281
1282
1283 float tmp_matchtp_t = 1.0 / tan(2.0 * atan(exp(-tmp_matchtp_eta)));
1284
1285 float delx = -tmp_matchtp_vx;
1286 float dely = -tmp_matchtp_vy;
1287
1288 float b_field = bFieldHandle.product()->inTesla(GlobalPoint(0, 0, 0)).z();
1289 float c_converted = CLHEP::c_light / 1.0E5;
1290 float r2_inv = my_tp->charge() * c_converted * b_field / tmp_matchtp_pt / 2.0;
1291
1292 float tmp_matchtp_x0p = delx - (1. / (2. * r2_inv) * sin(tmp_matchtp_phi));
1293 float tmp_matchtp_y0p = dely + (1. / (2. * r2_inv) * cos(tmp_matchtp_phi));
1294 float tmp_matchtp_rp = sqrt(tmp_matchtp_x0p * tmp_matchtp_x0p + tmp_matchtp_y0p * tmp_matchtp_y0p);
1295 tmp_matchtp_d0 = my_tp->charge() * tmp_matchtp_rp - (1. / (2. * r2_inv));
1296
1297 static double pi = M_PI;
1298 float delphi = tmp_matchtp_phi - atan2(-r2_inv * tmp_matchtp_x0p, r2_inv * tmp_matchtp_y0p);
1299 if (delphi < -pi)
1300 delphi += 2.0 * pi;
1301 if (delphi > pi)
1302 delphi -= 2.0 * pi;
1303 tmp_matchtp_z0 = tmp_matchtp_vz + tmp_matchtp_t * delphi / (2.0 * r2_inv);
1304
1305
1306 if (DebugMode) {
1307 edm::LogVerbatim("Tracklet") << "TP matched to track has pt = " << my_tp->p4().pt()
1308 << " eta = " << my_tp->momentum().eta() << " phi = " << my_tp->momentum().phi()
1309 << " z0 = " << my_tp->vertex().z() << " pdgid = " << my_tp->pdgId()
1310 << " lxy = " << tmp_matchtp_lxy;
1311 }
1312 }
1313
1314 m_trk_fake->push_back(myFake);
1315
1316 m_trk_matchtp_pdgid->push_back(tmp_matchtp_pdgid);
1317 m_trk_matchtp_pt->push_back(tmp_matchtp_pt);
1318 m_trk_matchtp_eta->push_back(tmp_matchtp_eta);
1319 m_trk_matchtp_phi->push_back(tmp_matchtp_phi);
1320 m_trk_matchtp_z0->push_back(tmp_matchtp_z0);
1321 m_trk_matchtp_lxy->push_back(tmp_matchtp_lxy);
1322 m_trk_matchtp_d0->push_back(tmp_matchtp_d0);
1323
1324
1325
1326
1327
1328 if (TrackingInJets) {
1329 if (DebugMode)
1330 edm::LogVerbatim("Tracklet") << "doing tracking in jets now";
1331
1332 int InJet = 0;
1333 int InJetHighpt = 0;
1334 int InJetVeryHighpt = 0;
1335
1336 for (int ij = 0; ij < (int)v_jets.size(); ij++) {
1337 float deta = tmp_trk_eta - (v_jets.at(ij)).eta();
1338 float dphi = tmp_trk_phi - (v_jets.at(ij)).phi();
1339 while (dphi > 3.14159)
1340 dphi = std::abs(2 * 3.14159 - dphi);
1341 float dR = sqrt(deta * deta + dphi * dphi);
1342
1343 if (dR < 0.4) {
1344 InJet = 1;
1345 if (v_jets_highpt.at(ij) == 1)
1346 InJetHighpt = 1;
1347 if (v_jets_vhighpt.at(ij) == 1)
1348 InJetVeryHighpt = 1;
1349 if (ij < NJETS)
1350 jets_trk_sumpt[ij] += tmp_trk_pt;
1351 }
1352 }
1353
1354 m_trk_injet->push_back(InJet);
1355 m_trk_injet_highpt->push_back(InJetHighpt);
1356 m_trk_injet_vhighpt->push_back(InJetVeryHighpt);
1357
1358 }
1359
1360
1361 const TTBV hitPattern((int)iterL1Track->hitPattern(), setup->numLayers());
1362 const double zT = iterL1Track->z0() + setup->chosenRofZ() * iterL1Track->tanL();
1363 const vector<int>& le = layerEncoding->layerEncoding(zT);
1364 vector<int> layers;
1365 layers.reserve(hitPattern.size());
1366 for (int layer : hitPattern.ids())
1367 layers.push_back(le[layer]);
1368 m_trk_layers->push_back(layers);
1369
1370 }
1371
1372 }
1373
1374
1375
1376
1377
1378 if (DebugMode)
1379 edm::LogVerbatim("Tracklet") << "\n Loop over tracking particles!";
1380
1381 int this_tp = 0;
1382 std::vector<TrackingParticle>::const_iterator iterTP;
1383 for (iterTP = TrackingParticleHandle->begin(); iterTP != TrackingParticleHandle->end(); ++iterTP) {
1384 edm::Ptr<TrackingParticle> tp_ptr(TrackingParticleHandle, this_tp);
1385 this_tp++;
1386
1387 int tmp_eventid = iterTP->eventId().event();
1388 if (MyProcess != 1 && tmp_eventid > 0)
1389 continue;
1390
1391 float tmp_tp_pt = iterTP->pt();
1392 float tmp_tp_charge = iterTP->charge();
1393 float tmp_tp_eta = iterTP->eta();
1394
1395 if (tmp_tp_pt < TP_minPt)
1396 continue;
1397 if (tmp_tp_charge == 0.)
1398 continue;
1399 if (std::abs(tmp_tp_eta) > TP_maxEta)
1400 continue;
1401
1402 float tmp_tp_phi = iterTP->phi();
1403 float tmp_tp_vz = iterTP->vz();
1404 float tmp_tp_vx = iterTP->vx();
1405 float tmp_tp_vy = iterTP->vy();
1406 int tmp_tp_pdgid = iterTP->pdgId();
1407 float tmp_tp_z0_prod = tmp_tp_vz;
1408 float tmp_tp_d0_prod = tmp_tp_vx * sin(tmp_tp_phi) - tmp_tp_vy * cos(tmp_tp_phi);
1409
1410
1411
1412
1413 float tmp_tp_t = 1.0 / tan(2.0 * atan(exp(-tmp_tp_eta)));
1414
1415 float delx = -tmp_tp_vx;
1416 float dely = -tmp_tp_vy;
1417
1418 float b_field = bFieldHandle.product()->inTesla(GlobalPoint(0, 0, 0)).z();
1419 float c_converted = CLHEP::c_light / 1.0E5;
1420 float r2_inv = tmp_tp_charge * c_converted * b_field / tmp_tp_pt / 2.0;
1421
1422 float tmp_tp_x0p = delx - (1. / (2. * r2_inv) * sin(tmp_tp_phi));
1423 float tmp_tp_y0p = dely + (1. / (2. * r2_inv) * cos(tmp_tp_phi));
1424 float tmp_tp_rp = sqrt(tmp_tp_x0p * tmp_tp_x0p + tmp_tp_y0p * tmp_tp_y0p);
1425 float tmp_tp_d0 = tmp_tp_charge * tmp_tp_rp - (1. / (2. * r2_inv));
1426
1427 static double pi = M_PI;
1428 float delphi = tmp_tp_phi - atan2(-r2_inv * tmp_tp_x0p, r2_inv * tmp_tp_y0p);
1429 if (delphi < -pi)
1430 delphi += 2.0 * pi;
1431 if (delphi > pi)
1432 delphi -= 2.0 * pi;
1433 float tmp_tp_z0 = tmp_tp_vz + tmp_tp_t * delphi / (2.0 * r2_inv);
1434
1435
1436 if (MyProcess == 13 && abs(tmp_tp_pdgid) != 13)
1437 continue;
1438 if (MyProcess == 11 && abs(tmp_tp_pdgid) != 11)
1439 continue;
1440 if ((MyProcess == 6 || MyProcess == 15 || MyProcess == 211) && abs(tmp_tp_pdgid) != 211)
1441 continue;
1442
1443 if (std::abs(tmp_tp_z0) > TP_maxZ0)
1444 continue;
1445
1446
1447 float lxy = sqrt(tmp_tp_vx * tmp_tp_vx + tmp_tp_vy * tmp_tp_vy);
1448 float tmp_tp_lxy = lxy;
1449 if (MyProcess == 6 && (lxy > 1.0))
1450 continue;
1451
1452 if (DebugMode)
1453 edm::LogVerbatim("Tracklet") << "Tracking particle, pt: " << tmp_tp_pt << " eta: " << tmp_tp_eta
1454 << " phi: " << tmp_tp_phi << " z0: " << tmp_tp_z0 << " d0: " << tmp_tp_d0
1455 << " z_prod: " << tmp_tp_z0_prod << " d_prod: " << tmp_tp_d0_prod
1456 << " pdgid: " << tmp_tp_pdgid << " eventID: " << iterTP->eventId().event()
1457 << " ttclusters " << MCTruthTTClusterHandle->findTTClusterRefs(tp_ptr).size()
1458 << " ttstubs " << MCTruthTTStubHandle->findTTStubRefs(tp_ptr).size() << " tttracks "
1459 << MCTruthTTTrackHandle->findTTTrackPtrs(tp_ptr).size();
1460
1461
1462
1463
1464 if (MCTruthTTClusterHandle->findTTClusterRefs(tp_ptr).empty()) {
1465 if (DebugMode)
1466 edm::LogVerbatim("Tracklet") << "No matching TTClusters for TP, continuing...";
1467 continue;
1468 }
1469
1470 std::vector<edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>>>
1471 theStubRefs = MCTruthTTStubHandle->findTTStubRefs(tp_ptr);
1472 int nStubTP = (int)theStubRefs.size();
1473
1474
1475 int hasStubInLayer[11] = {0};
1476 for (auto& theStubRef : theStubRefs) {
1477 DetId detid(theStubRef->getDetId());
1478
1479 int layer = -1;
1480 if (detid.subdetId() == StripSubdetector::TOB) {
1481 layer = static_cast<int>(tTopo->layer(detid)) - 1;
1482 } else if (detid.subdetId() == StripSubdetector::TID) {
1483 layer = static_cast<int>(tTopo->layer(detid)) + 5;
1484 }
1485
1486
1487
1488
1489 if (MCTruthTTStubHandle->findTrackingParticlePtr(theStubRef).isNull() && hasStubInLayer[layer] < 2)
1490 hasStubInLayer[layer] = 1;
1491 else
1492 hasStubInLayer[layer] = 2;
1493 }
1494
1495 int nStubLayerTP = 0;
1496 int nStubLayerTP_g = 0;
1497 for (int isum : hasStubInLayer) {
1498 if (isum >= 1)
1499 nStubLayerTP += 1;
1500 if (isum == 2)
1501 nStubLayerTP_g += 1;
1502 }
1503
1504 if (DebugMode)
1505 edm::LogVerbatim("Tracklet") << "TP is associated with " << nStubTP << " stubs, and has stubs in " << nStubLayerTP
1506 << " different layers/disks, and has GENUINE stubs in " << nStubLayerTP_g
1507 << " layers ";
1508
1509 if (TP_minNStub > 0) {
1510 if (DebugMode)
1511 edm::LogVerbatim("Tracklet") << "Only consider TPs with >= " << TP_minNStub << " stubs";
1512 if (nStubTP < TP_minNStub) {
1513 if (DebugMode)
1514 edm::LogVerbatim("Tracklet") << "TP fails minimum nbr stubs requirement! Continuing...";
1515 continue;
1516 }
1517 }
1518 if (TP_minNStubLayer > 0) {
1519 if (DebugMode)
1520 edm::LogVerbatim("Tracklet") << "Only consider TPs with stubs in >= " << TP_minNStubLayer << " layers/disks";
1521 if (nStubLayerTP < TP_minNStubLayer) {
1522 if (DebugMode)
1523 edm::LogVerbatim("Tracklet") << "TP fails stubs in minimum nbr of layers/disks requirement! Continuing...";
1524 continue;
1525 }
1526 }
1527
1528
1529
1530
1531 std::vector<edm::Ptr<TTTrack<Ref_Phase2TrackerDigi_>>> matchedTracks =
1532 MCTruthTTTrackHandle->findTTTrackPtrs(tp_ptr);
1533
1534 int nMatch = 0;
1535 int i_track = -1;
1536 float i_chi2dof = 99999;
1537
1538 if (!matchedTracks.empty()) {
1539 if (DebugMode && (matchedTracks.size() > 1))
1540 edm::LogVerbatim("Tracklet") << "TrackingParticle has more than one matched L1 track!";
1541
1542
1543
1544
1545
1546
1547 for (int it = 0; it < (int)matchedTracks.size(); it++) {
1548 bool tmp_trk_genuine = false;
1549 bool tmp_trk_loosegenuine = false;
1550 if (MCTruthTTTrackHandle->isGenuine(matchedTracks.at(it)))
1551 tmp_trk_genuine = true;
1552 if (MCTruthTTTrackHandle->isLooselyGenuine(matchedTracks.at(it)))
1553 tmp_trk_loosegenuine = true;
1554 if (!tmp_trk_loosegenuine)
1555 continue;
1556
1557 if (DebugMode) {
1558 if (MCTruthTTTrackHandle->findTrackingParticlePtr(matchedTracks.at(it)).isNull()) {
1559 edm::LogVerbatim("Tracklet") << "track matched to TP is NOT uniquely matched to a TP";
1560 } else {
1561 edm::Ptr<TrackingParticle> my_tp = MCTruthTTTrackHandle->findTrackingParticlePtr(matchedTracks.at(it));
1562 edm::LogVerbatim("Tracklet") << "TP matched to track matched to TP ... tp pt = " << my_tp->p4().pt()
1563 << " eta = " << my_tp->momentum().eta() << " phi = " << my_tp->momentum().phi()
1564 << " z0 = " << my_tp->vertex().z();
1565 }
1566 edm::LogVerbatim("Tracklet") << " ... matched L1 track has pt = " << matchedTracks.at(it)->momentum().perp()
1567 << " eta = " << matchedTracks.at(it)->momentum().eta()
1568 << " phi = " << matchedTracks.at(it)->momentum().phi()
1569 << " chi2 = " << matchedTracks.at(it)->chi2()
1570 << " consistency = " << matchedTracks.at(it)->stubPtConsistency()
1571 << " z0 = " << matchedTracks.at(it)->z0()
1572 << " nstub = " << matchedTracks.at(it)->getStubRefs().size();
1573 if (tmp_trk_genuine)
1574 edm::LogVerbatim("Tracklet") << " (genuine!) ";
1575 if (tmp_trk_loosegenuine)
1576 edm::LogVerbatim("Tracklet") << " (loose genuine!) ";
1577 }
1578
1579
1580
1581
1582
1583
1584 std::vector<edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>>>
1585 stubRefs = matchedTracks.at(it)->getStubRefs();
1586 int tmp_trk_nstub = stubRefs.size();
1587
1588 if (tmp_trk_nstub < L1Tk_minNStub)
1589 continue;
1590
1591
1592
1593
1594
1595
1596
1597
1598
1599
1600
1601
1602 float dmatch_pt = 999;
1603 float dmatch_eta = 999;
1604 float dmatch_phi = 999;
1605 int match_id = 999;
1606
1607 edm::Ptr<TrackingParticle> my_tp = MCTruthTTTrackHandle->findTrackingParticlePtr(matchedTracks.at(it));
1608 dmatch_pt = std::abs(my_tp->p4().pt() - tmp_tp_pt);
1609 dmatch_eta = std::abs(my_tp->p4().eta() - tmp_tp_eta);
1610 dmatch_phi = std::abs(my_tp->p4().phi() - tmp_tp_phi);
1611 match_id = my_tp->pdgId();
1612
1613 float tmp_trk_chi2dof = (matchedTracks.at(it)->chi2()) / (2 * tmp_trk_nstub - L1Tk_nPar);
1614
1615
1616 if (dmatch_pt < 0.1 && dmatch_eta < 0.1 && dmatch_phi < 0.1 && tmp_tp_pdgid == match_id && tmp_trk_genuine) {
1617 nMatch++;
1618 if (i_track < 0 || tmp_trk_chi2dof < i_chi2dof) {
1619 i_track = it;
1620 i_chi2dof = tmp_trk_chi2dof;
1621 }
1622 }
1623
1624 }
1625
1626 }
1627
1628
1629 float tmp_matchtrk_pt = -999;
1630 float tmp_matchtrk_eta = -999;
1631 float tmp_matchtrk_phi = -999;
1632 float tmp_matchtrk_z0 = -999;
1633 float tmp_matchtrk_d0 = -999;
1634 float tmp_matchtrk_chi2 = -999;
1635 float tmp_matchtrk_chi2_dof = -999;
1636 float tmp_matchtrk_chi2rphi = -999;
1637 float tmp_matchtrk_chi2rphi_dof = -999;
1638 float tmp_matchtrk_chi2rz = -999;
1639 float tmp_matchtrk_chi2rz_dof = -999;
1640 float tmp_matchtrk_bendchi2 = -999;
1641 float tmp_matchtrk_MVA1 = -999;
1642 int tmp_matchtrk_charge = -999;
1643 int tmp_matchtrk_nstub = -999;
1644 int tmp_matchtrk_dhits = -999;
1645 int tmp_matchtrk_lhits = -999;
1646 int tmp_matchtrk_seed = -999;
1647 int tmp_matchtrk_hitpattern = -999;
1648
1649 if (nMatch > 1 && DebugMode)
1650 edm::LogVerbatim("Tracklet") << "WARNING *** 2 or more matches to genuine L1 tracks ***";
1651
1652 if (nMatch > 0) {
1653 tmp_matchtrk_pt = matchedTracks.at(i_track)->momentum().perp();
1654 tmp_matchtrk_charge = (int)TMath::Sign(1, matchedTracks.at(i_track)->rInv());
1655 tmp_matchtrk_eta = matchedTracks.at(i_track)->momentum().eta();
1656 tmp_matchtrk_phi = matchedTracks.at(i_track)->momentum().phi();
1657 tmp_matchtrk_z0 = matchedTracks.at(i_track)->z0();
1658
1659 if (L1Tk_nPar == 5) {
1660 float tmp_matchtrk_x0 = matchedTracks.at(i_track)->POCA().x();
1661 float tmp_matchtrk_y0 = matchedTracks.at(i_track)->POCA().y();
1662 tmp_matchtrk_d0 = tmp_matchtrk_x0 * sin(tmp_matchtrk_phi) - tmp_matchtrk_y0 * cos(tmp_matchtrk_phi);
1663 }
1664
1665 tmp_matchtrk_chi2 = matchedTracks.at(i_track)->chi2();
1666 tmp_matchtrk_chi2rphi = matchedTracks.at(i_track)->chi2XY();
1667 tmp_matchtrk_chi2rz = matchedTracks.at(i_track)->chi2Z();
1668 tmp_matchtrk_bendchi2 = matchedTracks.at(i_track)->stubPtConsistency();
1669 tmp_matchtrk_MVA1 = matchedTracks.at(i_track)->trkMVA1();
1670 tmp_matchtrk_nstub = (int)matchedTracks.at(i_track)->getStubRefs().size();
1671 tmp_matchtrk_seed = (int)matchedTracks.at(i_track)->trackSeedType();
1672 tmp_matchtrk_hitpattern = (int)matchedTracks.at(i_track)->hitPattern();
1673
1674 int ndof = 2 * tmp_matchtrk_nstub - L1Tk_nPar;
1675 int ndofrphi = tmp_matchtrk_nstub - L1Tk_nPar + 2;
1676 int ndofrz = tmp_matchtrk_nstub - 2;
1677 tmp_matchtrk_chi2_dof = (float)tmp_matchtrk_chi2 / ndof;
1678 tmp_matchtrk_chi2rphi_dof = (float)tmp_matchtrk_chi2rphi / ndofrphi;
1679 tmp_matchtrk_chi2rz_dof = (float)tmp_matchtrk_chi2rz / ndofrz;
1680
1681
1682
1683
1684
1685 tmp_matchtrk_dhits = 0;
1686 tmp_matchtrk_lhits = 0;
1687
1688 std::vector<edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>>>
1689 stubRefs = matchedTracks.at(i_track)->getStubRefs();
1690 int tmp_nstub = stubRefs.size();
1691
1692 for (int is = 0; is < tmp_nstub; is++) {
1693 DetId detIdStub = theTrackerGeom->idToDet((stubRefs.at(is)->clusterRef(0))->getDetId())->geographicalId();
1694
1695
1696
1697
1698
1699
1700 int layer = -999999;
1701 if (detIdStub.subdetId() == StripSubdetector::TOB) {
1702 layer = static_cast<int>(tTopo->layer(detIdStub));
1703 tmp_matchtrk_lhits += pow(10, layer - 1);
1704 } else if (detIdStub.subdetId() == StripSubdetector::TID) {
1705 layer = static_cast<int>(tTopo->layer(detIdStub));
1706 tmp_matchtrk_dhits += pow(10, layer - 1);
1707 }
1708
1709
1710 }
1711 }
1712
1713 m_tp_pt->push_back(tmp_tp_pt);
1714 m_tp_eta->push_back(tmp_tp_eta);
1715 m_tp_phi->push_back(tmp_tp_phi);
1716 m_tp_lxy->push_back(tmp_tp_lxy);
1717 m_tp_z0->push_back(tmp_tp_z0);
1718 m_tp_d0->push_back(tmp_tp_d0);
1719 m_tp_z0_prod->push_back(tmp_tp_z0_prod);
1720 m_tp_d0_prod->push_back(tmp_tp_d0_prod);
1721 m_tp_pdgid->push_back(tmp_tp_pdgid);
1722 m_tp_nmatch->push_back(nMatch);
1723 m_tp_nstub->push_back(nStubTP);
1724 m_tp_eventid->push_back(tmp_eventid);
1725 m_tp_charge->push_back(tmp_tp_charge);
1726
1727 m_matchtrk_pt->push_back(tmp_matchtrk_pt);
1728 m_matchtrk_eta->push_back(tmp_matchtrk_eta);
1729 m_matchtrk_phi->push_back(tmp_matchtrk_phi);
1730 m_matchtrk_z0->push_back(tmp_matchtrk_z0);
1731 m_matchtrk_d0->push_back(tmp_matchtrk_d0);
1732 m_matchtrk_chi2->push_back(tmp_matchtrk_chi2);
1733 m_matchtrk_chi2rphi->push_back(tmp_matchtrk_chi2rphi);
1734 m_matchtrk_chi2rz->push_back(tmp_matchtrk_chi2rz);
1735 m_matchtrk_bendchi2->push_back(tmp_matchtrk_bendchi2);
1736 m_matchtrk_MVA1->push_back(tmp_matchtrk_MVA1);
1737 m_matchtrk_nstub->push_back(tmp_matchtrk_nstub);
1738 m_matchtrk_dhits->push_back(tmp_matchtrk_dhits);
1739 m_matchtrk_lhits->push_back(tmp_matchtrk_lhits);
1740 m_matchtrk_seed->push_back(tmp_matchtrk_seed);
1741 m_matchtrk_hitpattern->push_back(tmp_matchtrk_hitpattern);
1742 m_matchtrk_charge->push_back(tmp_matchtrk_charge);
1743 m_matchtrk_chi2_dof->push_back(tmp_matchtrk_chi2_dof);
1744 m_matchtrk_chi2rphi_dof->push_back(tmp_matchtrk_chi2rphi_dof);
1745 m_matchtrk_chi2rz_dof->push_back(tmp_matchtrk_chi2rz_dof);
1746
1747
1748
1749
1750
1751 if (TrackingInJets) {
1752 if (DebugMode)
1753 edm::LogVerbatim("Tracklet") << "check if TP/matched track is within jet";
1754
1755 int tp_InJet = 0;
1756 int matchtrk_InJet = 0;
1757 int tp_InJetHighpt = 0;
1758 int matchtrk_InJetHighpt = 0;
1759 int tp_InJetVeryHighpt = 0;
1760 int matchtrk_InJetVeryHighpt = 0;
1761
1762 for (int ij = 0; ij < (int)v_jets.size(); ij++) {
1763 float deta = tmp_tp_eta - (v_jets.at(ij)).eta();
1764 float dphi = tmp_tp_phi - (v_jets.at(ij)).phi();
1765 while (dphi > 3.14159)
1766 dphi = std::abs(2 * 3.14159 - dphi);
1767 float dR = sqrt(deta * deta + dphi * dphi);
1768 if (dR < 0.4) {
1769 tp_InJet = 1;
1770 if (v_jets_highpt.at(ij) == 1)
1771 tp_InJetHighpt = 1;
1772 if (v_jets_vhighpt.at(ij) == 1)
1773 tp_InJetVeryHighpt = 1;
1774 if (ij < NJETS)
1775 jets_tp_sumpt[ij] += tmp_tp_pt;
1776 }
1777
1778 if (nMatch > 0) {
1779 deta = tmp_matchtrk_eta - (v_jets.at(ij)).eta();
1780 dphi = tmp_matchtrk_phi - (v_jets.at(ij)).phi();
1781 while (dphi > 3.14159)
1782 dphi = std::abs(2 * 3.14159 - dphi);
1783 dR = sqrt(deta * deta + dphi * dphi);
1784 if (dR < 0.4) {
1785 matchtrk_InJet = 1;
1786 if (v_jets_highpt.at(ij) == 1)
1787 matchtrk_InJetHighpt = 1;
1788 if (v_jets_vhighpt.at(ij) == 1)
1789 matchtrk_InJetVeryHighpt = 1;
1790 if (ij < NJETS)
1791 jets_matchtrk_sumpt[ij] += tmp_matchtrk_pt;
1792 }
1793 }
1794 }
1795
1796 m_tp_injet->push_back(tp_InJet);
1797 m_tp_injet_highpt->push_back(tp_InJetHighpt);
1798 m_tp_injet_vhighpt->push_back(tp_InJetVeryHighpt);
1799 m_matchtrk_injet->push_back(matchtrk_InJet);
1800 m_matchtrk_injet_highpt->push_back(matchtrk_InJetHighpt);
1801 m_matchtrk_injet_vhighpt->push_back(matchtrk_InJetVeryHighpt);
1802
1803 }
1804
1805 }
1806
1807 if (TrackingInJets) {
1808 for (int ij = 0; ij < (int)v_jets.size(); ij++) {
1809 if (ij < NJETS) {
1810 m_jet_eta->push_back((v_jets.at(ij)).eta());
1811 m_jet_phi->push_back((v_jets.at(ij)).phi());
1812 m_jet_pt->push_back((v_jets.at(ij)).pt());
1813 m_jet_tp_sumpt->push_back(jets_tp_sumpt[ij]);
1814 m_jet_trk_sumpt->push_back(jets_trk_sumpt[ij]);
1815 m_jet_matchtrk_sumpt->push_back(jets_matchtrk_sumpt[ij]);
1816 }
1817 }
1818 }
1819
1820 eventTree->Fill();
1821
1822 }
1823
1824
1825
1826 DEFINE_FWK_MODULE(L1TrackNtupleMaker);