Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-06-03 00:12:17

0001 /*
0002 Track Quality Body file
0003 C.Brown & C.Savard 07/2020
0004 */
0005 
0006 #include "L1Trigger/TrackerTFP/interface/TrackQuality.h"
0007 #include "L1Trigger/TrackTrigger/interface/StubPtConsistency.h"
0008 
0009 #include <vector>
0010 #include <map>
0011 #include <string>
0012 #include "conifer.h"
0013 #include "ap_fixed.h"
0014 
0015 namespace trackerTFP {
0016 
0017   TrackQuality::TrackQuality(const ConfigTQ& iConfig, const DataFormats* dataFormats)
0018       : dataFormats_(dataFormats),
0019         model_(iConfig.model_),
0020         featureNames_(iConfig.featureNames_),
0021         baseShiftCot_(iConfig.baseShiftCot_),
0022         baseShiftZ0_(iConfig.baseShiftZ0_),
0023         baseShiftAPfixed_(iConfig.baseShiftAPfixed_),
0024         chi2rphiConv_(iConfig.chi2rphiConv_),
0025         chi2rzConv_(iConfig.chi2rzConv_),
0026         weightBinFraction_(iConfig.weightBinFraction_),
0027         dzTruncation_(iConfig.dzTruncation_),
0028         dphiTruncation_(iConfig.dphiTruncation_) {
0029     dataFormatsTQ_.reserve(+VariableTQ::end);
0030     fillDataFormats(iConfig);
0031   }
0032 
0033   // constructs TQ data formats
0034   template <VariableTQ v>
0035   void TrackQuality::fillDataFormats(const ConfigTQ& iConfig) {
0036     dataFormatsTQ_.emplace_back(makeDataFormat<v>(dataFormats_, iConfig));
0037     if constexpr (v + 1 != VariableTQ::end)
0038       fillDataFormats<v + 1>(iConfig);
0039   }
0040 
0041   // TQ MVA bin conversion LUT
0042   constexpr std::array<double, numBinsMVA_> TrackQuality::mvaPreSigBins() const {
0043     std::array<double, numBinsMVA_> lut = {};
0044     lut[0] = -16.;
0045     for (int i = 1; i < numBinsMVA_; i++)
0046       lut[i] = invSigmoid(TTTrack_TrackWord::tqMVABins[i]);
0047     return lut;
0048   }
0049 
0050   //
0051   template <class T>
0052   int TrackQuality::toBin(const T& bins, double d) const {
0053     int bin = 0;
0054     for (; bin < static_cast<int>(bins.size()) - 1; bin++)
0055       if (d < bins[bin + 1])
0056         break;
0057     return bin;
0058   }
0059 
0060   // Helper function to convert mvaPreSig to bin
0061   int TrackQuality::toBinMVA(double mva) const {
0062     static const std::array<double, numBinsMVA_> bins = mvaPreSigBins();
0063     return toBin(bins, mva);
0064   }
0065 
0066   // Helper function to convert chi2B to bin
0067   int TrackQuality::toBinChi2B(double chi2B) const {
0068     static const std::array<double, numBinsChi2B_> bins = TTTrack_TrackWord::bendChi2Bins;
0069     return toBin(bins, chi2B);
0070   }
0071 
0072   // Helper function to convert chi2rphi to bin
0073   int TrackQuality::toBinchi2rphi(double chi2rphi) const {
0074     static const std::array<double, numBinschi2rphi_> bins = TTTrack_TrackWord::chi2RPhiBins;
0075     double chi2 = chi2rphi * chi2rphiConv_;
0076     return toBin(bins, chi2);
0077   }
0078 
0079   // Helper function to convert chi2rz to bin
0080   int TrackQuality::toBinchi2rz(double chi2rz) const {
0081     static const std::array<double, numBinschi2rz_> bins = TTTrack_TrackWord::chi2RZBins;
0082     double chi2 = chi2rz * chi2rzConv_;
0083     return toBin(bins, chi2);
0084   }
0085 
0086   TrackQuality::Track::Track(const tt::FrameTrack& frameTrack, const tt::StreamStub& streamStub, const TrackQuality* tq)
0087       : frameTrack_(frameTrack), streamStub_(streamStub) {
0088     static const DataFormats* df = tq->dataFormats();
0089     static const tt::Setup* setup = df->setup();
0090     const TrackDR track(frameTrack, df);
0091     double trackchi2rphi(0.);
0092     double trackchi2rz(0.);
0093     TTBV hitPattern(0, streamStub.size());
0094     std::vector<TTStubRef> ttStubRefs;
0095     ttStubRefs.reserve(setup->numLayers());
0096     for (int layer = 0; layer < (int)streamStub.size(); layer++) {
0097       const tt::FrameStub& frameStub = streamStub[layer];
0098       if (frameStub.first.isNull())
0099         continue;
0100       const StubKF stub(frameStub, df);
0101       hitPattern.set(layer);
0102       ttStubRefs.push_back(frameStub.first);
0103       const double m20 = tq->format(VariableTQ::m20).digi(std::pow(stub.phi(), 2));
0104       const double m21 = tq->format(VariableTQ::m21).digi(std::pow(stub.z(), 2));
0105       const double invV0 = tq->format(VariableTQ::invV0).digi(1. / std::pow(stub.dPhi(), 2));
0106       const double invV1 = tq->format(VariableTQ::invV1).digi(1. / std::pow(stub.dZ(), 2));
0107       const double stubchi2rphi = tq->format(VariableTQ::chi2rphi).digi(m20 * invV0);
0108       const double stubchi2rz = tq->format(VariableTQ::chi2rz).digi(m21 * invV1);
0109       trackchi2rphi += stubchi2rphi;
0110       trackchi2rz += stubchi2rz;
0111     }
0112     if (trackchi2rphi > tq->range(VariableTQ::chi2rphi))
0113       trackchi2rphi = tq->range(VariableTQ::chi2rphi) - tq->base(VariableTQ::chi2rphi) / 2.;
0114     if (trackchi2rz > tq->range(VariableTQ::chi2rz))
0115       trackchi2rz = tq->range(VariableTQ::chi2rz) - tq->base(VariableTQ::chi2rz) / 2.;
0116     // calc bdt inputs
0117     const double cot = tq->scaleCot(df->format(Variable::cot, Process::dr).integer(track.cot()));
0118     const double z0 =
0119         tq->scaleZ0(df->format(Variable::zT, Process::kf).integer(track.zT() - setup->chosenRofZ() * track.cot()));
0120     const int nstub = hitPattern.count();
0121     const int n_missint = hitPattern.count(hitPattern.plEncode() + 1, setup->numLayers(), false);
0122     // use simulation for bendchi2
0123     const TTTrackRef& ttTrackRef = frameTrack.first;
0124     const int region = ttTrackRef->phiSector();
0125     const double aRinv = -.5 * track.inv2R();
0126     const double aphi =
0127         tt::deltaPhi(track.phiT() - track.inv2R() * setup->chosenRofPhi() + region * setup->baseRegion());
0128     const double aTanLambda = track.cot();
0129     const double az0 = track.zT() - track.cot() * setup->chosenRofZ();
0130     const double ad0 = ttTrackRef->d0();
0131     static constexpr double aChi2xyfit = 0.;
0132     static constexpr double aChi2zfit = 0.;
0133     static constexpr double trkMVA1 = 0.;
0134     static constexpr double trkMVA2 = 0.;
0135     static constexpr double trkMVA3 = 0.;
0136     static constexpr unsigned int aHitpattern = 0;
0137     const unsigned int nPar = ttTrackRef->nFitPars();
0138     static const double Bfield = setup->bField();
0139     TTTrack<Ref_Phase2TrackerDigi_> ttTrack(
0140         aRinv, aphi, aTanLambda, az0, ad0, aChi2xyfit, aChi2zfit, trkMVA1, trkMVA2, trkMVA3, aHitpattern, nPar, Bfield);
0141     ttTrack.setStubRefs(ttStubRefs);
0142     ttTrack.setStubPtConsistency(
0143         StubPtConsistency::getConsistency(ttTrack, setup->trackerGeometry(), setup->trackerTopology(), Bfield, nPar));
0144     const int chi2B = tq->toBinChi2B(ttTrack.chi2Bend());
0145     const int chi2rphi = tq->toBinchi2rphi(trackchi2rphi);
0146     const int chi2rz = tq->toBinchi2rz(trackchi2rz);
0147     // load in bdt
0148     conifer::BDT<ap_fixed<10, 5>, ap_fixed<10, 5>> bdt(tq->model().fullPath());
0149     // collect features and classify using bdt
0150     const std::vector<ap_fixed<10, 5>>& output =
0151         bdt.decision_function({cot, z0, chi2B, nstub, n_missint, chi2rphi, chi2rz});
0152     const float mva = output[0].to_float();
0153     // fill frame
0154     TTBV ttBV = hitPattern;
0155     ttBV += TTBV(tq->toBinMVA(mva), widthMVA_);
0156     tq->format(VariableTQ::chi2rphi).attach(trackchi2rphi, ttBV);
0157     tq->format(VariableTQ::chi2rz).attach(trackchi2rz, ttBV);
0158     frame_ = ttBV.bs();
0159   }
0160 
0161   template <>
0162   DataFormat makeDataFormat<VariableTQ::m20>(const DataFormats* dataFormats, const ConfigTQ& iConfig) {
0163     const DataFormat phi = makeDataFormat<Variable::phi, Process::kf>(dataFormats->setup());
0164     const int width = iConfig.widthM20_;
0165     const double base = std::pow(phi.base(), 2) * pow(2., width - phi.width());
0166     const double range = base * std::pow(2, width);
0167     return DataFormat(false, width, base, range);
0168   }
0169   template <>
0170   DataFormat makeDataFormat<VariableTQ::m21>(const DataFormats* dataFormats, const ConfigTQ& iConfig) {
0171     const DataFormat z = makeDataFormat<Variable::z, Process::gp>(dataFormats->setup());
0172     const int width = iConfig.widthM21_;
0173     const double base = std::pow(z.base(), 2) * std::pow(2., width - z.width());
0174     const double range = base * std::pow(2, width);
0175     return DataFormat(false, width, base, range);
0176   }
0177   template <>
0178   DataFormat makeDataFormat<VariableTQ::invV0>(const DataFormats* dataFormats, const ConfigTQ& iConfig) {
0179     const DataFormat dPhi = makeDataFormat<Variable::dPhi, Process::ctb>(dataFormats->setup());
0180     const int width = iConfig.widthInvV0_;
0181     const double range = 4.0 / std::pow(dPhi.base(), 2);
0182     const double base = range * std::pow(2, -width);
0183     return DataFormat(false, width, base, range);
0184   }
0185   template <>
0186   DataFormat makeDataFormat<VariableTQ::invV1>(const DataFormats* dataFormats, const ConfigTQ& iConfig) {
0187     const DataFormat dZ = makeDataFormat<Variable::dZ, Process::ctb>(dataFormats->setup());
0188     const int width = iConfig.widthInvV1_;
0189     const double range = 4.0 / std::pow(dZ.base(), 2);
0190     const double base = range * std::pow(2, -width);
0191     return DataFormat(false, width, base, range);
0192   }
0193   template <>
0194   DataFormat makeDataFormat<VariableTQ::chi2rphi>(const DataFormats* dataFormats, const ConfigTQ& iConfig) {
0195     const int shift = iConfig.baseShiftchi2rphi_;
0196     const int width = iConfig.widthchi2rphi_;
0197     const double base = std::pow(2., shift);
0198     const double range = base * std::pow(2, width);
0199     return DataFormat(false, width, base, range);
0200   }
0201   template <>
0202   DataFormat makeDataFormat<VariableTQ::chi2rz>(const DataFormats* dataFormats, const ConfigTQ& iConfig) {
0203     const int shift = iConfig.baseShiftchi2rz_;
0204     const int width = iConfig.widthchi2rz_;
0205     const double base = std::pow(2., shift);
0206     const double range = base * std::pow(2, width);
0207     return DataFormat(false, width, base, range);
0208   }
0209 
0210   // Controls the conversion between TTTrack features and ML model training features
0211   std::vector<float> TrackQuality::featureTransform(TTTrack<Ref_Phase2TrackerDigi_>& aTrack,
0212                                                     std::vector<std::string> const& featureNames) const {
0213     // List input features for MVA in proper order below, the current features options are
0214     // {"phi", "eta", "z0", "bendchi2_bin", "nstub", "nlaymiss_interior", "chi2rphi_bin",
0215     // "chi2rz_bin"}
0216     //
0217     // To use more features, they must be created here and added to feature_map below
0218     std::vector<float> transformedFeatures;
0219     // Define feature map, filled as features are generated
0220     std::map<std::string, float> feature_map;
0221     // -------- calculate feature variables --------
0222     // calculate number of missed interior layers from hitpattern
0223     int tmp_trk_hitpattern = aTrack.hitPattern();
0224     int nbits = std::floor(std::log2(tmp_trk_hitpattern)) + 1;
0225     int lay_i = 0;
0226     int tmp_trk_nlaymiss_interior = 0;
0227     bool seq = false;
0228     for (int i = 0; i < nbits; i++) {
0229       lay_i = ((1 << i) & tmp_trk_hitpattern) >> i;  //0 or 1 in ith bit (right to left)
0230 
0231       if (lay_i && !seq)
0232         seq = true;  //sequence starts when first 1 found
0233       if (!lay_i && seq)
0234         tmp_trk_nlaymiss_interior++;
0235     }
0236     // binned chi2 variables
0237     int tmp_trk_bendchi2_bin = aTrack.getBendChi2Bits();
0238     int tmp_trk_chi2rphi_bin = aTrack.getChi2RPhiBits();
0239     int tmp_trk_chi2rz_bin = aTrack.getChi2RZBits();
0240     // get the nstub
0241     std::vector<TTStubRef> stubRefs = aTrack.getStubRefs();
0242     int tmp_trk_nstub = stubRefs.size();
0243     // get other variables directly from TTTrack
0244     float tmp_trk_z0 = aTrack.z0();
0245     float tmp_trk_z0_scaled = tmp_trk_z0 / abs(aTrack.minZ0);
0246     float tmp_trk_phi = aTrack.phi();
0247     float tmp_trk_eta = aTrack.eta();
0248     float tmp_trk_tanl = aTrack.tanL();
0249     // -------- fill the feature map ---------
0250     feature_map["nstub"] = float(tmp_trk_nstub);
0251     feature_map["z0"] = tmp_trk_z0;
0252     feature_map["z0_scaled"] = tmp_trk_z0_scaled;
0253     feature_map["phi"] = tmp_trk_phi;
0254     feature_map["eta"] = tmp_trk_eta;
0255     feature_map["nlaymiss_interior"] = float(tmp_trk_nlaymiss_interior);
0256     feature_map["bendchi2_bin"] = tmp_trk_bendchi2_bin;
0257     feature_map["chi2rphi_bin"] = tmp_trk_chi2rphi_bin;
0258     feature_map["chi2rz_bin"] = tmp_trk_chi2rz_bin;
0259     feature_map["tanl"] = tmp_trk_tanl;
0260     // fill tensor with track params
0261     transformedFeatures.reserve(featureNames.size());
0262     for (const std::string& feature : featureNames)
0263       transformedFeatures.push_back(feature_map[feature]);
0264     return transformedFeatures;
0265   }
0266 
0267   // Passed by reference a track without MVA filled, method fills the track's MVA field
0268   void TrackQuality::setL1TrackQuality(TTTrack<Ref_Phase2TrackerDigi_>& aTrack) const {
0269     // load in bdt
0270     conifer::BDT<float, float> bdt(this->model_.fullPath());
0271     // collect features and classify using bdt
0272     std::vector<float> inputs = featureTransform(aTrack, this->featureNames_);
0273     std::vector<float> output = bdt.decision_function(inputs);
0274     aTrack.settrkMVA1(1. / (1. + exp(-output.at(0))));
0275   }
0276 
0277 }  // namespace trackerTFP