Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:22:07

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