Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-01 23:54:07

0001 #ifndef PHASE2GMT_TRACKCONVERTER
0002 #define PHASE2GMT_TRACKCONVERTER
0003 
0004 #include "L1Trigger/Phase2L1GMT/interface/ConvertedTTTrack.h"
0005 #include "DataFormats/L1TMuonPhase2/interface/TrackerMuon.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 
0008 namespace Phase2L1GMT {
0009 
0010   class TrackConverter {
0011   public:
0012     TrackConverter(const edm::ParameterSet& iConfig) : verbose_(iConfig.getParameter<int>("verbose")) {}
0013     ~TrackConverter() {}
0014 
0015     std::vector<ConvertedTTTrack> convertTracks(const std::vector<edm::Ptr<l1t::TrackerMuon::L1TTTrackType> >& tracks) {
0016       std::vector<ConvertedTTTrack> out;
0017       out.reserve(tracks.size());
0018       for (const auto& t : tracks)
0019         out.push_back(convert(t));
0020       return out;
0021     }
0022 
0023   private:
0024     int verbose_;
0025     typedef ap_uint<96> wordtype;
0026 
0027     uint generateQuality(const edm::Ptr<TTTrack<Ref_Phase2TrackerDigi_> >& track) { return 1; }
0028 
0029     uint ptLookup(uint absCurv) {
0030       for (auto i : ptShifts) {
0031         if (absCurv >= uint(i[0]) && absCurv < uint(i[1])) {
0032           if (i[2] < 0)
0033             return i[4];
0034           else
0035             return (absCurv >> i[2]) + i[3];
0036         }
0037       }
0038       return 0;
0039     }
0040 
0041     uint etaLookup(uint absTanL) {
0042       for (auto i : etaShifts) {
0043         if (absTanL >= uint(i[0]) && absTanL < uint(i[1])) {
0044           if (i[2] < 0)
0045             return i[4];
0046           else
0047             return (absTanL >> i[2]) + i[3];
0048         }
0049       }
0050       return 0;
0051     }
0052 
0053     ConvertedTTTrack convert(const edm::Ptr<TTTrack<Ref_Phase2TrackerDigi_> >& track) {
0054       uint charge = (track->rInv() < 0) ? 1 : 0;
0055       int curvature = track->rInv() * (1 << (BITSTTCURV - 1)) / maxCurv_;
0056       int phi = track->phi() * (1 << (BITSPHI - 1)) / (M_PI);
0057       int tanLambda = track->tanL() * (1 << (BITSTTTANL - 1)) / maxTanl_;
0058       int z0 = track->z0() * (1 << (BITSZ0 - 1)) / maxZ0_;
0059       int d0 = track->d0() * (1 << (BITSD0 - 1)) / maxD0_;
0060       //calculate pt
0061       uint absCurv = curvature > 0 ? (curvature) : (-curvature);
0062       uint pt = ptLUT[ptLookup(absCurv)];
0063       uint quality = generateQuality(track);
0064       uint absTanL = tanLambda > 0 ? (tanLambda) : (-tanLambda);
0065       uint absEta = etaLUT[etaLookup(absTanL)];
0066       int eta = tanLambda > 0 ? (absEta) : (-absEta);
0067 
0068       ap_int<BITSPHI> phiSec = ap_int<BITSPHI>(phi) -
0069                                ap_int<BITSPHI>((track->phiSector() * 40 * M_PI / 180.) * (1 << (BITSPHI - 1)) / (M_PI));
0070       ap_int<BITSPHI> phiCorrected = ap_int<BITSPHI>(phiSec + track->phiSector() * 910);
0071 
0072       wordtype word = 0;
0073       int bstart = 0;
0074       bstart = wordconcat<wordtype>(word, bstart, curvature, BITSTTCURV);
0075       bstart = wordconcat<wordtype>(word, bstart, phiSec, BITSTTPHI);
0076       bstart = wordconcat<wordtype>(word, bstart, tanLambda, BITSTTTANL);
0077       bstart = wordconcat<wordtype>(word, bstart, z0, BITSZ0);
0078       bstart = wordconcat<wordtype>(word, bstart, d0, BITSD0);
0079       bstart = wordconcat<wordtype>(word, bstart, uint(track->chi2()), 4);
0080 
0081       ConvertedTTTrack convertedTrack(charge, curvature, absEta, pt, eta, phiCorrected.to_int(), z0, d0, quality, word);
0082       convertedTrack.setOfflineQuantities(track->momentum().transverse(), track->eta(), track->phi());
0083       if (verbose_)
0084         convertedTrack.print();
0085       convertedTrack.setTrkPtr(track);
0086       return convertedTrack;
0087     }
0088   };
0089 }  // namespace Phase2L1GMT
0090 
0091 #endif