Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-10-14 01:44:19

0001 // This is a helper function that can be used to decode hitpattern, which is a 7-bit integer produced by the Kalman filter (KF).
0002 // Hitpattern is stored at TTTrack objects (DataFormats/L1TrackTrigger/interface/TTTrack.h)
0003 // It can be accessed via TTTrack: hitPattern()
0004 //
0005 // There are two classes declared in HitPatternHelper (hph) namesapce:
0006 // 1)Setup: This is used to produce a layermap and a collection of <tt::SensorModule> needed by HitPatternHelper.
0007 // 2)HitPatternHelper: This function returns more specific information (e.g. module type, layer id,...) about each stub on the TTTrack objects.
0008 // This function needs three variables from TTTrack: hitPattern(),tanL() and z0().
0009 // It makes predictions in two different ways depending on which version of the KF is deployed:
0010 //
0011 // Old KF (L1Trigger/TrackFindingTMTT/plugins/TMTrackProducer.cc) is a CMSSW emulation of KF firmware that
0012 // ignores truncation effect and is not bit/clock cycle accurate.
0013 // With this version of the KF, HitPatternHelper relys on a hard-coded layermap to deterimne layer IDs of each stub.
0014 //
0015 // New KF (L1Trigger/TrackerTFP/plugins/ProducerKF.cc) is a new CMSSW emulation of KF firmware that
0016 // considers truncaton effect and is bit/clock cycle accurate.
0017 // With this version of the KF, HitPatternHelper makes predictions based on the spatial coordinates of tracks and detector modules.
0018 //
0019 //  Created by J.Li on 1/23/21.
0020 //
0021 
0022 #ifndef L1Trigger_TrackTrigger_interface_HitPatternHelper_h
0023 #define L1Trigger_TrackTrigger_interface_HitPatternHelper_h
0024 
0025 #include "FWCore/Framework/interface/data_default_record_trait.h"
0026 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0027 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0028 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0029 #include "Geometry/CommonTopologies/interface/PixelGeomDetUnit.h"
0030 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0031 #include "L1Trigger/TrackTrigger/interface/HitPatternHelperRcd.h"
0032 #include "L1Trigger/TrackTrigger/interface/Setup.h"
0033 
0034 #include <bitset>
0035 #include <iostream>
0036 #include <vector>
0037 #include <utility>
0038 #include <map>
0039 
0040 namespace hph {
0041 
0042   //Class that stores configuration for HitPatternHelper
0043   class Setup {
0044   public:
0045     Setup() {}
0046     Setup(const edm::ParameterSet& iConfig, const tt::Setup& setupTT);
0047     ~Setup() {}
0048 
0049     bool hphDebug() const { return iConfig_.getParameter<bool>("hphDebug"); }
0050     bool useNewKF() const { return iConfig_.getParameter<bool>("useNewKF"); }
0051     double deltaTanL() const { return iConfig_.getParameter<double>("deltaTanL"); }
0052     double chosenRofZ() const { return setupTT_.chosenRofZ(); }
0053     std::vector<double> etaRegions() const { return setupTT_.boundarieEta(); }
0054     std::vector<tt::SensorModule> sensorModules() const { return setupTT_.sensorModules(); }
0055     std::map<int, std::map<int, std::vector<int>>> layermap() const { return layermap_; }
0056     int nKalmanLayers() const { return nKalmanLayers_; }
0057     static auto smallerID(std::pair<int, bool> lhs, std::pair<int, bool> rhs) { return lhs.first < rhs.first; }
0058     static auto equalID(std::pair<int, bool> lhs, std::pair<int, bool> rhs) { return lhs.first == rhs.first; }
0059 
0060   private:
0061     edm::ParameterSet iConfig_;
0062     const tt::Setup setupTT_;  // Helper class to store TrackTrigger configuration
0063     std::vector<std::pair<int, bool>>
0064         layerIds_;  // layer IDs (1~6->L1~L6;11~15->D1~D5) and whether or not they are from tracker barrel
0065                     // Only needed by Old KF
0066     std::map<int, std::map<int, std::vector<int>>> layermap_;  // Hard-coded layermap in Old KF
0067     int nEtaRegions_;                                          // # of eta regions
0068     int nKalmanLayers_;                                        // # of maximum KF layers allowed
0069   };                                                           // Only needed by Old KF
0070 
0071   //Class that returns decoded information from hitpattern
0072   class HitPatternHelper {
0073   public:
0074     HitPatternHelper() {}
0075     HitPatternHelper(const Setup* setup, int hitpattern, double cot, double z0);
0076     ~HitPatternHelper() {}
0077 
0078     int etaSector() { return etaSector_; }      //Eta sectors defined in KF
0079     int numExpLayer() { return numExpLayer_; }  //The number of layers KF expects
0080     int numMissingPS() {
0081       return numMissingPS_;
0082     }  //The number of PS layers that are missing. It includes layers that are missing:
0083        //1)before the innermost stub on the track,
0084        //2)after the outermost stub on the track.
0085     int numMissing2S() {
0086       return numMissing2S_;
0087     }  //The number of 2S layers that are missing. It includes the two types of layers mentioned above.
0088     int numPS() { return numPS_; }  //The number of PS layers are found in hitpattern
0089     int num2S() { return num2S_; }  //The number of 2S layers are found in hitpattern
0090     int numMissingInterior1() {
0091       return numMissingInterior1_;
0092     }  //The number of missing interior layers (using only hitpattern)
0093     int numMissingInterior2() {
0094       return numMissingInterior2_;
0095     }  //The number of missing interior layers (using hitpattern, layermap from Old KF and sensor modules)
0096     std::vector<int> binary() { return binary_; }  //11-bit hitmask needed by TrackQuality.cc (0~5->L1~L6;6~10->D1~D5)
0097     static auto smallerID(tt::SensorModule lhs, tt::SensorModule rhs) { return lhs.layerId() < rhs.layerId(); }
0098     static auto equalID(tt::SensorModule lhs, tt::SensorModule rhs) { return lhs.layerId() == rhs.layerId(); }
0099 
0100     int ReducedId(
0101         int layerId);  //Converts layer ID (1~6->L1~L6;11~15->D1~D5) to reduced layer ID (0~5->L1~L6;6~10->D1~D5)
0102     int findLayer(int layerId);  //Search for a layer ID from sensor modules
0103 
0104   private:
0105     int etaSector_;
0106     int hitpattern_;
0107     int numExpLayer_;
0108     int numMissingLayer_;
0109     int numMissingPS_;
0110     int numMissing2S_;
0111     int numPS_;
0112     int num2S_;
0113     int numMissingInterior1_;
0114     int numMissingInterior2_;
0115     double cot_;
0116     double z0_;
0117     const Setup* setup_;
0118     std::vector<tt::SensorModule> layers_;  //Sensor modules that particles are expected to hit
0119     std::vector<int> binary_;
0120     bool hphDebug_;
0121     bool useNewKF_;
0122     float chosenRofZ_;
0123     float deltaTanL_;  // Uncertainty added to tanL (cot) when layermap in new KF is determined
0124     std::vector<double> etaRegions_;
0125     int nKalmanLayers_;
0126     std::map<int, std::map<int, std::vector<int>>> layermap_;
0127   };
0128 
0129 }  // namespace hph
0130 
0131 EVENTSETUP_DATA_DEFAULT_RECORD(hph::Setup, hph::SetupRcd);
0132 
0133 #endif