src/Omtf/OmtfAngleConverter.cc

0001 /*
0002  * OmtfAngleConverter.cpp
0003  *
0004  *  Created on: Jan 14, 2019
0005  *      Author: kbunkow
0006  */
0007
0008 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/OmtfAngleConverter.h"
0009
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0012
0013 #include "Geometry/CSCGeometry/interface/CSCGeometry.h"
0014 #include "DataFormats/CSCDigi/interface/CSCConstants.h"
0015
0016 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0017 #include "L1Trigger/DTUtilities/interface/DTTrigGeom.h"
0018 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0019
0020 #include "DataFormats/CSCDigi/interface/CSCCorrelatedLCTDigi.h"
0021 #include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambPhDigi.h"
0022 #include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambThContainer.h"
0023 #include "DataFormats/RPCDigi/interface/RPCDigi.h"
0024
0025 namespace {
0026   template <typename T>
0027   int sgn(T val) {
0028     return (T(0) < val) - (val < T(0));
0029   }
0030
0031   //DT eta bins in the wheel +2
0032   const std::vector<float> bounds = {1.24, 1.14353, 1.09844, 1.05168, 1.00313, 0.952728, 0.90037, 0.8};
0033   //   0.8       -> 73
0034   //   0.85      -> 78
0035   //   0.9265    -> 85
0036   //   0.9779    -> 89.9 -> 90
0037   //   1.0274    -> 94.4 -> 94
0038   //   1.07506   -> 98.9 -> 99
0039   //   1.121     -> 103
0040   //   1.2       -> 110
0041   //   1.25      -> 115
0042   //
0043   // other (1.033) -> 1.033 -> 95
0044
0045   int etaVal2Bit(float eta) { return bounds.rend() - std::lower_bound(bounds.rbegin(), bounds.rend(), fabs(eta)); }
0046
0047   int etaBit2Code(unsigned int bit) {
0048     int code = 73;
0049     switch (bit) {
0050       case 0: {
0051         code = 115;
0052         break;
0053       }
0054       case 1: {
0055         code = 110;
0056         break;
0057       }
0058       case 2: {
0059         code = 103;
0060         break;
0061       }
0062       case 3: {
0063         code = 99;
0064         break;
0065       }
0066       case 4: {
0067         code = 94;
0068         break;
0069       }
0070       case 5: {
0071         code = 90;
0072         break;
0073       }
0074       case 6: {
0075         code = 85;
0076         break;
0077       }
0078       case 7: {
0079         code = 78;
0080         break;
0081       }
0082       case 8: {
0083         code = 73;
0084         break;
0085       }
0086       default: {
0087         code = 95;
0088         break;
0089       }
0090     }
0091     return code;
0092   }
0093
0094   int etaVal2Code(double etaVal) {
0095     int sign = sgn(etaVal);
0096     int bit = etaVal2Bit(fabs(etaVal));
0097     int code = etaBit2Code(bit);
0098     return sign * code;
0099   }
0100
0101   int etaKeyWG2Code(const CSCDetId &detId, uint16_t keyWG) {
0102     signed int etaCode = 121;
0103     if (detId.station() == 1 && detId.ring() == 2) {
0104       if (keyWG < 49)
0105         etaCode = 121;
0106       else if (keyWG <= 57)
0107         etaCode = etaBit2Code(0);
0108       else if (keyWG <= 63)
0109         etaCode = etaBit2Code(1);
0110     } else if (detId.station() == 1 && detId.ring() == 3) {
0111       if (keyWG <= 2)
0112         etaCode = etaBit2Code(2);
0113       else if (keyWG <= 8)
0114         etaCode = etaBit2Code(3);
0115       else if (keyWG <= 15)
0116         etaCode = etaBit2Code(4);
0117       else if (keyWG <= 23)
0118         etaCode = etaBit2Code(5);
0119       else if (keyWG <= 31)
0120         etaCode = etaBit2Code(6);
0121     } else if ((detId.station() == 2 || detId.station() == 3) && detId.ring() == 2) {
0122       if (keyWG < 24)
0123         etaCode = 121;
0124       else if (keyWG <= 29)
0125         etaCode = etaBit2Code(0);
0126       else if (keyWG <= 43)
0127         etaCode = etaBit2Code(1);
0128       else if (keyWG <= 49)
0129         etaCode = etaBit2Code(2);
0130       else if (keyWG <= 56)
0131         etaCode = etaBit2Code(3);
0132       else if (keyWG <= 63)
0133         etaCode = etaBit2Code(4);
0134     }
0135
0136     if (detId.endcap() == 2)
0137       etaCode *= -1;
0138     return etaCode;
0139   }
0140
0141 }  //namespace
0142
0143 OmtfAngleConverter::~OmtfAngleConverter() {}
0144
0145 ///////////////////////////////////////
0146 ///////////////////////////////////////
0147 int OmtfAngleConverter::getGlobalEta(const DTChamberId dTChamberId,
0148                                      const L1MuDTChambThContainer *dtThDigis,
0149                                      int bxNum) const {
0150   //const DTChamberId dTChamberId(aDigi.whNum(),aDigi.stNum(),aDigi.scNum()+1);
0151   DTTrigGeom trig_geom(_geodt->chamber(dTChamberId), false);
0152
0153   // super layer one is the theta superlayer in a DT chamber
0154   // station 4 does not have a theta super layer
0155   // the BTI index from the theta trigger is an OR of some BTI outputs
0156   // so, we choose the BTI that's in the middle of the group
0157   // as the BTI that we get theta from
0158   // TODO:::::>>> need to make sure this ordering doesn't flip under wheel sign
0159   const int NBTI_theta = ((dTChamberId.station() != 4) ? trig_geom.nCell(2) : trig_geom.nCell(3));
0160
0161   const L1MuDTChambThDigi *theta_segm =
0162       dtThDigis->chThetaSegm(dTChamberId.wheel(), dTChamberId.station(), dTChamberId.sector() - 1, bxNum);
0163
0164   int bti_group = -1;
0165   if (theta_segm) {
0166     for (unsigned int i = 0; i < 7; ++i)
0167       if (theta_segm->position(i) && bti_group < 0)
0168         bti_group = i;
0169       else if (theta_segm->position(i) && bti_group > -1)
0170         bti_group = 511;
0171   }
0172
0173   int iEta = 0;
0174   if (bti_group == 511)
0175     iEta = 95;
0176   else if (bti_group == -1 && dTChamberId.station() == 1)
0177     iEta = 92;
0178   else if (bti_group == -1 && dTChamberId.station() == 2)
0179     iEta = 79;
0180   else if (bti_group == -1 && dTChamberId.station() == 3)
0181     iEta = 75;
0182   else if (dTChamberId.station() != 4 && bti_group >= 0) {
0183     unsigned bti_actual = bti_group * NBTI_theta / 7 + NBTI_theta / 14 + 1;
0184     DTBtiId thetaBTI = DTBtiId(dTChamberId, 2, bti_actual);
0185     GlobalPoint theta_gp = trig_geom.CMSPosition(thetaBTI);
0186     iEta = etaVal2Code(fabs(theta_gp.eta()));
0187   }
0188   int signEta = sgn(dTChamberId.wheel());
0189   iEta *= signEta;
0190   return iEta;
0191 }
0192
0193 ///////////////////////////////////////
0194 ///////////////////////////////////////
0195 int OmtfAngleConverter::getGlobalEta(unsigned int rawid, const CSCCorrelatedLCTDigi &aDigi) const {
0196   ///Code taken from GeometryTranslator.
0197   ///Will be replaced by direct CSC phi local to global scale
0198   ///transformation as used in FPGA implementation
0199
0200   // a lot of this is transcription and consolidation of the CSC
0201   // global phi calculation code
0202   // this works directly with the geometry
0203   // rather than using the old phi luts
0204   const CSCDetId id(rawid);
0205   // we should change this to weak_ptrs at some point
0206   // requires introducing std::shared_ptrs to geometry
0207   auto chamb = _geocsc->chamber(id);
0208   auto layer_geom = chamb->layer(CSCConstants::KEY_ALCT_LAYER)->geometry();
0209   auto layer = chamb->layer(CSCConstants::KEY_ALCT_LAYER);
0210
0211   const uint16_t halfstrip = aDigi.getStrip();
0212   //const uint16_t pattern = aDigi.getPattern();
0213   const uint16_t keyWG = aDigi.getKeyWG();
0214   //const unsigned maxStrips = layer_geom->numberOfStrips();
0215
0216   // so we can extend this later
0217   // assume TMB2007 half-strips only as baseline
0218   double offset = 0.0;
0219   //K.B. CSCPatternLUT is removed since CMSSW_11_2_1, but it looks that this offset is effectively not needed here
0220   //offset = CSCPatternLUT::get2007Position(pattern);
0221
0222   const unsigned halfstrip_offs = unsigned(0.5 + halfstrip + offset);
0223   const unsigned strip = halfstrip_offs / 2 + 1;  // geom starts from 1
0224
0225   // the rough location of the hit at the ALCT key layer
0226   // we will refine this using the half strip information
0227   const LocalPoint coarse_lp = layer_geom->stripWireGroupIntersection(strip, keyWG);
0228   const GlobalPoint coarse_gp = layer->surface().toGlobal(coarse_lp);
0229
0230   // the strip width/4.0 gives the offset of the half-strip
0231   // center with respect to the strip center
0232   const double hs_offset = layer_geom->stripPhiPitch() / 4.0;
0233
0234   // determine handedness of the chamber
0235   const bool ccw = isCSCCounterClockwise(layer);
0236   // we need to subtract the offset of even half strips and add the odd ones
0237   const double phi_offset = ((halfstrip_offs % 2 ? 1 : -1) * (ccw ? -hs_offset : hs_offset));
0238
0239   // the global eta calculation uses the middle of the strip
0240   // so no need to increment it
0241   const GlobalPoint final_gp(
0242       GlobalPoint::Polar(coarse_gp.theta(), (coarse_gp.phi().value() + phi_offset), coarse_gp.mag()));
0243
0244   //  LogTrace("l1tOmtfEventPrint")<<id<<" st: " << id.station()<< "ri: "<<id.ring()<<" eta: " <<  final_gp.eta()
0245   //           <<" etaCode_simple: " <<  etaVal2Code( final_gp.eta() )<< " KW: "<<keyWG <<" etaKeyWG2Code: "<<etaKeyWG2Code(id,keyWG)<< std::endl;
0246   //  int station = (id.endcap()==1) ? id.station() : -id.station();
0247   //  LogTrace("l1tOmtfEventPrint")<<"ETA_CSC: " << station <<" "<<id.ring()<<" "<< final_gp.eta()<<" "<<keyWG <<" "<< etaKeyWG2Code(id,keyWG) << std::endl;
0248
0249   return etaKeyWG2Code(id, keyWG);
0250 }
0251 ///////////////////////////////////////
0252 ///////////////////////////////////////
0253 int OmtfAngleConverter::getGlobalEtaRpc(unsigned int rawid, const unsigned int &strip) const {
0254   const RPCDetId id(rawid);
0255   auto roll = _georpc->roll(id);
0256   const LocalPoint lp = roll->centreOfStrip((int)strip);
0257   const GlobalPoint gp = roll->toGlobal(lp);
0258
0259   return etaVal2Code(gp.eta());
0260 }
0261
0262 ///////////////////////////////////////
0263 ///////////////////////////////////////