Cbrown kfout

L1Trigger/TrackFindingTracklet/plugins/ProducerKFout.cc

@@ -32,6 +32,8 @@ namespace trackFindingTracklet {
   public:
     explicit ProducerKFout(const ParameterSet&);
     ~ProducerKFout() override {}
+    template<typename T>
+    int digitise(const vector<T> Bins, T Value, T factor = 1 );
 
   private:
     void beginRun(const Run&, const EventSetup&) override;
@@ -58,6 +60,25 @@ namespace trackFindingTracklet {
     const Setup* setup_;
     // helper class to extract structured data from TTDTC::Frames
     const DataFormats* dataFormats_;
+    // Cot bins used to convert from cot to TanL
+    vector<double> cotBins_;
+    // Factors used to convert between phi/Z at a radius T, modified to include scale conversions needed in calculation
+    double modChosenRofZ_ = 0;
+    double modChosenRofPhi_ = 0;
+    // Corrections to Phi depending on which phi sector
+    double BaseSectorCorr_ = 0;
+    double UnsignedBaseSector_ = 0;
+    // Bins for dPhi/dZ use to access weight LUT below
+    vector<double> dPhiBins_;
+    vector<double> dZBins_;
+    // LUT for weighting functions for chi2 calculation
+    vector<double> v0Bins_;
+    vector<double> v1Bins_;
+    // Bins for final Chi2 Packing
+    vector<double> chi2rphiBins_;
+    vector<double> chi2rzBins_;
+
+    int maxTracksPerEvent_;
   };
 
   ProducerKFout::ProducerKFout(const ParameterSet& iConfig) :
@@ -92,6 +113,37 @@ namespace trackFindingTracklet {
       setup_->checkHistory(iRun.processHistory());
     // helper class to extract structured data from TTDTC::Frames
     dataFormats_ = &iSetup.getData(esGetTokenDataFormats_);
+
+    for (int i = 0; i < setup_->numSectorsEta(); i++) cotBins_.push_back(setup_->sectorCot(i));
+
+    modChosenRofZ_ = setup_->chosenRofZ();
+    modChosenRofPhi_ = setup_->hybridChosenRofPhi();
+
+    UnsignedBaseSector_ = (M_PI / (double)(setup_->numRegions() * setup_->numSectorsPhi()) );
+
+    // Convert Integer bins to doubles for internal calculation
+    for(size_t i = 0; i < setup_->kfoutdPhiBins().size(); i++)
+      dPhiBins_.push_back((double)setup_->kfoutdPhiBins()[i] * dataFormats_->base(Variable::dPhi, Process::kfin));
+    for(size_t i = 0; i < setup_->kfoutdZBins().size(); i++) 
+      dZBins_.push_back((double)setup_->kfoutdZBins()[i] * dataFormats_->base(Variable::dZ, Process::kfin));
+    for(size_t i = 0; i < setup_->kfoutv0Bins().size(); i++) 
+      v0Bins_.push_back(setup_->kfoutv0Bins()[i]*dataFormats_->base(Variable::dPhi, Process::kfin));
+    for(size_t i = 0; i < setup_->kfoutv1Bins().size(); i++) 
+      v1Bins_.push_back(setup_->kfoutv1Bins()[i]*dataFormats_->base(Variable::dZ, Process::kfin));
+    for(size_t i = 0; i < setup_->kfoutchi2rphiBins().size(); i++)
+      chi2rphiBins_.push_back((double)setup_->kfoutchi2rphiBins()[i] *setup_->kfoutchi2rphiConv() * pow(dataFormats_->base(Variable::dPhi, Process::kfin),3));
+    for(size_t i = 0; i < setup_->kfoutchi2rzBins().size(); i++)
+      chi2rzBins_.push_back((double)setup_->kfoutchi2rzBins()[i] *setup_->kfoutchi2rzConv() * pow(dataFormats_->base(Variable::dZ, Process::kfin),3));
+    maxTracksPerEvent_ = setup_->numFramesIO() * setup_->clockRatio();
+  }
+
+
+  template<typename T>
+  int ProducerKFout::digitise(const vector<T> Bins, T Value, T factor ) {
+    for (int i = 0; i < (int)Bins.size(); i++){
+      if (Value > Bins[i]*factor && Value <= Bins[i+1]*factor) {return i;}
+    }
+    return -1;
   }
 
   void ProducerKFout::produce(Event& iEvent, const EventSetup& iSetup) {
@@ -109,19 +161,139 @@ namespace trackFindingTracklet {
       Handle<TTTrackRefMap> handleTTTrackRefMap;
       iEvent.getByToken<TTTrackRefMap>(edGetTokenTTTrackRefMap_, handleTTTrackRefMap);
       const TTTrackRefMap& ttTrackRefMap = *handleTTTrackRefMap.product();
-      // perform KFout emulation and fill accepted and lost
-      // StreamsTrack is a vector of StreamTrack which is a vector of FrameTrack which is a pair of an edm::Ref<TTTrack> and std::bitset<64>
-      // the std::bitset<64> are the frames of an emp link
-      // the edm::Ref<TTTrack> is used to meassure tracking efficiency
-      // your input streamsTracks contain edm::Ref<TTTrack> to KF input TTTracks
-      // use ttTrackRefMap to lookup edm::Ref<TTTrack> of KF output TTTracks, that allows us to meassure the tracking efficiency after the KFout block
-      // your output frames belong to either only one TTTrack or to two, in the later case chose any edm::Ref<TTTrack> of the two
-    }
+      // 18 Output Links (First Vector) each has a vector of tracks per event (second vector) each track is 3 32 bit TTBV partial tracks 
+      vector<vector<TTBV>> SortedPartialTracks(setup_->numRegions() * setup_->tfpNumChannel(),vector<TTBV>(0));   
+      StreamsTrack OutputStreamsTracks(setup_->numRegions() * setup_->tfpNumChannel());
+
+      for (int iLink = 0; iLink < (int)streamsTracks.size(); iLink++ ){
+
+        for (int iTrack = 0; iTrack < (int)streamsTracks[iLink].size(); iTrack++ ){
+          auto track = streamsTracks[iLink].at(iTrack);
+          TrackKF InTrack(track,dataFormats_);
+
+          double temp_z0 = InTrack.zT() - ((InTrack.cot() * setup_->chosenRofZ()));
+
+          // Correction to Phi calcuation depending if +ve/-ve phi sector
+          if (InTrack.sectorPhi() == 0) BaseSectorCorr_ = -UnsignedBaseSector_;
+          else BaseSectorCorr_ = UnsignedBaseSector_;
+
+          double temp_phi0 = InTrack.phiT() - ((InTrack.inv2R()) * setup_->hybridChosenRofPhi()) + BaseSectorCorr_;
+
+          double temp_tanL = cotBins_[InTrack.sectorEta()] + InTrack.cot();
+
+          TTBV HitPattern(0,setup_->numLayers());
+
+          double tempchi2rphi = 0;
+          double tempchi2rz   = 0;
+
+          for (int iStub = 0; iStub < setup_->numLayers() - 1; iStub++ ){
+            auto stub = streamsStubs[setup_->numLayers()*iLink+iStub].at(iTrack);
+            StubKF InStub(stub,dataFormats_,iStub);
+
+            if (stub.first.isNonnull()){
+              HitPattern.set(iStub);
+
+              double phiSquared = InStub.phi() * InStub.phi();
+              double zSquared   = InStub.z() * InStub.z();
+
+              double tempv0 = (v0Bins_[digitise(dPhiBins_, InStub.dPhi())]);
+              double tempv1 = (v1Bins_[digitise(dZBins_  , InStub.dZ())]);
+
+              double tempRphi = phiSquared * tempv0;
+              double tempRz   = zSquared * tempv1;
+
+              tempchi2rphi += tempRphi;
+              tempchi2rz   += tempRz;
+            }
+          }
+          // TODO extract TTTrack bit widths from TTTrack word pending update to the TTTrack_word class
+          TTBV TrackValid(1,1,false);
+          TTBV extraMVA(0,6,false);
+          TTBV TQMVA(0,3,false);
+          TTBV BendChi2(0,3,false); 
+          TTBV Chi2rphi(digitise(chi2rphiBins_,tempchi2rphi),4,false);
+          TTBV Chi2rz(digitise(chi2rzBins_,tempchi2rz),4,false);
+          TTBV D0(0,13,false);
+          TTBV z0(temp_z0 ,dataFormats_->base(Variable::zT,Process::kf) ,12,true);
+          TTBV TanL(temp_tanL,dataFormats_->base(Variable::cot,Process::kf),16,true);
+          TTBV phi0(temp_phi0,dataFormats_->base(Variable::phiT,Process::kf),12,true);
+          TTBV InvR(-1*InTrack.inv2R(),dataFormats_->base(Variable::inv2R,Process::kf) ,16,true );
+
+                              // 6      +   3   +   7        +  3       + 13
+          TTBV PartialTrack1((extraMVA + TQMVA + HitPattern + BendChi2 + D0),32,false);
+                             // 4       + 12    + 16
+          TTBV PartialTrack2((Chi2rz   + z0    + TanL),32,false);
+                             // 4       + 12    +  15  +    1
+          TTBV PartialTrack3((Chi2rphi + phi0  + InvR + TrackValid),32,false);
+
+          // Sort Tracks based on eta
+          if (iLink % 2 == 0){
+            if (InTrack.sectorEta() < (int)(setup_->numSectorsEta()/2)){
+                SortedPartialTracks[iLink].push_back(PartialTrack1);
+                SortedPartialTracks[iLink].push_back(PartialTrack2);
+                SortedPartialTracks[iLink].push_back(PartialTrack3);
+                OutputStreamsTracks[iLink].emplace_back(track);
+              }
+              else{
+                SortedPartialTracks[iLink+1].push_back(PartialTrack1);
+                SortedPartialTracks[iLink+1].push_back(PartialTrack2);
+                SortedPartialTracks[iLink+1].push_back(PartialTrack3);
+                OutputStreamsTracks[iLink+1].emplace_back(track);
+              }
+            }
+          else{
+            if (InTrack.sectorEta() < (int)(setup_->numSectorsEta()/2)){
+              SortedPartialTracks[iLink-1].push_back(PartialTrack1);
+              SortedPartialTracks[iLink-1].push_back(PartialTrack2);
+              SortedPartialTracks[iLink-1].push_back(PartialTrack3);
+              OutputStreamsTracks[iLink-1].emplace_back(track);
+            }
+            else{
+              SortedPartialTracks[iLink].push_back(PartialTrack1);
+              SortedPartialTracks[iLink].push_back(PartialTrack2);
+              SortedPartialTracks[iLink].push_back(PartialTrack3);
+              OutputStreamsTracks[iLink].emplace_back(track);
+            }
+
+          }
+        } // Tracks
+      } // Links
+
+      // Fill products and match up tracks
+      TTBV NullBitTrack(0,32,false);
+      for (int iLink = 0; iLink < (int)OutputStreamsTracks.size(); iLink++ ){
+        // Iterate through partial tracks
+        int numLinkTracks = (int)OutputStreamsTracks[iLink].size();
+        if (numLinkTracks > 0){
+          if ((numLinkTracks % 2 != 0)) { //If there is an odd number of tracks 
+            SortedPartialTracks[iLink].push_back(NullBitTrack);  //Pad out final set of bits
+            OutputStreamsTracks[iLink].emplace_back(OutputStreamsTracks[iLink][numLinkTracks]); //Pad out with final repeated track
+            numLinkTracks++;
+            } //If there is an odd number of tracks 
+          for (int iTrack = 0; iTrack < (int)(SortedPartialTracks[iLink].size()); iTrack++ ){  
+            if (iTrack % 2 == 1){
+              TTTrackRef TrackRef;
+              for (auto &it : ttTrackRefMap) {   //Iterate throguh ttTrackRefMap to find TTTrackRef Key by a TTTrack Value
+                if(it.second == OutputStreamsTracks[iLink][(int)(iTrack-1)/3].first) { 
+                  TrackRef = it.first;
+                } 
+              }
+              if ((int)iTrack/3 <= maxTracksPerEvent_){
+                accepted[iLink].emplace_back(std::make_pair(TrackRef,(SortedPartialTracks[iLink][iTrack].slice(32) + SortedPartialTracks[iLink][iTrack-1].slice(32)).bs()));
+              }
+              else{
+                lost[iLink].emplace_back(std::make_pair(TrackRef,(SortedPartialTracks[iLink][iTrack].slice(32) + SortedPartialTracks[iLink][iTrack-1].slice(32)).bs()));
+              }
+            }
+          }
+        }
+      }
+
+    } // Config Supported
     // store products
     iEvent.emplace(edPutTokenAccepted_, move(accepted));
     iEvent.emplace(edPutTokenLost_, move(lost));
   }
-
 } // namespace trackFindingTracklet
 
-DEFINE_FWK_MODULE(trackFindingTracklet::ProducerKFout);
+DEFINE_FWK_MODULE(trackFindingTracklet::ProducerKFout);

L1Trigger/TrackTrigger/interface/Setup.h

@@ -116,7 +116,7 @@ namespace tt {
     //
     const std::vector<SensorModule>& sensorModules() const { return sensorModules_; }
 
-    // Fimrware specific Parameter
+    // Firmware specific Parameter
 
     // width of the 'A' port of an DSP slice
     int widthDSPa() const { return widthDSPa_; }
@@ -152,6 +152,9 @@ namespace tt {
     double maxdPhi() const { return maxdPhi_; }
     // maximum representable stub z uncertainty
     double maxdZ() const { return maxdZ_; }
+    // firmware clock / 360 MHz
+    double clockRatio() const { return clockRatio_; }
+
 
     // Common track finding parameter
 
@@ -447,6 +450,25 @@ namespace tt {
     // search window of each track parameter in initial uncertainties
     double kfRangeFactor() const { return kfRangeFactor_; }
 
+    // Parameter specifying KalmanFilter Output Formatter
+
+    // Bins used to digitize dPhi for chi2 calculation
+    std::vector<int> kfoutdPhiBins() const { return kfoutdPhiBins_; }
+    // Bins used to digitize dZ for chi2 calculation
+    std::vector<int> kfoutdZBins() const { return kfoutdZBins_; }
+    // v0 weight Bins corresponding to dPhi Bins for chi2 calculation
+    std::vector<int> kfoutv0Bins() const { return kfoutv0Bins_; }
+    // v1 weight Bins corresponding to dZ Bins for chi2 calculation
+    std::vector<int> kfoutv1Bins() const { return kfoutv1Bins_; }
+    // Final Chi2rphi digitization TODO extract from TTTrack Word 
+    std::vector<double> kfoutchi2rphiBins() const { return kfoutchi2rphiBins_; }
+    // Final Chi2rz digitization TODO extract from TTTrack Word 
+    std::vector<double> kfoutchi2rzBins() const { return kfoutchi2rzBins_; }
+    // Conversion factor between dphi^2/weight and chi2rphi 
+    int kfoutchi2rphiConv() const { return kfoutchi2rphiConv_; }
+    // Conversion factor between dz^2/weight and chi2rz
+    int kfoutchi2rzConv() const { return kfoutchi2rzConv_; }
+
     // Parameter specifying DuplicateRemoval
 
     // internal memory depth
@@ -606,7 +628,7 @@ namespace tt {
     // required number of associated ps layers to a TP to consider it reconstruct-able
     int tpMinLayersPS_;
 
-    // Fimrware specific Parameter
+    // Firmware specific Parameter
     edm::ParameterSet pSetFW_;
     // width of the 'A' port of an DSP slice
     int widthDSPa_;
@@ -668,6 +690,8 @@ namespace tt {
     double mindZ_;
     // maximum representable stub z uncertainty
     double maxdZ_;
+    // firmware clock / 360 MHz
+    double clockRatio_;
 
     // Parameter specifying front-end
     edm::ParameterSet pSetFE_;
@@ -803,6 +827,25 @@ namespace tt {
     // search window of each track parameter in initial uncertainties
     double kfRangeFactor_;
 
+    // Parameter specifying KalmanFilter Output Formatter
+    edm::ParameterSet pSetKFOut_;
+    // Bins used to digitize dPhi for chi2 calculation
+    std::vector<int> kfoutdPhiBins_; 
+    // Bins used to digitize dZ for chi2 calculation
+    std::vector<int> kfoutdZBins_;
+    // v0 weight Bins corresponding to dPhi Bins for chi2 calculation
+    std::vector<int> kfoutv0Bins_;
+    // v1 weight Bins corresponding to dZ Bins for chi2 calculation
+    std::vector<int> kfoutv1Bins_;
+    // Final Chi2rphi digitization TODO extract from TTTrack Word 
+    std::vector<double> kfoutchi2rphiBins_;
+    // Final Chi2rz digitization TODO extract from TTTrack Word 
+    std::vector<double> kfoutchi2rzBins_;
+    // Conversion factor between dphi^2/weight and chi2rphi 
+    int kfoutchi2rphiConv_;
+    // Conversion factor between dz^2/weight and chi2rz
+    int kfoutchi2rzConv_;
+
     // Parameter specifying DuplicateRemoval
     edm::ParameterSet pSetDR_;
     // internal memory depth

L1Trigger/TrackTrigger/python/ProducerSetup_cfi.py

@@ -103,7 +103,8 @@
     #MaxdPhi             = cms.double(  0.01              ), # maximum representable stub phi uncertainty
     MaxdPhi             = cms.double(  0.02              ), # maximum representable stub phi uncertainty
     MindZ               = cms.double(  0.1               ), # minimum representable stub z uncertainty
-    MaxdZ               = cms.double( 30.                )  # maximum representable stub z uncertainty
+    MaxdZ               = cms.double( 30.                ), # maximum representable stub z uncertainty
+    ClockRatio          = cms.double( 0.666              )  # firmware clock / 360 MHz
   ),
 
   # Parmeter specifying front-end
@@ -191,6 +192,20 @@
     MaxLayers   = cms.int32 ( 4   )  # maximum number of  layers added to a track
   ),
 
+  # Parmeter specifying KalmanFilter Output Formatter
+  KalmanFilterOut = cms.PSet (
+    dPhiBins = cms.vint32( 0,10,35,70,100,130,165,200,225,260,300,330,360,512 ),        # Bins used to digitize dPhi for chi2 calculation
+    dZBins   = cms.vint32( 0,100,150,233,266,300,333,366,400,433,466,500,533,566,1024), # Bins used to digitize dZ for chi2 calculation
+    v0Bins   = cms.vint32(32716, 7139, 209, 82, 43, 26, 18, 13, 9, 7, 6, 5, 4),         # v0 weight Bins corresponding to dPhi Bins for chi2 calculation
+    v1Bins   = cms.vint32(63607, 486, 155, 121, 97, 77, 63, 52, 44, 39, 33, 29, 26, 0), # v1 weight Bins corresponding to dZ Bins for chi2 calculation
+
+    chi2rphiBins = cms.vdouble( 0, 0.25, 0.5, 1, 2, 3, 5, 7, 10, 20, 40, 100, 200, 500, 1000, 3000,6000 ), # Final Chi2rphi digitization TODO extract from TTTrack Word 
+    chi2rzBins   = cms.vdouble( 0, 0.25, 0.5, 1, 2, 3, 5, 7, 10, 20, 40, 100, 200, 500, 1000, 3000,6000 ), # Final Chi2rz digitization TODO extract from TTTrack Word 
+
+    chi2rphiConv = cms.int32 ( 537600 ), # Conversion factor between dphi^2/weight and chi2rphi
+    chi2rzConv   = cms.int32 ( 641024 ), # Conversion factor between dz^2/weight and chi2rz
+  ),
+
   # Parmeter specifying DuplicateRemoval
   DuplicateRemoval = cms.PSet (
     DepthMemory  = cms.int32( 16 ) # internal memory depth

L1Trigger/TrackTrigger/src/Setup.cc

@@ -118,6 +118,7 @@ namespace tt {
         maxdPhi_(pSetFW_.getParameter<double>("MaxdPhi")),
         mindZ_(pSetFW_.getParameter<double>("MindZ")),
         maxdZ_(pSetFW_.getParameter<double>("MaxdZ")),
+        clockRatio_(pSetFW_.getParameter<double>("ClockRatio")),
         // Parmeter specifying front-end
         pSetFE_(iConfig.getParameter<ParameterSet>("FrontEnd")),
         widthBend_(pSetFE_.getParameter<int>("WidthBend")),
@@ -187,6 +188,16 @@ namespace tt {
         kfMinLayers_(pSetKF_.getParameter<int>("MinLayers")),
         kfMaxLayers_(pSetKF_.getParameter<int>("MaxLayers")),
         kfRangeFactor_(pSetKF_.getParameter<double>("RangeFactor")),
+        // Parmeter specifying KalmanFilter Output Formatter
+        pSetKFOut_(iConfig.getParameter<ParameterSet>("KalmanFilterOut")),
+        kfoutdPhiBins_(pSetKFOut_.getParameter<vector<int>>("dPhiBins")),
+        kfoutdZBins_(pSetKFOut_.getParameter<vector<int>>("dZBins")),
+        kfoutv0Bins_(pSetKFOut_.getParameter<vector<int>>("v0Bins")),
+        kfoutv1Bins_(pSetKFOut_.getParameter<vector<int>>("v1Bins")),
+        kfoutchi2rphiBins_(pSetKFOut_.getParameter<vector<double>>("chi2rphiBins")),
+        kfoutchi2rzBins_(pSetKFOut_.getParameter<vector<double>>("chi2rzBins")),
+        kfoutchi2rphiConv_(pSetKFOut_.getParameter<int>("chi2rphiConv")),
+        kfoutchi2rzConv_(pSetKFOut_.getParameter<int>("chi2rzConv")),
         // Parmeter specifying DuplicateRemoval
         pSetDR_(iConfig.getParameter<ParameterSet>("DuplicateRemoval")),
         drDepthMemory_(pSetDR_.getParameter<int>("DepthMemory")) {