Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-07 22:33:25

0001 #include "L1Trigger/TrackerDTC/interface/Stub.h"
0002 
0003 #include <cmath>
0004 #include <iterator>
0005 #include <algorithm>
0006 #include <utility>
0007 
0008 using namespace edm;
0009 using namespace std;
0010 
0011 namespace trackerDTC {
0012 
0013   Stub::Stub(const ParameterSet& iConfig, const Setup& setup, SensorModule* sm, const TTStubRef& ttStubRef)
0014       : setup_(&setup), sm_(sm), ttStubRef_(ttStubRef), hybrid_(iConfig.getParameter<bool>("UseHybrid")), valid_(true) {
0015     regions_.reserve(setup.numOverlappingRegions());
0016     // get stub local coordinates
0017     const MeasurementPoint& mp = ttStubRef->clusterRef(0)->findAverageLocalCoordinatesCentered();
0018 
0019     // convert to uniformed local coordinates
0020 
0021     // column number in pitch units
0022     col_ = (int)floor(pow(-1, sm->signCol()) * (mp.y() - sm->numColumns() / 2) / setup.baseCol());
0023     // row number in half pitch units
0024     row_ = (int)floor(pow(-1, sm->signRow()) * (mp.x() - sm->numRows() / 2) / setup.baseRow());
0025     // bend number in quarter pitch units
0026     bend_ = (int)floor(pow(-1, sm->signBend()) * (ttStubRef->bendBE()) / setup.baseBend());
0027     // reduced row number for look up
0028     rowLUT_ = (int)floor((double)row_ / pow(2., setup.widthRow() - setup.dtcWidthRowLUT()));
0029     // sub row number inside reduced row number
0030     rowSub_ = row_ - (rowLUT_ + .5) * pow(2, setup.widthRow() - setup.dtcWidthRowLUT());
0031 
0032     // convert local to global coordinates
0033 
0034     const double y = (col_ + .5) * setup.baseCol() * sm->pitchCol();
0035     // radius of a column of strips/pixel in cm
0036     d_ = sm->r() + y * sm->sin();
0037     // stub z in cm
0038     z_ = digi(sm->z() + y * sm->cos(), setup.baseZ());
0039 
0040     const double x0 = rowLUT_ * setup.baseRow() * setup.dtcNumMergedRows() * sm->pitchRow();
0041     const double x1 = (rowLUT_ + 1) * setup.baseRow() * setup.dtcNumMergedRows() * sm->pitchRow();
0042     const double x = (rowLUT_ + .5) * setup.baseRow() * setup.dtcNumMergedRows() * sm->pitchRow();
0043     // stub r in cm
0044     r_ = sqrt(d_ * d_ + x * x);
0045 
0046     const double phi0 = sm->phi() + atan2(x0, d_);
0047     const double phi1 = sm->phi() + atan2(x1, d_);
0048     const double c = (phi0 + phi1) / 2.;
0049     const double m = (phi1 - phi0) / setup.dtcNumMergedRows();
0050 
0051     // intercept of linearized stub phi in rad
0052     c_ = digi(c, setup.basePhi());
0053     // slope of linearized stub phi in rad / strip
0054     m_ = digi(m, setup.dtcBaseM());
0055 
0056     if (hybrid_) {
0057       if (abs(z_ / r_) > setup.hybridMaxCot())
0058         // did not pass eta cut
0059         valid_ = false;
0060     } else {
0061       // extrapolated z at radius T assuming z0=0
0062       const double zT = setup.chosenRofZ() * z_ / r_;
0063       // extrapolated z0 window at radius T
0064       const double dZT = setup.beamWindowZ() * abs(1. - setup.chosenRofZ() / r_);
0065       double zTMin = zT - dZT;
0066       double zTMax = zT + dZT;
0067       if (zTMin >= setup.maxZT() || zTMax < -setup.maxZT())
0068         // did not pass "eta" cut
0069         valid_ = false;
0070       else {
0071         zTMin = max(zTMin, -setup.maxZT());
0072         zTMax = min(zTMax, setup.maxZT());
0073       }
0074       // range of stub cot(theta)
0075       cot_ = {zTMin / setup.chosenRofZ(), zTMax / setup.chosenRofZ()};
0076     }
0077 
0078     // stub r w.r.t. chosenRofPhi in cm
0079     static const double chosenRofPhi = hybrid_ ? setup.hybridChosenRofPhi() : setup.chosenRofPhi();
0080     r_ = digi(r_ - chosenRofPhi, setup.baseR());
0081 
0082     // radial (cylindrical) component of sensor separation
0083     const double dr = sm->sep() / (sm->cos() - sm->sin() * z_ / d_);
0084     // converts bend into qOverPt in 1/cm
0085     const double qOverPtOverBend = sm->pitchRow() / dr / d_;
0086     // qOverPt in 1/cm
0087     const double qOverPt = bend_ * setup.baseBend() * qOverPtOverBend;
0088     // qOverPt uncertainty in 1/cm
0089     const double dQoverPt = setup.bendCut() * qOverPtOverBend;
0090     const double minPt = hybrid_ ? setup.hybridMinPt() : setup.minPt();
0091     const double maxQoverPt = setup.invPtToDphi() / minPt - setup.dtcBaseQoverPt() / 2.;
0092     double qOverPtMin = digi(qOverPt - dQoverPt, setup.dtcBaseQoverPt());
0093     double qOverPtMax = digi(qOverPt + dQoverPt, setup.dtcBaseQoverPt());
0094     if (qOverPtMin >= maxQoverPt || qOverPtMax < -maxQoverPt)
0095       // did not pass pt cut
0096       valid_ = false;
0097     else {
0098       qOverPtMin = max(qOverPtMin, -maxQoverPt);
0099       qOverPtMax = min(qOverPtMax, maxQoverPt);
0100     }
0101     // range of stub qOverPt in 1/cm
0102     qOverPt_ = {qOverPtMin, qOverPtMax};
0103 
0104     // stub phi w.r.t. detector region centre in rad
0105     phi_ = c_ + rowSub_ * m_;
0106 
0107     // range of stub extrapolated phi to radius chosenRofPhi in rad
0108     phiT_.first = phi_ + r_ * qOverPt_.first;
0109     phiT_.second = phi_ + r_ * qOverPt_.second;
0110     if (phiT_.first > phiT_.second)
0111       swap(phiT_.first, phiT_.second);
0112 
0113     if (phiT_.first < 0.)
0114       regions_.push_back(0);
0115     if (phiT_.second >= 0.)
0116       regions_.push_back(1);
0117 
0118     // apply data format specific manipulations
0119     if (!hybrid_)
0120       return;
0121 
0122     // stub r w.r.t. an offset in cm
0123     r_ -= sm->offsetR() - chosenRofPhi;
0124     // stub z w.r.t. an offset in cm
0125     z_ -= sm->offsetZ();
0126     if (sm->type() == SensorModule::Disk2S) {
0127       // encoded r
0128       r_ = sm->encodedR() + (sm->side() ? -col_ : (col_ + sm->numColumns() / 2));
0129       r_ = (r_ + 0.5) * setup.hybridBaseR(sm->type());
0130     }
0131 
0132     // encode bend
0133     const vector<double>& encodingBend = setup.encodingBend(sm->windowSize(), sm->psModule());
0134     const auto pos = find(encodingBend.begin(), encodingBend.end(), abs(ttStubRef->bendBE()));
0135     const int uBend = distance(encodingBend.begin(), pos);
0136     bend_ = pow(-1, signbit(bend_)) * uBend;
0137   }
0138 
0139   // returns bit accurate representation of Stub
0140   TTDTC::BV Stub::frame(int region) const { return hybrid_ ? formatHybrid(region) : formatTMTT(region); }
0141 
0142   // returns true if stub belongs to region
0143   bool Stub::inRegion(int region) const { return find(regions_.begin(), regions_.end(), region) != regions_.end(); }
0144 
0145   // truncates double precision to f/w integer equivalent
0146   double Stub::digi(double value, double precision) const { return (floor(value / precision) + .5) * precision; }
0147 
0148   // returns 64 bit stub in hybrid data format
0149   TTDTC::BV Stub::formatHybrid(int region) const {
0150     const SensorModule::Type type = sm_->type();
0151     // stub phi w.r.t. processing region centre in rad
0152     const double phi = phi_ - (region - .5) * setup_->baseRegion() +
0153                        0.5 * setup_->hybridBasePhi(type) * (1 << setup_->hybridWidthPhi(type));
0154 
0155     // convert stub variables into bit vectors
0156     const TTBV hwR(r_, setup_->hybridBaseR(type), setup_->hybridWidthR(type), true);
0157     const TTBV hwPhi(phi, setup_->hybridBasePhi(type), setup_->hybridWidthPhi(type), true);
0158     const TTBV hwZ(z_, setup_->hybridBaseZ(type), setup_->hybridWidthZ(type), true);
0159     const TTBV hwAlpha(row_, setup_->hybridBaseAlpha(type), setup_->hybridWidthAlpha(type), true);
0160     const TTBV hwBend(bend_, setup_->hybridWidthBend(type), true);
0161     const TTBV hwLayer(sm_->encodedLayerId(), setup_->hybridWidthLayer());
0162     const TTBV hwGap(0, setup_->hybridNumUnusedBits(type));
0163     const TTBV hwValid(1, 1);
0164     // assemble final bitset
0165     return TTDTC::BV(hwGap.str() + hwR.str() + hwZ.str() + hwPhi.str() + hwAlpha.str() + hwBend.str() + hwLayer.str() +
0166                      hwValid.str());
0167   }
0168 
0169   TTDTC::BV Stub::formatTMTT(int region) const {
0170     int layerM = sm_->layerId();
0171     // convert unique layer id [1-6,11-15] into reduced layer id [0-6]
0172     // a fiducial track may not cross more then 7 detector layers, for stubs from a given track the reduced layer id is actually unique
0173     int layer(-1);
0174     if (layerM == 1)
0175       layer = 0;
0176     else if (layerM == 2)
0177       layer = 1;
0178     else if (layerM == 6 || layerM == 11)
0179       layer = 2;
0180     else if (layerM == 5 || layerM == 12)
0181       layer = 3;
0182     else if (layerM == 4 || layerM == 13)
0183       layer = 4;
0184     else if (layerM == 14)
0185       layer = 5;
0186     else if (layerM == 3 || layerM == 15)
0187       layer = 6;
0188     // assign stub to phi sectors within a processing region, to be generalized
0189     TTBV sectorsPhi(0, setup_->numOverlappingRegions() * setup_->numSectorsPhi());
0190     if (phiT_.first < 0.) {
0191       if (phiT_.first < -setup_->baseSector())
0192         sectorsPhi.set(0);
0193       else
0194         sectorsPhi.set(1);
0195       if (phiT_.second < 0. && phiT_.second >= -setup_->baseSector())
0196         sectorsPhi.set(1);
0197     }
0198     if (phiT_.second >= 0.) {
0199       if (phiT_.second < setup_->baseSector())
0200         sectorsPhi.set(2);
0201       else
0202         sectorsPhi.set(3);
0203       if (phiT_.first >= 0. && phiT_.first < setup_->baseSector())
0204         sectorsPhi.set(2);
0205     }
0206     // assign stub to eta sectors within a processing region
0207     pair<int, int> setcorEta({0, setup_->numSectorsEta() - 1});
0208     for (int bin = 0; bin < setup_->numSectorsEta(); bin++)
0209       if (asinh(cot_.first) < setup_->boundarieEta(bin + 1)) {
0210         setcorEta.first = bin;
0211         break;
0212       }
0213     for (int bin = setcorEta.first; bin < setup_->numSectorsEta(); bin++)
0214       if (asinh(cot_.second) < setup_->boundarieEta(bin + 1)) {
0215         setcorEta.second = bin;
0216         break;
0217       }
0218     // stub phi w.r.t. processing region centre in rad
0219     const double phi = phi_ - (region - .5) * setup_->baseRegion();
0220     // convert stub variables into bit vectors
0221     const TTBV hwValid(1, 1);
0222     const TTBV hwGap(0, setup_->dtcNumUnusedBits());
0223     const TTBV hwLayer(layer, setup_->widthLayer());
0224     const TTBV hwSectorEtaMin(setcorEta.first, setup_->widthSectorEta());
0225     const TTBV hwSectorEtaMax(setcorEta.second, setup_->widthSectorEta());
0226     const TTBV hwR(r_, setup_->baseR(), setup_->widthR(), true);
0227     const TTBV hwPhi(phi, setup_->basePhi(), setup_->widthPhiDTC(), true);
0228     const TTBV hwZ(z_, setup_->baseZ(), setup_->widthZ(), true);
0229     const TTBV hwQoverPtMin(qOverPt_.first, setup_->htBaseQoverPt(), setup_->htWidthQoverPt(), true);
0230     const TTBV hwQoverPtMax(qOverPt_.second, setup_->htBaseQoverPt(), setup_->htWidthQoverPt(), true);
0231     TTBV hwSectorPhis(0, setup_->numSectorsPhi());
0232     for (int sectorPhi = 0; sectorPhi < setup_->numSectorsPhi(); sectorPhi++)
0233       hwSectorPhis[sectorPhi] = sectorsPhi[region * setup_->numSectorsPhi() + sectorPhi];
0234     // assemble final bitsetTTDTC::BV(hwGap.str() + hwValid.str() + hwR.str() + hwPhi.str() + hwZ.str() + hwQoverPtMin.str() +
0235     return TTDTC::BV(hwGap.str() + hwValid.str() + hwR.str() + hwPhi.str() + hwZ.str() + hwQoverPtMin.str() +
0236                      hwQoverPtMax.str() + hwSectorEtaMin.str() + hwSectorEtaMax.str() + hwSectorPhis.str() +
0237                      hwLayer.str());
0238   }
0239 
0240 }  // namespace trackerDTC