TTStubBuilder.cc

Go to the documentation of this file.

#include "L1Trigger/TrackTrigger/plugins/TTStubBuilder.h"


template <>
void TTStubBuilder<Ref_Phase2TrackerDigi_>::updateStubs(
    const edm::OrphanHandle<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>>& clusterHandle,
    const edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>& inputEDstubs,
    edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>& outputEDstubs) const {
  for (const auto& module : inputEDstubs) {
    DetId thisStackedDetId = module.id();
    typename edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>::FastFiller outputFiller(outputEDstubs,
                                                                                           thisStackedDetId);

    DetId lowerDetid = thisStackedDetId + 1;
    DetId upperDetid = thisStackedDetId + 2;

    edmNew::DetSet<TTCluster<Ref_Phase2TrackerDigi_>> lowerClusters = (*clusterHandle)[lowerDetid];
    edmNew::DetSet<TTCluster<Ref_Phase2TrackerDigi_>> upperClusters = (*clusterHandle)[upperDetid];

    edmNew::DetSet<TTStub<Ref_Phase2TrackerDigi_>> theseStubs = inputEDstubs[thisStackedDetId];

    typename edmNew::DetSet<TTCluster<Ref_Phase2TrackerDigi_>>::const_iterator clusterIter;
    for (const auto& stub : theseStubs) {
      TTStub<Ref_Phase2TrackerDigi_> tempTTStub(stub.getDetId());

      const edm::Ref<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>, TTCluster<Ref_Phase2TrackerDigi_>>&
          lowerClusterToBeReplaced = stub.clusterRef(0);
      const edm::Ref<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>, TTCluster<Ref_Phase2TrackerDigi_>>&
          upperClusterToBeReplaced = stub.clusterRef(1);

      bool lowerOK = false;
      bool upperOK = false;

      for (clusterIter = lowerClusters.begin(); clusterIter != lowerClusters.end() && !lowerOK; ++clusterIter) {
        if (clusterIter->getHits() == lowerClusterToBeReplaced->getHits()) {
          tempTTStub.addClusterRef(edmNew::makeRefTo(clusterHandle, clusterIter));
          lowerOK = true;
        }
      }

      for (clusterIter = upperClusters.begin(); clusterIter != upperClusters.end() && !upperOK; ++clusterIter) {
        if (clusterIter->getHits() == upperClusterToBeReplaced->getHits()) {
          tempTTStub.addClusterRef(edmNew::makeRefTo(clusterHandle, clusterIter));
          upperOK = true;
        }
      }

      if (!lowerOK || !upperOK)
        continue;

      tempTTStub.setRawBend(2. * stub.rawBend());
      tempTTStub.setBendOffset(2. * stub.bendOffset());
      tempTTStub.setBendBE(stub.bendBE());
      tempTTStub.setModuleTypePS(stub.moduleTypePS());

      outputFiller.push_back(tempTTStub);

    }  
  }    
}


template <>
void TTStubBuilder<Ref_Phase2TrackerDigi_>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
  //Retrieve tracker topology from geometry
  edm::ESHandle<TrackerTopology> tTopoHandle = iSetup.getHandle(tTopoToken);
  const TrackerTopology* const tTopo = tTopoHandle.product();
  edm::ESHandle<TrackerGeometry> tGeomHandle = iSetup.getHandle(tGeomToken);
  const TrackerGeometry* const theTrackerGeom = tGeomHandle.product();

  auto ttClusterDSVForOutputAcc = std::make_unique<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>>();
  auto ttClusterDSVForOutputRej = std::make_unique<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>>();
  auto ttStubDSVForOutputAccTemp = std::make_unique<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>>();
  auto ttStubDSVForOutputRejTemp = std::make_unique<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>>();
  auto ttStubDSVForOutputAcc = std::make_unique<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>>();
  auto ttStubDSVForOutputRej = std::make_unique<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>>();

  edm::Handle<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>> clusterHandle;
  iEvent.getByToken(clustersToken, clusterHandle);

  int nmod = -1;

  // Loop over all the tracker elements

  for (const auto& gd : theTrackerGeom->dets()) {
    DetId detid = (*gd).geographicalId();
    if (detid.subdetId() != StripSubdetector::TOB && detid.subdetId() != StripSubdetector::TID)
      continue;  // only run on OT
    if (!tTopo->isLower(detid))
      continue;  // loop on the stacks: choose the lower sensor
    DetId lowerDetid = detid;
    DetId upperDetid = tTopo->partnerDetId(detid);
    DetId stackDetid = tTopo->stack(detid);
    bool isPS = (theTrackerGeom->getDetectorType(stackDetid) == TrackerGeometry::ModuleType::Ph2PSP);

    bool is10G_PS = false;

    // Determine if this module is a 10G transmission scheme module
    //
    // TO FIX: take this info from Tracker -> DTC cabling map.

    if (detid.subdetId() == StripSubdetector::TOB) {
      if (tTopo->layer(detid) <= high_rate_max_layer)
        is10G_PS = true;
    } else if (detid.subdetId() == StripSubdetector::TID) {
      if (tTopo->tidRing(detid) <= high_rate_max_ring[tTopo->tidWheel(detid) - 1])
        is10G_PS = true;
    }

    ++nmod;

    unsigned int maxStubs;
    std::vector<std::pair<unsigned int, double>> bendMap;

    if (clusterHandle->find(lowerDetid) == clusterHandle->end() ||
        clusterHandle->find(upperDetid) == clusterHandle->end())
      continue;

    edmNew::DetSet<TTCluster<Ref_Phase2TrackerDigi_>> lowerClusters = (*clusterHandle)[lowerDetid];
    edmNew::DetSet<TTCluster<Ref_Phase2TrackerDigi_>> upperClusters = (*clusterHandle)[upperDetid];

    if (lowerClusters.empty() || upperClusters.empty())
      continue;

    std::vector<TTCluster<Ref_Phase2TrackerDigi_>> tempClusLowerAcc;
    std::vector<TTCluster<Ref_Phase2TrackerDigi_>> tempClusLowerRej;
    std::vector<TTCluster<Ref_Phase2TrackerDigi_>> tempClusUpperAcc;
    std::vector<TTCluster<Ref_Phase2TrackerDigi_>> tempClusUpperRej;
    std::vector<TTStub<Ref_Phase2TrackerDigi_>> tempStubAcc;
    std::vector<TTStub<Ref_Phase2TrackerDigi_>> tempStubRej;
    tempClusLowerAcc.clear();
    tempClusLowerRej.clear();
    tempClusUpperAcc.clear();
    tempClusUpperRej.clear();
    tempStubAcc.clear();
    tempStubRej.clear();

    const int chipSize = isPS ? 120 : 127;
    // No. of macro pixels along local y in each half of 2S module.
    constexpr int numMacroPixels = 16;

    for (auto lowerClusterIter = lowerClusters.begin(); lowerClusterIter != lowerClusters.end(); ++lowerClusterIter) {
      std::vector<TTStub<Ref_Phase2TrackerDigi_>> tempOutput;

      for (auto upperClusterIter = upperClusters.begin(); upperClusterIter != upperClusters.end(); ++upperClusterIter) {
        TTStub<Ref_Phase2TrackerDigi_> tempTTStub(stackDetid);
        tempTTStub.addClusterRef(edmNew::makeRefTo(clusterHandle, lowerClusterIter));
        tempTTStub.addClusterRef(edmNew::makeRefTo(clusterHandle, upperClusterIter));
        tempTTStub.setModuleTypePS(isPS);

        bool thisConfirmation = false;
        int thisDisplacement = 999999;
        int thisOffset = 0;
        float thisHardBend = 0;

        theStubFindingAlgoHandle->PatternHitCorrelation(
            thisConfirmation, thisDisplacement, thisOffset, thisHardBend, tempTTStub);
        // Removed real offset.  Ivan Reid 10/2019

        if (thisConfirmation) {
          tempTTStub.setRawBend(thisDisplacement);
          tempTTStub.setBendOffset(thisOffset);
          tempTTStub.setBendBE(thisHardBend);
          tempOutput.push_back(tempTTStub);
        }  
      }    

      if (ForbidMultipleStubs && tempOutput.size() > 1 + static_cast<unsigned int>(isPS)) {
        std::sort(tempOutput.begin(), tempOutput.end(), TTStubBuilder<Ref_Phase2TrackerDigi_>::SortStubsBend);

        typename std::vector<TTStub<Ref_Phase2TrackerDigi_>>::iterator tempIter = tempOutput.begin();
        ++tempIter;
        if (isPS)
          ++tempIter;  // PS module case

        tempOutput.erase(tempIter, tempOutput.end());
      }



      for (auto& tempTTStub : tempOutput) {
        if (not applyFE)  // No dynamic inefficiencies
        {
          tempClusLowerAcc.push_back(*(tempTTStub.clusterRef(0)));
          tempClusUpperAcc.push_back(*(tempTTStub.clusterRef(1)));
          tempStubAcc.push_back(tempTTStub);
        } else {
          bool FEreject = false;

          MeasurementPoint mp0 = tempTTStub.clusterRef(0)->findAverageLocalCoordinates();
          int seg = static_cast<int>(mp0.y());  // Identifies which half of module
          if (isPS)
            seg = seg / numMacroPixels;
          int chip = 1000 * nmod + 10 * int(tempTTStub.innerClusterPosition() / chipSize) + seg;
          int CIC_chip = 10 * nmod + seg;

          // First look is the stub is passing trough the very front end (CBC/MPA)
          maxStubs = isPS ? maxStubs_PS : maxStubs_2S;

          if (isPS)  // MPA
          {
            if (moduleStubs_MPA[chip] < int(maxStubs)) {
              moduleStubs_MPA[chip]++;
            } else {
              FEreject = true;
            }
          } else  // CBC
          {
            if (moduleStubs_CBC[chip] < int(maxStubs)) {
              moduleStubs_CBC[chip]++;
            } else {
              FEreject = true;
            }
          }

          // End of the MPA/CBC loop

          // If the stub has been already thrown out, there is no reason to include it into the CIC stream
          // We put it in the rejected container, flagged with offset to indicate reason.

          if (FEreject) {
            tempTTStub.setRawBend(CBCFailOffset + 2. * tempTTStub.rawBend());
            tempTTStub.setBendOffset(CBCFailOffset + 2. * tempTTStub.bendOffset());
            tempClusLowerRej.push_back(*(tempTTStub.clusterRef(0)));
            tempClusUpperRej.push_back(*(tempTTStub.clusterRef(1)));
            tempStubRej.push_back(tempTTStub);
            continue;
          }

          maxStubs = isPS ? maxStubs_PS_CIC_5 : maxStubs_2S_CIC_5;

          if (is10G_PS)
            maxStubs = maxStubs_PS_CIC_10;

          bool CIC_reject = true;

          moduleStubs_CIC[CIC_chip].push_back(tempTTStub);  //We temporarily add the new stub

          if (moduleStubs_CIC[CIC_chip].size() <= maxStubs) {
            tempClusLowerAcc.push_back(*(tempTTStub.clusterRef(0)));
            tempClusUpperAcc.push_back(*(tempTTStub.clusterRef(1)));
            tempStubAcc.push_back(tempTTStub);  // The stub is kept

          } else {

            bendMap.clear();
            bendMap.reserve(moduleStubs_CIC[CIC_chip].size());

            for (unsigned int i = 0; i < moduleStubs_CIC[CIC_chip].size(); ++i) {
              bendMap.emplace_back(i, moduleStubs_CIC[CIC_chip].at(i).bendFE());
            }

            std::sort(bendMap.begin(), bendMap.end(), TTStubBuilder<Ref_Phase2TrackerDigi_>::SortStubBendPairs);

            // bendMap contains link over all the stubs included in moduleStubs_CIC[CIC_chip]

            for (unsigned int i = 0; i < maxStubs; ++i) {
              // The stub we have just added is amongst those with smallest |bend| so keep it.
              if (bendMap[i].first == moduleStubs_CIC[CIC_chip].size() - 1) {
                CIC_reject = false;
              }
            }

            if (CIC_reject)  // The stub added does not pass the cut
            {
              tempTTStub.setRawBend(CICFailOffset + 2. * tempTTStub.rawBend());
              tempTTStub.setBendOffset(CICFailOffset + 2. * tempTTStub.bendOffset());
              tempClusLowerRej.push_back(*(tempTTStub.clusterRef(0)));
              tempClusUpperRej.push_back(*(tempTTStub.clusterRef(1)));
              tempStubRej.push_back(tempTTStub);
            } else {
              tempClusLowerAcc.push_back(*(tempTTStub.clusterRef(0)));
              tempClusUpperAcc.push_back(*(tempTTStub.clusterRef(1)));
              tempStubAcc.push_back(tempTTStub);  // The stub is added
            }
          }
        }  
      }    
    }      

    if (not tempClusLowerAcc.empty())
      this->fill(*ttClusterDSVForOutputAcc, lowerDetid, tempClusLowerAcc);
    if (not tempClusLowerRej.empty())
      this->fill(*ttClusterDSVForOutputRej, lowerDetid, tempClusLowerRej);
    if (not tempClusUpperAcc.empty())
      this->fill(*ttClusterDSVForOutputAcc, upperDetid, tempClusUpperAcc);
    if (not tempClusUpperRej.empty())
      this->fill(*ttClusterDSVForOutputRej, upperDetid, tempClusUpperRej);
    if (not tempStubAcc.empty())
      this->fill(*ttStubDSVForOutputAccTemp, stackDetid, tempStubAcc);
    if (not tempStubRej.empty())
      this->fill(*ttStubDSVForOutputRejTemp, stackDetid, tempStubRej);

  }  

  // Store the subset of clusters associated to stubs.

  edm::OrphanHandle<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>> ttClusterAccHandle =
      iEvent.put(std::move(ttClusterDSVForOutputAcc), "ClusterAccepted");
  edm::OrphanHandle<edmNew::DetSetVector<TTCluster<Ref_Phase2TrackerDigi_>>> ttClusterRejHandle =
      iEvent.put(std::move(ttClusterDSVForOutputRej), "ClusterRejected");


  this->updateStubs(ttClusterAccHandle, *ttStubDSVForOutputAccTemp, *ttStubDSVForOutputAcc);
  this->updateStubs(ttClusterRejHandle, *ttStubDSVForOutputRejTemp, *ttStubDSVForOutputRej);

  iEvent.put(std::move(ttStubDSVForOutputAcc), "StubAccepted");
  iEvent.put(std::move(ttStubDSVForOutputRej), "StubRejected");

  ++ievt;
  if (ievt % 8 == 0)
    moduleStubs_CIC.clear();  // Everything is cleared up after 8BX
  if (ievt % 2 == 0)
    moduleStubs_MPA.clear();  // Everything is cleared up after 2BX
  moduleStubs_CBC.clear();    // Everything is cleared up after everyBX
}