L1TMuonEndCap/src/PtAssignmentEngineDxy.cc

0001 #include "L1Trigger/L1TMuonEndCap/interface/PtAssignmentEngineDxy.h"
0002
0003 #include <cassert>
0004 #include <iostream>
0005 #include <sstream>
0006
0007 #include "helper.h"  // assert_no_abort
0008
0009 PtAssignmentEngineDxy::PtAssignmentEngineDxy() : graphDefDxy_(nullptr), sessionDxy_(nullptr) {}
0010
0011 PtAssignmentEngineDxy::~PtAssignmentEngineDxy() {
0012   if (sessionDxy_ != nullptr) {
0013     tensorflow::closeSession(sessionDxy_);
0014   }
0015   delete graphDefDxy_;
0016 }
0017
0018 void PtAssignmentEngineDxy::configure(int verbose, const std::string pbFileNameDxy) {
0019   verbose_ = verbose;
0020
0021   pbFileNameDxy_ = pbFileNameDxy;
0022   std::string pbFilePathDxy_ = "L1Trigger/L1TMuon/data/emtf_luts/" + pbFileNameDxy_;
0023
0024   inputNameDxy_ = "input1";
0025   outputNamesDxy_ = {"Identity"};
0026
0027   if (graphDefDxy_ == nullptr) {
0028     graphDefDxy_ = tensorflow::loadGraphDef(edm::FileInPath(pbFilePathDxy_).fullPath());
0029   }
0030   emtf_assert(graphDefDxy_ != nullptr);
0031
0032   if (sessionDxy_ == nullptr) {
0033     sessionDxy_ = tensorflow::createSession(graphDefDxy_);
0034   }
0035
0036   emtf_assert(sessionDxy_ != nullptr);
0037 }
0038
0039 const PtAssignmentEngineAux2017& PtAssignmentEngineDxy::aux() const {
0040   static const PtAssignmentEngineAux2017 instance;
0041   return instance;
0042 }
0043
0044 void PtAssignmentEngineDxy::calculate_pt_dxy(const EMTFTrack& track,
0045                                              emtf::Feature& feature,
0046                                              emtf::Prediction& prediction) const {
0047   // This is called for each track instead of for entire track collection as was done in Phase-2 implementation
0048   preprocessing_dxy(track, feature);
0049   call_tensorflow_dxy(feature, prediction);
0050   return;
0051 }
0052
0053 void PtAssignmentEngineDxy::preprocessing_dxy(const EMTFTrack& track, emtf::Feature& feature) const {
0054   // Mimic Phase-1 EMTF input calculations
0055   // 6 delta Phis: S1-S2, S1-S3, S1-S4, S2-S3, S2-S4, S3-S4
0056   // 6 delta Thetas: S1-S2, S1-S3, S1-S4, S2-S3, S2-S4, S3-S4
0057   // 4 bends : set to zero if no CSC hit and thus RPC hit is used
0058   // 1 FR bit: for ME1 only
0059   // 1 Ring bit: for ME1 only
0060   // 1 track Theta taken from stub coordinate in ME2, ME3, ME4 (in this priority)
0061   // 4 RPC bits indicating if ME or RE hit was used in each station (S1, S2, S3, S4)
0062   // Total: 23 variables
0063   std::array<float, 6> x_dphi;
0064   std::array<float, 6> x_dtheta;
0065   std::array<float, 4> x_bend_emtf;
0066   std::array<float, 1> x_fr_emtf;
0067   std::array<float, 1> x_trk_theta;
0068   std::array<float, 1> x_me11ring;
0069   std::array<float, 4> x_rpcbit;
0070
0071   // Initialize to zeros
0072   x_dphi.fill(0);
0073   x_dtheta.fill(0);
0074   //
0075   x_bend_emtf.fill(0);
0076   x_fr_emtf.fill(0);
0077   x_trk_theta.fill(0);
0078   x_me11ring.fill(0);
0079   x_rpcbit.fill(0);
0080
0081   EMTFPtLUT data = track.PtLUT();
0082
0083   const int invalid_dtheta = 127;
0084   const int invalid_dphi = 8191;
0085
0086   // // Variables to extract from the PtLUT
0087   int dPhi_12, dPhi_13, dPhi_14, dPhi_23, dPhi_24, dPhi_34;
0088   int dTh_12, dTh_13, dTh_14, dTh_23, dTh_24, dTh_34;
0089   int fr_1;
0090   int bend_1, bend_2, bend_3, bend_4;
0091   int rpc_1, rpc_2, rpc_3, rpc_4;
0092   int St1_ring2 = data.st1_ring2;
0093
0094   int pat1 = -99, pat2 = -99, pat3 = -99, pat4 = -99;
0095
0096   // // Which stations have hits
0097   int st1 = (track.Mode() >= 8);
0098   int st2 = ((track.Mode() % 8) >= 4);
0099   int st3 = ((track.Mode() % 4) >= 2);
0100   int st4 = ((track.Mode() % 2) == 1);
0101
0102   // Get valid pattern values
0103   if (st1)
0104     pat1 = data.cpattern[0];
0105   if (st2)
0106     pat2 = data.cpattern[1];
0107   if (st3)
0108     pat3 = data.cpattern[2];
0109   if (st4)
0110     pat4 = data.cpattern[3];
0111
0112   // F/R bit
0113   fr_1 = data.fr[0];
0114
0115   // RPC hit in station
0116   rpc_1 = (st1 ? (pat1 == 0) : 0);
0117   rpc_2 = (st2 ? (pat2 == 0) : 0);
0118   rpc_3 = (st3 ? (pat3 == 0) : 0);
0119   rpc_4 = (st4 ? (pat4 == 0) : 0);
0120
0121   // Calculate bends from patterns
0122   bend_1 = aux().calcBendFromPattern(pat1, track.Endcap());
0123   bend_2 = aux().calcBendFromPattern(pat2, track.Endcap());
0124   bend_3 = aux().calcBendFromPattern(pat3, track.Endcap());
0125   bend_4 = aux().calcBendFromPattern(pat4, track.Endcap());
0126
0127   // Invalid bend value is 0 in the NN
0128   if (bend_1 == -99)
0129     bend_1 = 0;
0130   if (bend_2 == -99)
0131     bend_2 = 0;
0132   if (bend_3 == -99)
0133     bend_3 = 0;
0134   if (bend_4 == -99)
0135     bend_4 = 0;
0136
0137   // In the emulator RPCs get assigned abs(bend) = 5. This needs to be 0 for the NN.
0138   if (std::abs(bend_1) == 5 && rpc_1 == 1)
0139     bend_1 = 0;
0140   if (std::abs(bend_2) == 5 && rpc_2 == 1)
0141     bend_2 = 0;
0142   if (std::abs(bend_3) == 5 && rpc_3 == 1)
0143     bend_3 = 0;
0144   if (std::abs(bend_4) == 5 && rpc_4 == 1)
0145     bend_4 = 0;
0146
0147   // Calculate delta phi
0148   dPhi_12 = (data.delta_ph[0] != invalid_dphi) ? data.delta_ph[0] * (data.sign_ph[0] ? 1 : -1) : 0;
0149   dPhi_13 = (data.delta_ph[1] != invalid_dphi) ? data.delta_ph[1] * (data.sign_ph[1] ? 1 : -1) : 0;
0150   dPhi_14 = (data.delta_ph[2] != invalid_dphi) ? data.delta_ph[2] * (data.sign_ph[2] ? 1 : -1) : 0;
0151   dPhi_23 = (data.delta_ph[3] != invalid_dphi) ? data.delta_ph[3] * (data.sign_ph[3] ? 1 : -1) : 0;
0152   dPhi_24 = (data.delta_ph[4] != invalid_dphi) ? data.delta_ph[4] * (data.sign_ph[4] ? 1 : -1) : 0;
0153   dPhi_34 = (data.delta_ph[5] != invalid_dphi) ? data.delta_ph[5] * (data.sign_ph[5] ? 1 : -1) : 0;
0154
0155   // Calculate delta theta
0156   dTh_12 = (data.delta_th[0] != invalid_dtheta) ? data.delta_th[0] * (data.sign_th[0] ? 1 : -1) : 0;
0157   dTh_13 = (data.delta_th[1] != invalid_dtheta) ? data.delta_th[1] * (data.sign_th[1] ? 1 : -1) : 0;
0158   dTh_14 = (data.delta_th[2] != invalid_dtheta) ? data.delta_th[2] * (data.sign_th[2] ? 1 : -1) : 0;
0159   dTh_23 = (data.delta_th[3] != invalid_dtheta) ? data.delta_th[3] * (data.sign_th[3] ? 1 : -1) : 0;
0160   dTh_24 = (data.delta_th[4] != invalid_dtheta) ? data.delta_th[4] * (data.sign_th[4] ? 1 : -1) : 0;
0161   dTh_34 = (data.delta_th[5] != invalid_dtheta) ? data.delta_th[5] * (data.sign_th[5] ? 1 : -1) : 0;
0162
0163   // Set dPhi and dTheta values to 0 if there was no hit in the station
0164   if (!st1) {
0165     dPhi_12 = 0;
0166     dPhi_13 = 0;
0167     dPhi_14 = 0;
0168
0169     dTh_12 = 0;
0170     dTh_13 = 0;
0171     dTh_14 = 0;
0172   }
0173   if (!st2) {
0174     dPhi_12 = 0;
0175     dPhi_23 = 0;
0176     dPhi_24 = 0;
0177
0178     dTh_12 = 0;
0179     dTh_23 = 0;
0180     dTh_24 = 0;
0181   }
0182   if (!st3) {
0183     dPhi_13 = 0;
0184     dPhi_23 = 0;
0185     dPhi_34 = 0;
0186
0187     dTh_13 = 0;
0188     dTh_23 = 0;
0189     dTh_34 = 0;
0190   }
0191   if (!st4) {
0192     dPhi_14 = 0;
0193     dPhi_24 = 0;
0194     dPhi_34 = 0;
0195
0196     dTh_14 = 0;
0197     dTh_24 = 0;
0198     dTh_34 = 0;
0199   }
0200
0201   // Set NN inputs
0202
0203   // NN was trained with the wrong sign convention. TO BE CHANGED LATER!
0204   x_dphi[0] = dPhi_12;
0205   x_dphi[1] = dPhi_13;
0206   x_dphi[2] = dPhi_14;
0207   x_dphi[3] = dPhi_23;
0208   x_dphi[4] = dPhi_24;
0209   x_dphi[5] = dPhi_34;
0210
0211   // NN was trained with the wrong sign convention. TO BE CHANGED LATER!
0212   x_dtheta[0] = dTh_12;
0213   x_dtheta[1] = dTh_13;
0214   x_dtheta[2] = dTh_14;
0215   x_dtheta[3] = dTh_23;
0216   x_dtheta[4] = dTh_24;
0217   x_dtheta[5] = dTh_34;
0218
0219   // NN was trained with the wrong sign convention. TO BE CHANGED LATER!
0220   x_bend_emtf[0] = bend_1;
0221   x_bend_emtf[1] = bend_2;
0222   x_bend_emtf[2] = bend_3;
0223   x_bend_emtf[3] = bend_4;
0224
0225   x_fr_emtf[0] = fr_1;
0226   x_trk_theta[0] = track.Theta_fp();
0227   x_me11ring[0] = St1_ring2;
0228
0229   x_rpcbit[0] = rpc_1;
0230   x_rpcbit[1] = rpc_2;
0231   x_rpcbit[2] = rpc_3;
0232   x_rpcbit[3] = rpc_4;
0233
0234   feature = {{x_dphi[0],      x_dphi[1],      x_dphi[2],      x_dphi[3],      x_dphi[4],    x_dphi[5],
0235               x_dtheta[0],    x_dtheta[1],    x_dtheta[2],    x_dtheta[3],    x_dtheta[4],  x_dtheta[5],
0236               x_bend_emtf[0], x_bend_emtf[1], x_bend_emtf[2], x_bend_emtf[3], x_fr_emtf[0], x_trk_theta[0],
0237               x_me11ring[0],  x_rpcbit[0],    x_rpcbit[1],    x_rpcbit[2],    x_rpcbit[3]}};
0238   return;
0239 }
0240
0241 void PtAssignmentEngineDxy::call_tensorflow_dxy(const emtf::Feature& feature, emtf::Prediction& prediction) const {
0242   tensorflow::Tensor input(tensorflow::DT_FLOAT, {1, emtf::NUM_FEATURES});
0243   std::vector<tensorflow::Tensor> outputs;
0244   emtf_assert(feature.size() == emtf::NUM_FEATURES);
0245
0246   float* d = input.flat<float>().data();
0247   std::copy(feature.begin(), feature.end(), d);
0248   tensorflow::run(sessionDxy_, {{inputNameDxy_, input}}, outputNamesDxy_, &outputs);
0249   emtf_assert(outputs.size() == 1);
0250   emtf_assert(prediction.size() == emtf::NUM_PREDICTIONS);
0251
0252   const float reg_pt_scale = 100.0;  // a scale factor applied to regression during training
0253   const float reg_dxy_scale = 1.0;   // a scale factor applied to regression during training
0254
0255   prediction.at(0) = outputs[0].matrix<float>()(0, 0);
0256   prediction.at(1) = outputs[0].matrix<float>()(0, 1);
0257
0258   // Remove scale factor used during training
0259   prediction.at(0) /= reg_pt_scale;
0260   prediction.at(1) /= reg_dxy_scale;
0261   return;
0262 }