TrackFindingTracklet/plugins/ProducerKFout.cc

0001 #include "FWCore/Framework/interface/stream/EDProducer.h"
0002 #include "FWCore/Framework/interface/Run.h"
0003 #include "FWCore/Framework/interface/EventSetup.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/Utilities/interface/EDGetToken.h"
0007 #include "FWCore/Utilities/interface/EDPutToken.h"
0008 #include "FWCore/Utilities/interface/ESGetToken.h"
0009 #include "FWCore/Utilities/interface/InputTag.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "DataFormats/Common/interface/Handle.h"
0012
0013 #include "L1Trigger/TrackTrigger/interface/Setup.h"
0014 #include "L1Trigger/TrackerTFP/interface/DataFormats.h"
0015 #include "L1Trigger/TrackTrigger/interface/L1TrackQuality.h"
0016
0017 #include <string>
0018 #include <numeric>
0019
0020 using namespace std;
0021 using namespace edm;
0022 using namespace trackerTFP;
0023 using namespace tt;
0024
0025 namespace trklet {
0026
0027   /*! \class  trklet::ProducerKFout
0028    *  \brief  Converts KF output into TFP output
0029    *  A bit accurate emulation of the track transformation, the
0030    *  eta routing and splitting of the 96-bit track words into 64-bit
0031    *  packets. Also run is a bit accurate emulation of the track quality
0032    *  BDT, whose output is also added to the track word.
0033    *  \author Christopher Brown
0034    *  \date   2021, Aug
0035    */
0036   class ProducerKFout : public stream::EDProducer<> {
0037   public:
0038     explicit ProducerKFout(const ParameterSet&);
0039     ~ProducerKFout() override {}
0040
0041   private:
0042     void beginRun(const Run&, const EventSetup&) override;
0043     void produce(Event&, const EventSetup&) override;
0044     void endJob() {}
0045
0046     // ED input token of kf stubs
0047     EDGetTokenT<StreamsStub> edGetTokenStubs_;
0048     // ED input token of kf tracks
0049     EDGetTokenT<StreamsTrack> edGetTokenTracks_;
0050     // ED input token of kf input to kf output TTTrack map
0051     EDGetTokenT<TTTrackRefMap> edGetTokenTTTrackRefMap_;
0052     // ED output token for accepted kfout tracks
0053     EDPutTokenT<StreamsTrack> edPutTokenAccepted_;
0054     // ED output token for truncated kfout tracks
0055     EDPutTokenT<StreamsTrack> edPutTokenLost_;
0056     // Setup token
0057     ESGetToken<Setup, SetupRcd> esGetTokenSetup_;
0058     // DataFormats token
0059     ESGetToken<DataFormats, DataFormatsRcd> esGetTokenDataFormats_;
0060     // configuration
0061     ParameterSet iConfig_;
0062     // helper class to store configurations
0063     const Setup* setup_;
0064     // helper class to extract structured data from tt::Frames
0065     const DataFormats* dataFormats_;
0066     // Bins for dPhi/dZ use to create weight LUT
0067     vector<double> dPhiBins_;
0068     vector<double> dZBins_;
0069
0070     std::unique_ptr<L1TrackQuality> trackQualityModel_;
0071     vector<int> tqBins_;
0072     double tqTanlScale_;
0073     double tqZ0Scale_;
0074     static constexpr double ap_fixed_rescale = 32.0;
0075
0076     // For convenience and keeping readable code, accessed many times
0077     int numWorkers_;
0078
0079     int partialTrackWordBits_;
0080
0081     // Helper function to convert floating value to bin
0082     template <typename T>
0083     unsigned int digitise(const T& bins, double value, double factor) {
0084       unsigned int bin = 0;
0085       for (unsigned int i = 0; i < bins.size() - 1; i++) {
0086         if (value * factor > bins[i] && value * factor <= bins[i + 1])
0087           break;
0088         bin++;
0089       }
0090       return bin;
0091     }
0092   };
0093
0094   ProducerKFout::ProducerKFout(const ParameterSet& iConfig) : iConfig_(iConfig) {
0095     const string& labelKF = iConfig.getParameter<string>("LabelKF");
0096     const string& labelAS = iConfig.getParameter<string>("LabelAS");
0097     const string& branchStubs = iConfig.getParameter<string>("BranchAcceptedStubs");
0098     const string& branchTracks = iConfig.getParameter<string>("BranchAcceptedTracks");
0099     const string& branchLost = iConfig.getParameter<string>("BranchLostTracks");
0100     // book in- and output ED products
0101     edGetTokenStubs_ = consumes<StreamsStub>(InputTag(labelKF, branchStubs));
0102     edGetTokenTracks_ = consumes<StreamsTrack>(InputTag(labelKF, branchTracks));
0103     edGetTokenTTTrackRefMap_ = consumes<TTTrackRefMap>(InputTag(labelAS, branchTracks));
0104     edPutTokenAccepted_ = produces<StreamsTrack>(branchTracks);
0105     edPutTokenLost_ = produces<StreamsTrack>(branchLost);
0106     // book ES products
0107     esGetTokenSetup_ = esConsumes<Setup, SetupRcd, Transition::BeginRun>();
0108     esGetTokenDataFormats_ = esConsumes<DataFormats, DataFormatsRcd, Transition::BeginRun>();
0109     // initial ES products
0110     setup_ = nullptr;
0111     dataFormats_ = nullptr;
0112
0113     trackQualityModel_ = std::make_unique<L1TrackQuality>(iConfig.getParameter<edm::ParameterSet>("TrackQualityPSet"));
0114     edm::ParameterSet trackQualityPSset = iConfig.getParameter<edm::ParameterSet>("TrackQualityPSet");
0115     tqBins_ = trackQualityPSset.getParameter<vector<int>>("tqemu_bins");
0116     tqTanlScale_ = trackQualityPSset.getParameter<double>("tqemu_TanlScale");
0117     tqZ0Scale_ = trackQualityPSset.getParameter<double>("tqemu_Z0Scale");
0118   }
0119
0120   void ProducerKFout::beginRun(const Run& iRun, const EventSetup& iSetup) {
0121     // helper class to store configurations
0122     setup_ = &iSetup.getData(esGetTokenSetup_);
0123     if (!setup_->configurationSupported())
0124       return;
0125     // check process history if desired
0126     if (iConfig_.getParameter<bool>("CheckHistory"))
0127       setup_->checkHistory(iRun.processHistory());
0128     // helper class to extract structured data from tt::Frames
0129     dataFormats_ = &iSetup.getData(esGetTokenDataFormats_);
0130
0131     // Calculate 1/dz**2 and 1/dphi**2 bins for v0 and v1 weightings
0132
0133     float temp_dphi = 0.0;
0134     float temp_dz = 0.0;
0135     for (int i = 0;
0136          i < pow(2, dataFormats_->width(Variable::dPhi, Process::kfin)) / pow(2, setup_->weightBinFraction());
0137          i++) {
0138       temp_dphi =
0139           pow(dataFormats_->base(Variable::dPhi, Process::kfin) * (i + 1) * pow(2, setup_->weightBinFraction()), -2);
0140       temp_dphi = temp_dphi / setup_->dphiTruncation();
0141       temp_dphi = std::floor(temp_dphi);
0142       dPhiBins_.push_back(temp_dphi * setup_->dphiTruncation());
0143     }
0144     for (int i = 0; i < pow(2, dataFormats_->width(Variable::dZ, Process::kfin)) / pow(2, setup_->weightBinFraction());
0145          i++) {
0146       temp_dz =
0147           pow(dataFormats_->base(Variable::dZ, Process::kfin) * (i + 1) * pow(2, setup_->weightBinFraction()), -2);
0148       temp_dz = temp_dz * setup_->dzTruncation();
0149       temp_dz = std::ceil(temp_dz);
0150       dZBins_.push_back(temp_dz / setup_->dzTruncation());
0151     }
0152     numWorkers_ = setup_->kfNumWorker();
0153     partialTrackWordBits_ = TTBV::S_ / 2;
0154   }
0155
0156   void ProducerKFout::produce(Event& iEvent, const EventSetup& iSetup) {
0157     // empty KFout product
0158     StreamsTrack accepted(setup_->numRegions() * setup_->tfpNumChannel());
0159     StreamsTrack lost(setup_->numRegions() * setup_->tfpNumChannel());
0160     // read in KF Product and produce KFout product
0161     if (setup_->configurationSupported()) {
0162       Handle<StreamsStub> handleStubs;
0163       iEvent.getByToken<StreamsStub>(edGetTokenStubs_, handleStubs);
0164       const StreamsStub& streamsStubs = *handleStubs.product();
0165       Handle<StreamsTrack> handleTracks;
0166       iEvent.getByToken<StreamsTrack>(edGetTokenTracks_, handleTracks);
0167       const StreamsTrack& streamsTracks = *handleTracks.product();
0168       Handle<TTTrackRefMap> handleTTTrackRefMap;
0169       iEvent.getByToken<TTTrackRefMap>(edGetTokenTTTrackRefMap_, handleTTTrackRefMap);
0170       const TTTrackRefMap& ttTrackRefMap = *handleTTTrackRefMap.product();
0171       // 18 Output Links (First Vector) each has a vector of tracks per event (second vector) each track is 3 32 bit TTBV partial tracks
0172       vector<vector<TTBV>> sortedPartialTracks(setup_->numRegions() * setup_->tfpNumChannel(), vector<TTBV>(0));
0173
0174       TrackKFOutSAPtrCollectionss inTrackStreams;
0175       TrackKFOutSAPtrCollectionss outTrackStreams;
0176
0177       // Setup empty collections for input tracks to be routed
0178       for (int iRegion = 0; iRegion < setup_->numRegions(); iRegion++) {
0179         TrackKFOutSAPtrCollections temp_collection;
0180         for (int iLink = 0; iLink < setup_->tfpNumChannel(); iLink++) {
0181           TrackKFOutSAPtrCollection temp;
0182           for (int iTrack = 0; iTrack < setup_->numFramesIO(); iTrack++)
0183             temp.emplace_back(std::make_shared<TrackKFOut>());
0184           temp_collection.push_back(temp);
0185         }
0186         outTrackStreams.push_back(temp_collection);
0187       }
0188
0189       // Setup empty collections for output tracks from routing
0190       for (int iRegion = 0; iRegion < setup_->numRegions(); iRegion++) {
0191         TrackKFOutSAPtrCollections temp_collection;
0192         for (int iLink = 0; iLink < numWorkers_; iLink++) {
0193           TrackKFOutSAPtrCollection temp;
0194           for (int iTrack = 0; iTrack < setup_->numFramesIO(); iTrack++)
0195             temp.emplace_back(std::make_shared<TrackKFOut>());
0196           temp_collection.push_back(temp);
0197         }
0198         inTrackStreams.push_back(temp_collection);
0199       }
0200
0201       StreamsTrack outputStreamsTracks(setup_->numRegions() * setup_->tfpNumChannel());
0202
0203       // Setup containers for track quality
0204       float tempTQMVAPreSig = 0.0;
0205       // Due to ap_fixed implementation in CMSSW this 10,5 must be specified at compile time, TODO make this a changeable parameter
0206       std::vector<ap_fixed<10, 5>> trackQuality_inputs = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
0207
0208       for (int iLink = 0; iLink < (int)streamsTracks.size(); iLink++) {
0209         for (int iTrack = 0; iTrack < (int)streamsTracks[iLink].size(); iTrack++) {
0210           const auto& track = streamsTracks[iLink].at(iTrack);
0211           TrackKF inTrack(track, dataFormats_);
0212
0213           double temp_z0 = inTrack.zT() - ((inTrack.cot() * setup_->chosenRofZ()));
0214
0215           // Correction to Phi calcuation depending if +ve/-ve phi sector
0216           const double baseSectorCorr = inTrack.sectorPhi() ? -setup_->baseSector() : setup_->baseSector();
0217
0218           double temp_phi0 = inTrack.phiT() - ((inTrack.inv2R()) * setup_->hybridChosenRofPhi()) + baseSectorCorr;
0219
0220           double temp_tanL = inTrack.cotGlobal();
0221
0222           TTBV hitPattern(0, setup_->numLayers());
0223
0224           double tempchi2rphi = 0;
0225           double tempchi2rz = 0;
0226
0227           int temp_nstub = 0;
0228           int temp_ninterior = 0;
0229           bool counter = false;
0230
0231           for (int iStub = 0; iStub < setup_->numLayers() - 1; iStub++) {
0232             const auto& stub = streamsStubs[setup_->numLayers() * iLink + iStub].at(iTrack);
0233             StubKF inStub(stub, dataFormats_, iStub);
0234
0235             if (!stub.first.isNonnull()) {
0236               if (counter)
0237                 temp_ninterior += 1;
0238               continue;
0239             }
0240
0241             counter = true;
0242
0243             hitPattern.set(iStub);
0244             temp_nstub += 1;
0245             double phiSquared = pow(inStub.phi(), 2);
0246             double zSquared = pow(inStub.z(), 2);
0247
0248             double tempv0 = dPhiBins_[(inStub.dPhi() / (dataFormats_->base(Variable::dPhi, Process::kfin) *
0249                                                         pow(2, setup_->weightBinFraction())))];
0250             double tempv1 = dZBins_[(
0251                 inStub.dZ() / (dataFormats_->base(Variable::dZ, Process::kfin) * pow(2, setup_->weightBinFraction())))];
0252
0253             double tempRphi = phiSquared * tempv0;
0254             double tempRz = zSquared * tempv1;
0255
0256             tempchi2rphi += tempRphi;
0257             tempchi2rz += tempRz;
0258           }  // Iterate over track stubs
0259
0260           // Create bit vectors for eacch output, including digitisation of chi2
0261           // TODO implement extraMVA, bendChi2, d0
0262           TTBV trackValid(1, TTTrack_TrackWord::TrackBitWidths::kValidSize, false);
0263           TTBV extraMVA(0, TTTrack_TrackWord::TrackBitWidths::kMVAOtherSize, false);
0264           TTBV bendChi2(0, TTTrack_TrackWord::TrackBitWidths::kBendChi2Size, false);
0265           TTBV chi2rphi(digitise(TTTrack_TrackWord::chi2RPhiBins, tempchi2rphi, (double)setup_->kfoutchi2rphiConv()),
0266                         TTTrack_TrackWord::TrackBitWidths::kChi2RPhiSize,
0267                         false);
0268           TTBV chi2rz(digitise(TTTrack_TrackWord::chi2RZBins, tempchi2rz, (double)setup_->kfoutchi2rzConv()),
0269                       TTTrack_TrackWord::TrackBitWidths::kChi2RZSize,
0270                       false);
0271           TTBV d0(0, TTTrack_TrackWord::TrackBitWidths::kD0Size, false);
0272           TTBV z0(
0273               temp_z0, dataFormats_->base(Variable::zT, Process::kf), TTTrack_TrackWord::TrackBitWidths::kZ0Size, true);
0274           TTBV tanL(temp_tanL,
0275                     dataFormats_->base(Variable::cot, Process::kf),
0276                     TTTrack_TrackWord::TrackBitWidths::kTanlSize,
0277                     true);
0278           TTBV phi0(temp_phi0,
0279                     dataFormats_->base(Variable::phiT, Process::kf),
0280                     TTTrack_TrackWord::TrackBitWidths::kPhiSize,
0281                     true);
0282           TTBV invR(-inTrack.inv2R(),
0283                     dataFormats_->base(Variable::inv2R, Process::kf),
0284                     TTTrack_TrackWord::TrackBitWidths::kRinvSize + 1,
0285                     true);
0286           invR.resize(TTTrack_TrackWord::TrackBitWidths::kRinvSize);
0287
0288           // Create input vector for BDT
0289           trackQuality_inputs = {
0290               (std::trunc(tanL.val() / tqTanlScale_)) / ap_fixed_rescale,
0291               (std::trunc(z0.val() / tqZ0Scale_)) / ap_fixed_rescale,
0292               0,
0293               temp_nstub,
0294               temp_ninterior,
0295               digitise(TTTrack_TrackWord::chi2RPhiBins, tempchi2rphi, (double)setup_->kfoutchi2rphiConv()),
0296               digitise(TTTrack_TrackWord::chi2RZBins, tempchi2rz, (double)setup_->kfoutchi2rzConv())};
0297
0298           // Run BDT emulation and package output into 3 bits
0299           // output needs sigmoid transformation applied
0300           tempTQMVAPreSig = trackQualityModel_->runEmulatedTQ(trackQuality_inputs);
0301           TTBV tqMVA(digitise(L1TrackQuality::getTqMVAPreSigBins(), tempTQMVAPreSig, 1.0),
0302                      TTTrack_TrackWord::TrackBitWidths::kMVAQualitySize,
0303                      false);
0304
0305           // Build 32 bit partial tracks for outputting in 64 bit packets
0306           //                  12 +  3       +  7         +  3    +  6
0307           TTBV partialTrack3((d0 + bendChi2 + hitPattern + tqMVA + extraMVA), partialTrackWordBits_, false);
0308           //                  16   + 12    + 4
0309           TTBV partialTrack2((tanL + z0 + chi2rz), partialTrackWordBits_, false);
0310           //                    1        + 15   +  12 +    4
0311           TTBV partialTrack1((trackValid + invR + phi0 + chi2rphi), partialTrackWordBits_, false);
0312
0313           int sortKey = (inTrack.sectorEta() < (int)(setup_->numSectorsEta() / 2)) ? 0 : 1;
0314           // Set correct bit to valid for track valid
0315           TrackKFOut temp_track(partialTrack1.set((partialTrackWordBits_ - 1)),
0316                                 partialTrack2,
0317                                 partialTrack3,
0318                                 sortKey,
0319                                 track,
0320                                 iTrack,
0321                                 iLink,
0322                                 true);
0323
0324           inTrackStreams[iLink / setup_->kfNumWorker()][iLink % setup_->kfNumWorker()][iTrack] =
0325               (std::make_shared<TrackKFOut>(temp_track));
0326         }  // Iterate over Tracks
0327       }    // Iterate over Links
0328            // Route Tracks in eta based on their sort key
0329       for (int iRegion = 0; iRegion < setup_->numRegions(); iRegion++) {
0330         int buffered_tracks[] = {0, 0};
0331         for (int iTrack = 0; iTrack < setup_->numFramesIO() * ((double)TTBV::S_ / TTTrack_TrackWord::kTrackWordSize);
0332              iTrack++) {
0333           for (int iWorker = 0; iWorker < setup_->kfNumWorker(); iWorker++) {
0334             for (int iLink = 0; iLink < setup_->tfpNumChannel(); iLink++) {
0335               if ((inTrackStreams[iRegion][iWorker][iTrack]->sortKey() == iLink) &&
0336                   (inTrackStreams[iRegion][iWorker][iTrack]->dataValid() == true)) {
0337                 outTrackStreams[iRegion][iLink][buffered_tracks[iLink]] = inTrackStreams[iRegion][iWorker][iTrack];
0338                 buffered_tracks[iLink] = buffered_tracks[iLink] + 1;
0339               }
0340             }
0341           }
0342         }
0343       }
0344
0345       // Pack output of router onto each link, with correct partial tracks in correct places
0346       for (int iRegion = 0; iRegion < setup_->numRegions(); iRegion++) {
0347         for (int iLink = 0; iLink < setup_->tfpNumChannel(); iLink++) {
0348           for (int iTrack = 0; iTrack < (int)outTrackStreams[iRegion][iLink].size(); iTrack++) {
0349             sortedPartialTracks[2 * iRegion + iLink].push_back(
0350                 outTrackStreams[iRegion][iLink][iTrack]->PartialTrack1());
0351             sortedPartialTracks[2 * iRegion + iLink].push_back(
0352                 outTrackStreams[iRegion][iLink][iTrack]->PartialTrack2());
0353             sortedPartialTracks[2 * iRegion + iLink].push_back(
0354                 outTrackStreams[iRegion][iLink][iTrack]->PartialTrack3());
0355             outputStreamsTracks[2 * iRegion + iLink].emplace_back(outTrackStreams[iRegion][iLink][iTrack]->track());
0356           }
0357         }
0358       }
0359       // Fill products and match up tracks
0360       // store products
0361       const TTBV nullBitTrack(0, partialTrackWordBits_, false);
0362       for (int iLink = 0; iLink < (int)outputStreamsTracks.size(); iLink++) {
0363         // Iterate through partial tracks
0364         int numLinkTracks = (int)outputStreamsTracks[iLink].size();
0365         if (numLinkTracks == 0)
0366           continue;  // Don't fill links if no tracks
0367         if ((numLinkTracks % 2 != 0)) {
0368           sortedPartialTracks[iLink].push_back(nullBitTrack);  //Pad out final set of bits
0369           outputStreamsTracks[iLink].emplace_back(
0370               outputStreamsTracks[iLink][numLinkTracks++]);  //Pad out with final repeated track
0371         }                                                    //If there is an odd number of tracks
0372         for (int iTrack = 0; iTrack < (int)(sortedPartialTracks[iLink].size()); iTrack++) {
0373           if (iTrack % 2 != 1)  // Write to links every other partial track, 3 partial tracks per full TTTrack
0374             continue;
0375           TTTrackRef trackRef;
0376           for (auto& it : ttTrackRefMap) {  //Iterate through ttTrackRefMap to find TTTrackRef Key by a TTTrack Value
0377             if (it.second == outputStreamsTracks[iLink][(int)(iTrack - 1) / 3].first)
0378               trackRef = it.first;
0379           }
0380           if ((int)iTrack / 3 <= setup_->numFramesIO() * ((double)TTBV::S_ / TTTrack_TrackWord::kTrackWordSize))
0381             accepted[iLink].emplace_back(
0382                 std::make_pair(trackRef,
0383                                (sortedPartialTracks[iLink][iTrack - 1].slice(partialTrackWordBits_) +
0384                                 sortedPartialTracks[iLink][iTrack].slice(partialTrackWordBits_))
0385                                    .bs()));
0386           else
0387             lost[iLink].emplace_back(
0388                 std::make_pair(trackRef,
0389                                (sortedPartialTracks[iLink][iTrack - 1].slice(partialTrackWordBits_) +
0390                                 sortedPartialTracks[iLink][iTrack].slice(partialTrackWordBits_))
0391                                    .bs()));
0392         }  //Iterate through sorted partial tracks
0393       }    // Iterate through links
0394     }      // Config Supported
0395     // store products
0396     iEvent.emplace(edPutTokenAccepted_, std::move(accepted));
0397     iEvent.emplace(edPutTokenLost_, std::move(lost));
0398   }
0399 }  // namespace trklet
0400
0401 DEFINE_FWK_MODULE(trklet::ProducerKFout);