Project CMSSW displayed by LXR CMSSW_14_2_X_2024-09-06-2300/​L1Trigger/​DTTriggerPhase2/​src/​GlobalCoordsObtainer.cc	Source navigation Diff markup Identifier search General search
	Version: CMSSW_14_2_X_2024-09-06-2300 CMSSW_14_2_X_2024-09-05-2300 CMSSW_14_2_X_2024-09-04-2300 CMSSW_14_2_X_2024-09-03-2300 CMSSW_14_2_X_2024-09-01-2300 Revert   Change

File indexing completed on 2024-09-04 04:35:15

 #include "L1Trigger/DTTriggerPhase2/interface/GlobalCoordsObtainer.h"
 
 using namespace edm;
 using namespace cmsdt;
 using namespace std;
 
 // ============================================================================
 // Constructors and destructor
 // ============================================================================
 GlobalCoordsObtainer::GlobalCoordsObtainer(const ParameterSet& pset) {
   global_coords_filename_ = pset.getParameter<edm::FileInPath>("global_coords_filename");
   std::ifstream ifin3(global_coords_filename_.fullPath());
 
   if (ifin3.fail()) {
     throw cms::Exception("Missing Input File")
         << "GlobalCoordsObtainer::GlobalCoordsObtainer() -  Cannot find " << global_coords_filename_.fullPath();
   }
 
   int wh, st, se, sl;
   double perp, x_phi0;
   std::string line;
 
   global_constant_per_sl sl1_constants;
   global_constant_per_sl sl3_constants;
 
   while (ifin3.good()) {
     ifin3 >> wh >> st >> se >> sl >> perp >> x_phi0;
 
     if (sl == 1) {
       sl1_constants.perp = perp;
       sl1_constants.x_phi0 = x_phi0;
     } else {
       sl3_constants.perp = perp;
       sl3_constants.x_phi0 = x_phi0;
 
       DTChamberId ChId(wh, st, se);
       global_constants.push_back({ChId.rawId(), sl1_constants, sl3_constants});
     }
   }
 
   int maxdrift;
   maxdrift_filename_ = pset.getParameter<edm::FileInPath>("maxdrift_filename");
   std::ifstream ifind(maxdrift_filename_.fullPath());
   if (ifind.fail()) {
     throw cms::Exception("Missing Input File")
         << "MPSLFilter::MPSLFilter() -  Cannot find " << maxdrift_filename_.fullPath();
   }
   while (ifind.good()) {
     ifind >> wh >> st >> se >> maxdrift;
     maxdriftinfo_[wh][st][se] = maxdrift;
   }
 }
 
 GlobalCoordsObtainer::~GlobalCoordsObtainer() {}
 
 void GlobalCoordsObtainer::generate_luts() {
   for (auto& global_constant : global_constants) {
     int sgn = 1;
     DTChamberId ChId(global_constant.chid);
     // typical hasPosRF function
     if (ChId.wheel() > 0 || (ChId.wheel() == 0 && ChId.sector() % 4 > 1)) {
       sgn = -1;
     }
 
     max_drift_tdc = maxdriftinfo_[ChId.wheel() + 2][ChId.station() - 1][ChId.sector() - 1];
     int SEMICELL_IN_TDC_COUNTS = max_drift_tdc;  //495 en un mundo peor y mas oscuro
     int SL1_SEMICELL_OFFSET = 96;
 
     //Read atan function from LUTs
     auto phi =
         calc_atan_lut(11,
                       6,
                       21. / SEMICELL_IN_TDC_COUNTS / ((global_constant.sl1.perp + global_constant.sl3.perp) / .2),
                       (global_constant.sl1.x_phi0 - 2.1 * SL1_SEMICELL_OFFSET) /
                           ((global_constant.sl1.perp + global_constant.sl3.perp) / 2),
                       1. / std::pow(2, 17),
                       9,
                       4,
                       11,
                       21,
                       sgn);
 
     auto phib =
         calc_atan_lut(9, 5, 21. / SEMICELL_IN_TDC_COUNTS / (6.5 * 16), 0., 4. / std::pow(2, 13), 9, 3, 10, 16, sgn);
 
     luts[global_constant.chid] = {phi, phib};
   }
 }
 
 std::map<int, lut_value> GlobalCoordsObtainer::calc_atan_lut(int msb_num,
                                                              int lsb_num,
                                                              double in_res,
                                                              double abscissa_0,
                                                              double out_res,
                                                              int a_extra_bits,
                                                              int b_extra_bits,
                                                              int a_size,
                                                              int b_size,
                                                              int sgn) {
   // Calculates the coefficients needed to calculate the arc-tan function in fw
   // by doing a piece-wise linear approximation.
   // In fw, input (x) and output (y) are integers, these conversions are needed
   //   t = x*in_res - abscissa_0
   //   phi = arctan(t)
   //   y = phi/out_res
   //   => y = arctan(x*in_res - abcissa_0)/out_res
   // The linear function is approximated as
   //   y = a*x_lsb + b
   // Where a, b = func(x_msb) are the coefficients stored in the lut
 
   // a is stored as unsigned, b as signed, with their respective sizes a_size, b_size,
   // previously shifted left by a_extra_bits and b_extra_bits, respectively
 
   long int a_min = -std::pow(2, (a_size - 1));
   long int a_max = std::pow(2, (a_size - 1)) - 1;
   long int b_min = -std::pow(2, (b_size - 1));
   long int b_max = std::pow(2, (b_size - 1)) - 1;
 
   std::map<int, lut_value> lut;
 
   for (long int x_msb = -(long int)std::pow(2, msb_num - 1); x_msb < (long int)std::pow(2, msb_num - 1); x_msb++) {
     int x1 = ((x_msb) << lsb_num);
     int x2 = ((x_msb + 1) << lsb_num) - 1;
 
     double t1 = x1 * in_res - abscissa_0;
     double t2 = x2 * in_res - abscissa_0;
 
     double phi1 = sgn * atan(t1);
     double phi2 = sgn * atan(t2);
 
     double y1 = phi1 / out_res;
     double y2 = phi2 / out_res;
 
     // we want to find a, b so that the error in the extremes is the same as the error in the center
     // so the error in the extremes will be the same, so the "a" is determined by those points
     double a = (y2 - y1) / (x2 - x1);
 
     // by derivating the error function and equaling to 0, you get this is the point
     // towards the interval center with the highest error
     // err_f = y - (a*x+b) = sgn*arctan(x*in_res - abcissa_0)/out_res - (a*x+b)
     // d(err_f)/dx = sgn*(1/(1+(x*in_res - abcissa_0)^2))*in_res/out_res - a
     // d(err_f)/dx = 0 => x_max_err = (sqrt(in_res/out_res/a-1) + abscissa_0)/in_res
     // There is sign ambiguity in the sqrt operation. The sqrt has units of t (adimensional).
     // It is resolved by setting the sqrt to have the same sign as (t1+t2)/2
 
     double t_max_err = sqrt(sgn * in_res / out_res / a - 1);
     if ((t1 + t2) < 0) {
       t_max_err *= -1;
     }
 
     double x_max_err = (t_max_err + abscissa_0) / in_res;
     double phi_max_err = sgn * atan(t_max_err);
     double y_max_err = phi_max_err / out_res;
 
     // once you have the point of max error, the "b" parameter is chosen as the average between
     // those two numbers, which makes the error at the center be equal in absolute value
     // to the error in the extremes
     // units: rad
 
     double b = (y1 + y_max_err - a * (x_max_err - x1)) / 2;
 
     // increase b in 1/2 of y_lsb, so that fw truncate operation on the of the result
     // is equivalent to a round function instead of a floor function
     b += 0.5;
 
     // shift left and round
     long int a_int = (long int)(round(a * (pow(2, a_extra_bits))));
     long int b_int = (long int)(round(b * (pow(2, b_extra_bits))));
 
     // tuck a, b constants into the bit size of the output (un)signed integer
     std::vector<long int> as = {a_min, a_int, a_max};
     std::vector<long int> bs = {b_min, b_int, b_max};
 
     std::sort(as.begin(), as.end());
     std::sort(bs.begin(), bs.end());
 
     a_int = as.at(1);
     b_int = bs.at(1);
 
     // // convert a, b to two's complement
     // auto a_signed = a_int % (long int) (pow(2, a_size));
     // auto b_signed = b_int % (long int) (pow(2, b_size));
 
     // convert x_msb to two's complement signed
     int index = to_two_comp(x_msb, msb_num);
     lut[index] = {a_int, b_int};
   }
   return lut;
 }
 
 std::vector<double> GlobalCoordsObtainer::get_global_coordinates(uint32_t chid, int sl, int x, int tanpsi) {
   // Depending on the type of primitive (SL1, SL3 or correlated), choose the
   // appropriate input data (x, tanpsi) from the input primitive data structure
   // and the corresponding phi-lut from the 3 available options
 
   auto phi_lut = &luts[chid].phi;  //One parameter to rule them all
   auto phib_lut = &luts[chid].phib;
 
   // x and slope are given in two's complement in fw
   x = to_two_comp(x, X_SIZE);
   tanpsi = to_two_comp(tanpsi, TANPSI_SIZE);
 
   // Slice x and tanpsi
   // Both x and tanpsi are represented in vhdl as signed, this means their values
   // are stored as two's complement.
 
   // The MSB part is going to be used to index the luts and obtain a and b parameters
   // Converting the upper part of the signed to an integer (with sign).
 
   int x_msb = x >> (X_SIZE - PHI_LUT_ADDR_WIDTH);
   x_msb = from_two_comp(x_msb, PHI_LUT_ADDR_WIDTH);
 
   int tanpsi_msb = tanpsi >> (TANPSI_SIZE - PHIB_LUT_ADDR_WIDTH);
   tanpsi_msb = from_two_comp(tanpsi_msb, PHIB_LUT_ADDR_WIDTH);
 
   // The LSB part can be sliced right away because it must yield a positive integer
   int x_lsb = x & (int)(std::pow(2, (X_SIZE - PHI_LUT_ADDR_WIDTH)) - 1);
   int tanpsi_lsb = tanpsi & (int)(std::pow(2, (TANPSI_SIZE - PHIB_LUT_ADDR_WIDTH)) - 1);
 
   // Index the luts wiht the MSB parts already calculated
   auto phi_lut_q = phi_lut->at(to_two_comp(x_msb, PHI_LUT_ADDR_WIDTH));
   auto phib_lut_q = phib_lut->at(to_two_comp(tanpsi_msb, PHIB_LUT_ADDR_WIDTH));
 
   // Separate this data into the coefficients a and b
   auto phi_lut_a = phi_lut_q.a;
   auto phi_lut_b = phi_lut_q.b;
   auto phib_lut_a = phib_lut_q.a;
   auto phib_lut_b = phib_lut_q.b;
 
   // Do the linear piece-wise operations
   // At this point all variables that can be negative have already been converted
   // so will yield negative values when necessary
   int phi_uncut = (phi_lut_b << PHI_B_SHL_BITS) + x_lsb * phi_lut_a;
   int psi_uncut = (phib_lut_b << PHIB_B_SHL_BITS) + tanpsi_lsb * phib_lut_a;
 
   // Trim phi to its final size
   int phi = (phi_uncut >> PHI_MULT_SHR_BITS);
 
   // Calculate phi_bending from the uncut version of phi and psi, and the trim it to size
   int phib_uncut = psi_uncut - (phi_uncut >> (PHI_PHIB_RES_DIFF_BITS + PHI_MULT_SHR_BITS - PHIB_MULT_SHR_BITS));
   int phib = (phib_uncut >> PHIB_MULT_SHR_BITS);
 
   double phi_f = (double)phi * PHI_SIZE;
   double phib_f = (double)phib * PHIB_SIZE;
 
   return std::vector({phi_f, phib_f});
 }

This page was automatically generated by the 2.2.1 LXR engine. The LXR team