TrackFindingTracklet/src/KFin.cc

0001 #include "L1Trigger/TrackFindingTracklet/interface/KFin.h"
0002
0003 #include <vector>
0004 #include <numeric>
0005 #include <algorithm>
0006
0007 using namespace std;
0008 using namespace edm;
0009 using namespace tt;
0010 using namespace trackerTFP;
0011
0012 namespace trklet {
0013
0014   KFin::KFin(const ParameterSet& iConfig,
0015              const Setup* setup,
0016              const DataFormats* dataFormats,
0017              const LayerEncoding* layerEncoding,
0018              const ChannelAssignment* channelAssignment,
0019              int region)
0020       : enableTruncation_(iConfig.getParameter<bool>("EnableTruncation")),
0021         setup_(setup),
0022         dataFormats_(dataFormats),
0023         layerEncoding_(layerEncoding),
0024         channelAssignment_(channelAssignment),
0025         region_(region),
0026         input_(channelAssignment_->numNodesDR()) {}
0027
0028   // read in and organize input tracks and stubs
0029   void KFin::consume(const StreamsTrack& streamsTrack, const StreamsStub& streamsStub) {
0030     const int offsetTrack = region_ * channelAssignment_->numNodesDR();
0031     auto nonNullTrack = [](int sum, const FrameTrack& frame) { return sum + (frame.first.isNonnull() ? 1 : 0); };
0032     auto nonNullStub = [](int sum, const FrameStub& frame) { return sum + (frame.first.isNonnull() ? 1 : 0); };
0033     // count tracks and stubs and reserve corresponding vectors
0034     int sizeTracks(0);
0035     int sizeStubs(0);
0036     for (int channel = 0; channel < channelAssignment_->numNodesDR(); channel++) {
0037       const int streamTrackId = offsetTrack + channel;
0038       const int offsetStub = streamTrackId * setup_->numLayers();
0039       const StreamTrack& streamTrack = streamsTrack[streamTrackId];
0040       input_[channel].reserve(streamTrack.size());
0041       sizeTracks += accumulate(streamTrack.begin(), streamTrack.end(), 0, nonNullTrack);
0042       for (int layer = 0; layer < setup_->numLayers(); layer++) {
0043         const StreamStub& streamStub = streamsStub[offsetStub + layer];
0044         sizeStubs += accumulate(streamStub.begin(), streamStub.end(), 0, nonNullStub);
0045       }
0046     }
0047     tracks_.reserve(sizeTracks);
0048     stubs_.reserve(sizeStubs);
0049     // transform input data into handy structs
0050     for (int channel = 0; channel < channelAssignment_->numNodesDR(); channel++) {
0051       vector<Track*>& input = input_[channel];
0052       const int streamTrackId = offsetTrack + channel;
0053       const int offsetStub = streamTrackId * setup_->numLayers();
0054       const StreamTrack& streamTrack = streamsTrack[streamTrackId];
0055       for (int frame = 0; frame < (int)streamTrack.size(); frame++) {
0056         const FrameTrack& frameTrack = streamTrack[frame];
0057         if (frameTrack.first.isNull()) {
0058           input.push_back(nullptr);
0059           continue;
0060         }
0061         vector<Stub*> stubs;
0062         stubs.reserve(setup_->numLayers());
0063         for (int layer = 0; layer < setup_->numLayers(); layer++) {
0064           const FrameStub& frameStub = streamsStub[offsetStub + layer][frame];
0065           if (frameStub.first.isNull())
0066             continue;
0067           TTBV ttBV = frameStub.second;
0068           const TTBV zBV(ttBV, dataFormats_->format(Variable::z, Process::kfin).width(), 0, true);
0069           ttBV >>= dataFormats_->format(Variable::z, Process::kfin).width();
0070           const TTBV phiBV(ttBV, dataFormats_->format(Variable::phi, Process::kfin).width(), 0, true);
0071           ttBV >>= dataFormats_->format(Variable::phi, Process::kfin).width();
0072           const TTBV rBV(ttBV, dataFormats_->format(Variable::r, Process::kfin).width(), 0, true);
0073           ttBV >>= dataFormats_->format(Variable::r, Process::kfin).width();
0074           const TTBV layerIdBV(ttBV, channelAssignment_->widthLayerId(), 0);
0075           ttBV >>= channelAssignment_->widthPSTilt();
0076           const TTBV tiltBV(ttBV, channelAssignment_->widthPSTilt(), 0);
0077           const double r = rBV.val(dataFormats_->base(Variable::r, Process::kfin));
0078           const double phi = phiBV.val(dataFormats_->base(Variable::phi, Process::kfin));
0079           const double z = zBV.val(dataFormats_->base(Variable::z, Process::kfin));
0080           stubs_.emplace_back(frameStub.first, r, phi, z, layerIdBV.val(), tiltBV.val(), layer);
0081           stubs.push_back(&stubs_.back());
0082         }
0083         TTBV ttBV = frameTrack.second;
0084         const TTBV cotBV(ttBV, dataFormats_->format(Variable::cot, Process::kfin).width(), 0, true);
0085         ttBV >>= dataFormats_->format(Variable::cot, Process::kfin).width();
0086         const TTBV zTBV(ttBV, dataFormats_->format(Variable::zT, Process::kfin).width(), 0, true);
0087         ttBV >>= dataFormats_->format(Variable::zT, Process::kfin).width();
0088         const TTBV phiTBV(ttBV, dataFormats_->format(Variable::phiT, Process::kfin).width(), 0, true);
0089         ttBV >>= dataFormats_->format(Variable::phiT, Process::kfin).width();
0090         const TTBV inv2RBV(ttBV, dataFormats_->format(Variable::inv2R, Process::kfin).width(), 0, true);
0091         ttBV >>= dataFormats_->format(Variable::inv2R, Process::kfin).width();
0092         const TTBV sectorEtaBV(ttBV, dataFormats_->format(Variable::sectorEta, Process::kfin).width(), 0);
0093         ttBV >>= dataFormats_->format(Variable::sectorEta, Process::kfin).width();
0094         const TTBV sectorPhiBV(ttBV, dataFormats_->format(Variable::sectorPhi, Process::kfin).width(), 0);
0095         const double cot = cotBV.val(dataFormats_->base(Variable::cot, Process::kfin));
0096         const double zT = zTBV.val(dataFormats_->base(Variable::zT, Process::kfin));
0097         const double inv2R = inv2RBV.val(dataFormats_->base(Variable::inv2R, Process::kfin));
0098         const int sectorEta = sectorEtaBV.val();
0099         const int zTu = dataFormats_->format(Variable::zT, Process::kfin).toUnsigned(zT);
0100         const int cotu = dataFormats_->format(Variable::cot, Process::kfin).toUnsigned(cot);
0101         const TTBV maybe = layerEncoding_->maybePattern(sectorEta, zTu, cotu);
0102         const FrameTrack frameT(frameTrack.first,
0103                                 Frame("1" + maybe.str() + sectorPhiBV.str() + sectorEtaBV.str() + phiTBV.str() +
0104                                       inv2RBV.str() + zTBV.str() + cotBV.str()));
0105         tracks_.emplace_back(frameT, stubs, cot, zT, inv2R, sectorEtaBV.val());
0106         input.push_back(&tracks_.back());
0107       }
0108       // remove all gaps between end and last track
0109       for (auto it = input.end(); it != input.begin();)
0110         it = (*--it) ? input.begin() : input.erase(it);
0111     }
0112   }
0113
0114   // fill output products
0115   void KFin::produce(StreamsStub& accpetedStubs,
0116                      StreamsTrack& acceptedTracks,
0117                      StreamsStub& lostStubs,
0118                      StreamsTrack& lostTracks) {
0119     // calculate stub uncertainties
0120     static constexpr int usedMSBpitchOverRaddr = 1;
0121     static const double baseRlut =
0122         dataFormats_->base(Variable::r, Process::kfin) *
0123         pow(2, dataFormats_->width(Variable::r, Process::zht) - setup_->widthAddrBRAM18() + usedMSBpitchOverRaddr);
0124     static const double baseRinvR = dataFormats_->base(Variable::r, Process::kfin) *
0125                                     pow(2, dataFormats_->width(Variable::r, Process::zht) - setup_->widthAddrBRAM18());
0126     static const double basePhi =
0127         dataFormats_->base(Variable::inv2R, Process::kfin) * dataFormats_->base(Variable::r, Process::kfin);
0128     static const double baseInvR =
0129         pow(2.,
0130             ceil(log2(dataFormats_->base(Variable::r, Process::kfin) / setup_->tbInnerRadius())) -
0131                 setup_->widthDSPbu()) /
0132         dataFormats_->base(Variable::r, Process::kfin);
0133     static const double maxCot = sinh(setup_->maxEta()) + setup_->beamWindowZ() / setup_->chosenRofZ();
0134     static constexpr int usedMSBCotLutaddr = 3;
0135     static const double baseCotLut = pow(2., ceil(log2(maxCot)) - setup_->widthAddrBRAM18() + usedMSBCotLutaddr);
0136     static const double baseCot = dataFormats_->base(Variable::cot, Process::kfin);
0137     static const double baseZ = dataFormats_->base(Variable::z, Process::kfin);
0138     static const double baseR = dataFormats_->base(Variable::r, Process::kfin);
0139     for (const Track& track : tracks_) {
0140       const int sectorEta = track.sectorEta_;
0141       const double inv2R = abs(track.inv2R_);
0142       for (Stub* stub : track.stubs_) {
0143         const bool barrel = setup_->barrel(stub->ttStubRef_);
0144         const bool ps = barrel ? setup_->psModule(stub->ttStubRef_) : stub->psTilt_;
0145         const bool tilt = barrel ? (ps && !stub->psTilt_) : false;
0146         const double length = ps ? setup_->lengthPS() : setup_->length2S();
0147         const double pitch = ps ? setup_->pitchPS() : setup_->pitch2S();
0148         const double pitchOverR = digi(pitch / (digi(stub->r_, baseRlut) + dataFormats_->chosenRofPhi()), basePhi);
0149         const double r = digi(stub->r_, baseRinvR) + dataFormats_->chosenRofPhi();
0150         const double sumdz = track.zT_ + stub->z_;
0151         const double dZ = digi(sumdz - digi(setup_->chosenRofZ(), baseR) * track.cot_, baseCot * baseR);
0152         const double sumcot = track.cot_ + digi(setup_->sectorCot(sectorEta), baseCot);
0153         const double cot = digi(abs(dZ * digi(1. / r, baseInvR) + sumcot), baseCotLut);
0154         double lengthZ = length;
0155         double lengthR = 0.;
0156         if (!barrel) {
0157           lengthZ = length * cot;
0158           lengthR = length;
0159         } else if (tilt) {
0160           lengthZ = length * abs(setup_->tiltApproxSlope() * cot + setup_->tiltApproxIntercept());
0161           lengthR = setup_->tiltUncertaintyR();
0162         }
0163         const double scat = digi(setup_->scattering(), baseR);
0164         stub->dZ_ = lengthZ + baseZ;
0165         stub->dPhi_ = (scat + digi(lengthR, baseR)) * inv2R + pitchOverR;
0166         stub->dPhi_ = digi(stub->dPhi_, basePhi) + basePhi;
0167       }
0168     }
0169     // store helper
0170     auto frameTrack = [](Track* track) { return track->frame_; };
0171     auto frameStub = [this](Track* track, int layer) {
0172       auto equal = [layer](Stub* stub) { return stub->channel_ == layer; };
0173       const auto it = find_if(track->stubs_.begin(), track->stubs_.end(), equal);
0174       if (it == track->stubs_.end())
0175         return FrameStub();
0176       Stub* stub = *it;
0177       const TTBV r(dataFormats_->format(Variable::r, Process::kfin).ttBV(stub->r_));
0178       const TTBV phi(dataFormats_->format(Variable::phi, Process::kfin).ttBV(stub->phi_));
0179       const TTBV z(dataFormats_->format(Variable::z, Process::kfin).ttBV(stub->z_));
0180       const TTBV dPhi(dataFormats_->format(Variable::dPhi, Process::kfin).ttBV(stub->dPhi_));
0181       const TTBV dZ(dataFormats_->format(Variable::dZ, Process::kfin).ttBV(stub->dZ_));
0182       return FrameStub(stub->ttStubRef_, Frame("1" + r.str() + phi.str() + z.str() + dPhi.str() + dZ.str()));
0183     };
0184     // merge number of nodes DR to number of Nodes KF and store result
0185     static const int nMux = channelAssignment_->numNodesDR() / setup_->kfNumWorker();
0186     const int offsetTrack = region_ * setup_->kfNumWorker();
0187     for (int nodeKF = 0; nodeKF < setup_->kfNumWorker(); nodeKF++) {
0188       const int offset = nodeKF * nMux;
0189       deque<Track*> accepted;
0190       deque<Track*> lost;
0191       vector<deque<Track*>> stacks(nMux);
0192       vector<deque<Track*>> inputs(nMux);
0193       for (int channel = 0; channel < nMux; channel++) {
0194         const vector<Track*>& input = input_[offset + channel];
0195         inputs[channel] = deque<Track*>(input.begin(), input.end());
0196       }
0197       // clock accurate firmware emulation, each while trip describes one clock tick, one stub in and one stub out per tick
0198       while (!all_of(inputs.begin(), inputs.end(), [](const deque<Track*>& tracks) { return tracks.empty(); }) or
0199              !all_of(stacks.begin(), stacks.end(), [](const deque<Track*>& tracks) { return tracks.empty(); })) {
0200         // fill input fifos
0201         for (int channel = 0; channel < nMux; channel++) {
0202           deque<Track*>& stack = stacks[channel];
0203           Track* track = pop_front(inputs[channel]);
0204           if (track)
0205             stack.push_back(track);
0206         }
0207         // merge input fifos to one stream, prioritizing higher input channel over lower channel
0208         bool nothingToRoute(true);
0209         for (int channel = nMux - 1; channel >= 0; channel--) {
0210           Track* track = pop_front(stacks[channel]);
0211           if (track) {
0212             nothingToRoute = false;
0213             accepted.push_back(track);
0214             break;
0215           }
0216         }
0217         if (nothingToRoute)
0218           accepted.push_back(nullptr);
0219       }
0220       // truncate if desired
0221       if (enableTruncation_ && (int)accepted.size() > setup_->numFrames()) {
0222         const auto limit = next(accepted.begin(), setup_->numFrames());
0223         copy_if(limit, accepted.end(), back_inserter(lost), [](const Track* track) { return track; });
0224         accepted.erase(limit, accepted.end());
0225       }
0226       // remove all gaps between end and last track
0227       for (auto it = accepted.end(); it != accepted.begin();)
0228         it = (*--it) ? accepted.begin() : accepted.erase(it);
0229       // fill products StreamsStub& accpetedStubs, StreamsTrack& acceptedTracks, StreamsStub& lostStubs, StreamsTrack& lostTracks
0230       const int channelTrack = offsetTrack + nodeKF;
0231       const int offsetStub = channelTrack * setup_->numLayers();
0232       // fill lost tracks and stubs without gaps
0233       lostTracks[channelTrack].reserve(lost.size());
0234       for (int layer = 0; layer < setup_->numLayers(); layer++)
0235         lostStubs[offsetStub + layer].reserve(lost.size());
0236       for (Track* track : lost) {
0237         lostTracks[channelTrack].emplace_back(frameTrack(track));
0238         for (int layer = 0; layer < setup_->numLayers(); layer++)
0239           lostStubs[offsetStub + layer].emplace_back(frameStub(track, layer));
0240       }
0241       // fill accepted tracks and stubs with gaps
0242       acceptedTracks[channelTrack].reserve(accepted.size());
0243       for (int layer = 0; layer < setup_->numLayers(); layer++)
0244         accpetedStubs[offsetStub + layer].reserve(accepted.size());
0245       for (Track* track : accepted) {
0246         if (!track) {  // fill gap
0247           acceptedTracks[channelTrack].emplace_back(FrameTrack());
0248           for (int layer = 0; layer < setup_->numLayers(); layer++)
0249             accpetedStubs[offsetStub + layer].emplace_back(FrameStub());
0250           continue;
0251         }
0252         acceptedTracks[channelTrack].emplace_back(frameTrack(track));
0253         for (int layer = 0; layer < setup_->numLayers(); layer++)
0254           accpetedStubs[offsetStub + layer].emplace_back(frameStub(track, layer));
0255       }
0256     }
0257   }
0258
0259   // remove and return first element of deque, returns nullptr if empty
0260   template <class T>
0261   T* KFin::pop_front(deque<T*>& ts) const {
0262     T* t = nullptr;
0263     if (!ts.empty()) {
0264       t = ts.front();
0265       ts.pop_front();
0266     }
0267     return t;
0268   }
0269
0270 }  // namespace trklet