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