HcalAlgos/src/HcalDbHardcode.cc

0001 //
0002 // F.Ratnikov (UMd), Dec 14, 2005
0003 //
0004 #include <vector>
0005 #include <string>
0006 #include <cmath>
0007 #include <sstream>
0008 #include <iostream>
0009 #include <memory>
0010 #include <cassert>
0011
0012 #include "CLHEP/Random/RandGauss.h"
0013 #include "CalibCalorimetry/HcalAlgos/interface/HcalDbHardcode.h"
0014 #include "CalibFormats/HcalObjects/interface/HcalSiPMType.h"
0015 #include "DataFormats/HcalDigi/interface/HcalQIENum.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017
0018 HcalDbHardcode::HcalDbHardcode()
0019     :                                              //"generic" set of conditions
0020       theDefaultParameters_(3.0,                   //pedestal
0021                             0.5,                   //pedestal width
0022                             {0.2, 0.2},            //gains
0023                             {0.0, 0.0},            //gain widths
0024                             0,                     //ZS threshold
0025                             0,                     //QIE type
0026                             {0.0, 0.0, 0.0, 0.0},  //QIE offsets
0027                             {0.9, 0.9, 0.9, 0.9},  //QIE slopes
0028                             125,                   //MC shape
0029                             105,                   //Reco shape
0030                             0.0,                   //photoelectronsToAnalog
0031                             {0.0},                 //dark current
0032                             {0.0},                 //noise correlation
0033                             0.0,                   //PF noise threshold
0034                             0.1                    //PF seed threshold
0035                             ),
0036       setHB_(false),
0037       setHE_(false),
0038       setHF_(false),
0039       setHO_(false),
0040       setHBUpgrade_(false),
0041       setHEUpgrade_(false),
0042       setHFUpgrade_(false),
0043       useHBUpgrade_(false),
0044       useHEUpgrade_(false),
0045       useHOUpgrade_(true),
0046       useHFUpgrade_(false),
0047       testHFQIE10_(false),
0048       testHEPlan1_(false) {}
0049
0050 const HcalHardcodeParameters& HcalDbHardcode::getParameters(HcalGenericDetId fId) const {
0051   if (fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel) {
0052     if (useHBUpgrade_ && setHBUpgrade_)
0053       return theHBUpgradeParameters_;
0054     else if (!useHBUpgrade_ && setHB_)
0055       return theHBParameters_;
0056     else
0057       return theDefaultParameters_;
0058   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenEndcap) {
0059     bool b_isHEPlan1 = testHEPlan1_ ? isHEPlan1(fId) : false;
0060     if ((useHEUpgrade_ || b_isHEPlan1) && setHEUpgrade_)
0061       return theHEUpgradeParameters_;
0062     else if (!useHEUpgrade_ && !b_isHEPlan1 && setHE_)
0063       return theHEParameters_;
0064     else
0065       return theDefaultParameters_;
0066   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenForward) {
0067     if (useHFUpgrade_ && setHFUpgrade_)
0068       return theHFUpgradeParameters_;
0069     else if (testHFQIE10_ && fId.isHcalDetId()) {
0070       HcalDetId hid(fId);
0071       //special mixed case for HF 2016
0072       if (hid.iphi() == 39 && hid.zside() == 1 &&
0073           (hid.depth() >= 3 || (hid.depth() == 2 && (hid.ieta() == 30 || hid.ieta() == 34))) && setHFUpgrade_)
0074         return theHFUpgradeParameters_;
0075       else if (setHF_)
0076         return theHFParameters_;
0077       else
0078         return theDefaultParameters_;
0079     } else if (!useHFUpgrade_ && setHF_)
0080       return theHFParameters_;
0081     else
0082       return theDefaultParameters_;
0083   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenOuter) {
0084     if (setHO_)
0085       return theHOParameters_;
0086     else
0087       return theDefaultParameters_;
0088   } else
0089     return theDefaultParameters_;
0090 }
0091
0092 const int HcalDbHardcode::getGainIndex(HcalGenericDetId fId) const {
0093   int index = 0;
0094   if (fId.genericSubdet() == HcalGenericDetId::HcalGenOuter) {
0095     HcalDetId hid(fId);
0096     if ((hid.ieta() > -5) && (hid.ieta() < 5))
0097       index = 0;
0098     else
0099       index = 1;
0100   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenForward) {
0101     HcalDetId hid(fId);
0102     if (hid.depth() % 2 == 1)
0103       index = 0;  //depths 1,3
0104     else if (hid.depth() % 2 == 0)
0105       index = 1;  //depths 2,4
0106   }
0107   return index;
0108 }
0109
0110 HcalPedestal HcalDbHardcode::makePedestal(
0111     HcalGenericDetId fId, bool fSmear, bool eff, const HcalTopology* topo, double intlumi) {
0112   HcalPedestalWidth width = makePedestalWidth(fId, eff, topo, intlumi);
0113   float value0 = getParameters(fId).pedestal();
0114   if (eff) {
0115     //account for dark current + crosstalk
0116     auto sipmpar = makeHardcodeSiPMParameter(fId, topo, intlumi);
0117     auto sipmchar = makeHardcodeSiPMCharacteristics();
0118     value0 += sipmpar.getDarkCurrent() * 25. / (1. - sipmchar->getCrossTalk(sipmpar.getType()));
0119   }
0120   float value[4] = {value0, value0, value0, value0};
0121   if (fSmear) {
0122     for (int i = 0; i < 4; i++) {
0123       value[i] = 0.0f;
0124       while (value[i] <= 0.0f)
0125         // ignore correlations, assume 10K pedestal run
0126         value[i] = value0 + (float)CLHEP::RandGauss::shoot(0.0, width.getWidth(i) / 100.);
0127     }
0128   }
0129   HcalPedestal result(fId.rawId(), value[0], value[1], value[2], value[3]);
0130   return result;
0131 }
0132
0133 HcalPedestalWidth HcalDbHardcode::makePedestalWidth(HcalGenericDetId fId,
0134                                                     bool eff,
0135                                                     const HcalTopology* topo,
0136                                                     double intlumi) {
0137   float value = getParameters(fId).pedestalWidth();
0138   float width2 = value * value;
0139   // everything in fC
0140
0141   if (eff) {
0142     //account for dark current + crosstalk
0143     auto sipmpar = makeHardcodeSiPMParameter(fId, topo, intlumi);
0144     auto sipmchar = makeHardcodeSiPMCharacteristics();
0145     //add in quadrature
0146     width2 += sipmpar.getDarkCurrent() * 25. / std::pow(1 - sipmchar->getCrossTalk(sipmpar.getType()), 3) *
0147               sipmpar.getFCByPE();
0148   }
0149
0150   HcalPedestalWidth result(fId.rawId());
0151   for (int i = 0; i < 4; i++) {
0152     for (int j = 0; j < 4; j++) {
0153       result.setSigma(i, j, 0.0);
0154     }
0155     result.setSigma(i, i, width2);
0156   }
0157   return result;
0158 }
0159
0160 HcalGain HcalDbHardcode::makeGain(HcalGenericDetId fId, bool fSmear) const {  // GeV/fC
0161   HcalGainWidth width = makeGainWidth(fId);
0162   float value0 = getParameters(fId).gain(getGainIndex(fId));
0163   float value[4] = {value0, value0, value0, value0};
0164   if (fSmear) {
0165     for (int i = 0; i < 4; i++) {
0166       value[i] = 0.0f;
0167       while (value[i] <= 0.0f)
0168         value[i] = value0 + (float)CLHEP::RandGauss::shoot(0.0, width.getValue(i));
0169     }
0170   }
0171   HcalGain result(fId.rawId(), value[0], value[1], value[2], value[3]);
0172   return result;
0173 }
0174
0175 HcalGainWidth HcalDbHardcode::makeGainWidth(HcalGenericDetId fId) const {  // GeV/fC
0176   float value = getParameters(fId).gainWidth(getGainIndex(fId));
0177   HcalGainWidth result(fId.rawId(), value, value, value, value);
0178   return result;
0179 }
0180
0181 HcalPFCut HcalDbHardcode::makePFCut(HcalGenericDetId fId, double intLumi, bool noHE) const {  // GeV
0182
0183   // assign default dummy parameters
0184   float value0 = getParameters(fId).noiseThreshold();
0185   float value1 = getParameters(fId).seedThreshold();
0186
0187   if (noHE && fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel) {  // HB Phase2
0188
0189     // aging-independent (lumi=0) intialization for Phase2 using "HBphase1" numbers
0190     const float cuts_Phase1[] = {0.1, 0.2, 0.3, 0.3};
0191     const float seeds_Phase1[] = {0.125, 0.25, 0.35, 0.35};
0192
0193     HcalDetId hid(fId);
0194     int depth_m1 = hid.depth() - 1;
0195     value0 = cuts_Phase1[depth_m1];
0196     value1 = seeds_Phase1[depth_m1];
0197
0198     // lumi-dependent stuff for Phase2 from SLHCUpgradeSimulations/Configuration/python/aging.py
0199     const double lumis[] = {300, 1000, 3000, 4500};  // integrated lumi points
0200     // row by row initialization
0201     const float cuts[4][4] = {{0.4, 0.5, 0.6, 0.6}, {0.8, 1.2, 1.2, 1.2}, {1.0, 2.0, 2.0, 2.0}, {1.25, 2.5, 2.5, 2.5}};
0202     const float seeds[4][4] = {
0203         {0.5, 0.625, 0.75, 0.75}, {1.0, 1.5, 1.5, 1.5}, {1.25, 2.5, 2.5, 2.5}, {1.5, 3.0, 3.0, 3.0}};
0204     const double eps = 1.e-6;
0205
0206     for (int i = 0; i < 4; i++) {
0207       if (std::abs(intLumi - lumis[i]) < eps) {
0208         value0 = cuts[i][depth_m1];
0209         value1 = seeds[i][depth_m1];
0210       }
0211     }
0212   }
0213
0214   HcalPFCut result(fId.rawId(), value0, value1);
0215   return result;
0216 }
0217
0218 HcalZSThreshold HcalDbHardcode::makeZSThreshold(HcalGenericDetId fId) const {
0219   int value = getParameters(fId).zsThreshold();
0220   HcalZSThreshold result(fId.rawId(), value);
0221   return result;
0222 }
0223
0224 HcalQIECoder HcalDbHardcode::makeQIECoder(HcalGenericDetId fId) const {
0225   HcalQIECoder result(fId.rawId());
0226   // slope in ADC/fC
0227   const HcalHardcodeParameters& param(getParameters(fId));
0228   for (unsigned range = 0; range < 4; range++) {
0229     for (unsigned capid = 0; capid < 4; capid++) {
0230       result.setOffset(capid, range, param.qieOffset(range));
0231       result.setSlope(capid, range, param.qieSlope(range));
0232     }
0233   }
0234
0235   return result;
0236 }
0237
0238 HcalQIEType HcalDbHardcode::makeQIEType(HcalGenericDetId fId) const {
0239   HcalQIENum qieType = (HcalQIENum)(getParameters(fId).qieType());
0240   HcalQIEType result(fId.rawId(), qieType);
0241   return result;
0242 }
0243
0244 HcalCalibrationQIECoder HcalDbHardcode::makeCalibrationQIECoder(HcalGenericDetId fId) const {
0245   HcalCalibrationQIECoder result(fId.rawId());
0246   float lowEdges[64];
0247   for (int i = 0; i < 64; i++) {
0248     lowEdges[i] = -1.5 + i;
0249   }
0250   result.setMinCharges(lowEdges);
0251   return result;
0252 }
0253
0254 HcalQIEShape HcalDbHardcode::makeQIEShape() const { return HcalQIEShape(); }
0255
0256 HcalMCParam HcalDbHardcode::makeMCParam(HcalGenericDetId fId) const {
0257   int r1bit[5];
0258   r1bit[0] = 9;  //  [0,9]
0259   int syncPhase = 0;
0260   r1bit[1] = 1;
0261   int binOfMaximum = 0;
0262   r1bit[2] = 4;
0263   float phase = -25.0f;                   // [-25.0,25.0]
0264   float Xphase = (phase + 32.0f) * 4.0f;  // never change this line
0265                                           // (offset 50nsec,  0.25ns step)
0266   int Iphase = Xphase;
0267   r1bit[3] = 8;  // [0,256] offset 50ns, .25ns step
0268   int timeSmearing = 0;
0269   r1bit[4] = 1;  //  bool
0270
0271   const HcalHardcodeParameters& hparam(getParameters(fId));
0272   int pulseShapeID = hparam.mcShape();  // a0
0273
0274   if (fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel) {
0275     syncPhase = 1;                    // a1  bool
0276     binOfMaximum = 5;                 // a2
0277     phase = 5.0f;                     // a3  [-25.0,25.0]
0278     Xphase = (phase + 32.0f) * 4.0f;  // never change this line
0279                                       // (offset 50nsec,  0.25ns step)
0280     Iphase = Xphase;
0281     timeSmearing = 1;  // a4
0282
0283   }
0284
0285   else if (fId.genericSubdet() == HcalGenericDetId::HcalGenEndcap) {
0286     syncPhase = 1;                    // a1  bool
0287     binOfMaximum = 5;                 // a2
0288     phase = 5.0f;                     // a3  [-25.0,25.0]
0289     Xphase = (phase + 32.0f) * 4.0f;  // never change this line
0290                                       // (offset 50nsec,  0.25ns step)
0291     Iphase = Xphase;
0292     timeSmearing = 1;  // a4
0293
0294   }
0295
0296   else if (fId.genericSubdet() == HcalGenericDetId::HcalGenOuter) {
0297     syncPhase = 1;                    // a1  bool
0298     binOfMaximum = 5;                 // a2
0299     phase = 5.0f;                     // a3  [-25.0,25.0]
0300     Xphase = (phase + 32.0f) * 4.0f;  // never change this line
0301                                       // (offset 50nsec,  0.25ns step)
0302     Iphase = Xphase;
0303     timeSmearing = 0;  // a4
0304
0305     HcalDetId cell = HcalDetId(fId);
0306     if (cell.ieta() == 1 && cell.iphi() == 1)
0307       pulseShapeID = 125;
0308
0309   }
0310
0311   else if (fId.genericSubdet() == HcalGenericDetId::HcalGenForward) {
0312     syncPhase = 1;                    // a1  bool
0313     binOfMaximum = 3;                 // a2
0314     phase = 14.0f;                    // a3  [-25.0,25.0]
0315     Xphase = (phase + 32.0f) * 4.0f;  // never change this line
0316                                       // (offset 50nsec,  0.25ns step)
0317     Iphase = Xphase;
0318     timeSmearing = 0;  // a4
0319
0320   }
0321
0322   else if (fId.genericSubdet() == HcalGenericDetId::HcalGenZDC) {
0323     pulseShapeID = 401;               // a0
0324     syncPhase = 1;                    // a1  bool
0325     binOfMaximum = 5;                 // a2
0326     phase = -4.0f;                    // a3  [-25.0,25.0]
0327     Xphase = (phase + 32.0f) * 4.0f;  // never change this line
0328                                       // (offset 50nsec,  0.25ns step)
0329     Iphase = Xphase;
0330     timeSmearing = 0;  // a4
0331   }
0332
0333   int rshift[7];
0334   rshift[0] = 0;
0335   for (int k = 0; k < 5; k++) {
0336     rshift[k + 1] = rshift[k] + r1bit[k];
0337   }
0338
0339   int packingScheme = 1;  // a5
0340   unsigned int param = pulseShapeID | syncPhase << rshift[1] | (binOfMaximum << rshift[2]) | (Iphase << rshift[3]) |
0341                        (timeSmearing << rshift[4] | packingScheme << 27);
0342
0343   HcalMCParam result(fId.rawId(), param);
0344   return result;
0345 }
0346
0347 HcalRecoParam HcalDbHardcode::makeRecoParam(HcalGenericDetId fId) const {
0348   // Mostly comes from S.Kunori's macro
0349   int p1bit[6];
0350
0351   // param1
0352   int containmentCorrectionFlag = 0;
0353   p1bit[0] = 1;  // bool
0354   int containmentCorrectionPreSample = 0;
0355   p1bit[1] = 1;                         // bool
0356   float phase = 13.0;                   // [-25.0,25.0]
0357   float Xphase = (phase + 32.0) * 4.0;  //never change this line
0358                                         // (offset 50nsec,  0.25ns step)
0359   int Iphase = Xphase;
0360   p1bit[2] = 8;  // [0,256]
0361                  // (offset 50ns, 0.25ns step
0362   int firstSample = 4;
0363   p1bit[3] = 4;  // [0,9]
0364   int samplesToAdd = 2;
0365   p1bit[4] = 4;  // [0,9]
0366   p1bit[5] = 9;  // [0,9]
0367
0368   const HcalHardcodeParameters& hparam(getParameters(fId));
0369   int pulseShapeID = hparam.recoShape();  // a5
0370
0371   int q2bit[10];
0372   //  param2.
0373   int useLeakCorrection = 0;
0374   q2bit[0] = 1;  // bool
0375   int LeakCorrectionID = 0;
0376   q2bit[1] = 4;  // [0,15]
0377   int correctForTimeslew = 0;
0378   q2bit[2] = 1;  // bool
0379   int timeCorrectionID = 0;
0380   q2bit[3] = 4;  // [0,15]
0381   int correctTiming = 0;
0382   q2bit[4] = 1;  // bool
0383   int firstAuxTS = 0;
0384   q2bit[5] = 4;  // [0,15]
0385   int specialCaseID = 0;
0386   q2bit[6] = 4;  // [0,15]
0387   int noiseFlaggingID = 0;
0388   q2bit[7] = 4;  // [0,15]
0389   int pileupCleaningID = 0;
0390   q2bit[8] = 4;  // [0,15]
0391   int packingScheme = 1;
0392   q2bit[9] = 4;
0393
0394   if ((fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel) ||
0395       (fId.genericSubdet() == HcalGenericDetId::HcalGenEndcap)) {
0396     //  param1.
0397     containmentCorrectionFlag = 1;       // p0
0398     containmentCorrectionPreSample = 0;  // p1
0399     float phase = 6.0;
0400     float Xphase = (phase + 32.0) * 4.0;  // never change this line
0401                                           //(offset 50nsec, 0.25ns step)
0402     Iphase = Xphase;                      // p2
0403     firstSample = 4;                      // p3
0404     samplesToAdd = 2;                     // p4
0405
0406     //  param2.
0407     useLeakCorrection = 0;   //  q0
0408     LeakCorrectionID = 0;    //  q1
0409     correctForTimeslew = 1;  //  q2
0410     timeCorrectionID = 0;    //  q3
0411     correctTiming = 1;       //  q4
0412     firstAuxTS = 4;          //  q5
0413     specialCaseID = 0;       //  q6
0414     noiseFlaggingID = 1;     //  q7
0415     pileupCleaningID = 0;    //  q8
0416   }
0417
0418   else if (fId.genericSubdet() == HcalGenericDetId::HcalGenOuter) {
0419     //  param1.
0420     containmentCorrectionFlag = 1;       // p0
0421     containmentCorrectionPreSample = 0;  // p1
0422     float phase = 13.0;
0423     float Xphase = (phase + 32.0) * 4.0;  // never change this line
0424                                           // (offset 50nsec,  0.25ns step)
0425     Iphase = Xphase;                      // p2
0426     firstSample = 4;                      // p3
0427     samplesToAdd = 4;                     // p4
0428
0429     //  param2.
0430     useLeakCorrection = 0;   //  q0
0431     LeakCorrectionID = 0;    //  q1
0432     correctForTimeslew = 1;  //  q2
0433     timeCorrectionID = 0;    //  q3
0434     correctTiming = 1;       //  q4
0435     firstAuxTS = 4;          //  q5
0436     specialCaseID = 0;       //  q6
0437     noiseFlaggingID = 1;     //  q7
0438     pileupCleaningID = 0;    //  q8
0439
0440   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenForward) {
0441     //  param1.
0442     containmentCorrectionFlag = 0;       // p0
0443     containmentCorrectionPreSample = 0;  // p1
0444     float phase = 13.0;
0445     float Xphase = (phase + 32.0) * 4.0;  // never change this line
0446                                           // (offset 50nsec,  0.25ns step)
0447     Iphase = Xphase;                      // p2
0448     firstSample = 2;                      // p3
0449     samplesToAdd = 1;                     // p4
0450
0451     //  param2.
0452     useLeakCorrection = 0;   //  q0
0453     LeakCorrectionID = 0;    //  q1
0454     correctForTimeslew = 0;  //  q2
0455     timeCorrectionID = 0;    //  q3
0456     correctTiming = 1;       //  q4
0457     firstAuxTS = 1;          //  q5
0458     specialCaseID = 0;       //  q6
0459     noiseFlaggingID = 1;     //  q7
0460     pileupCleaningID = 0;    //  q8
0461   }
0462
0463   // Packing parameters in two words
0464
0465   int p1shift[7];
0466   p1shift[0] = 0;
0467   for (int k = 0; k < 6; k++) {
0468     int j = k + 1;
0469     p1shift[j] = p1shift[k] + p1bit[k];
0470     //     cout<<"  j= "<<j<<"  shift "<< p1shift[j]<<endl;
0471   }
0472   int param1 = 0;
0473   param1 = containmentCorrectionFlag | (containmentCorrectionPreSample << p1shift[1]) | (Iphase << p1shift[2]) |
0474            (firstSample << p1shift[3]) | (samplesToAdd << p1shift[4]) | (pulseShapeID << p1shift[5]);
0475
0476   int q2shift[10];
0477   q2shift[0] = 0;
0478   for (int k = 0; k < 9; k++) {
0479     int j = k + 1;
0480     q2shift[j] = q2shift[k] + q2bit[k];
0481     //  cout<<"  j= "<<j<<"  shift "<< q2shift[j]<<endl;
0482   }
0483   int param2 = 0;
0484   param2 = useLeakCorrection | (LeakCorrectionID << q2shift[1]) | (correctForTimeslew << q2shift[2]) |
0485            (timeCorrectionID << q2shift[3]) | (correctTiming << q2shift[4]) | (firstAuxTS << q2shift[5]) |
0486            (specialCaseID << q2shift[6]) | (noiseFlaggingID << q2shift[7]) | (pileupCleaningID << q2shift[8]) |
0487            (packingScheme << q2shift[9]);
0488
0489   HcalRecoParam result(fId.rawId(), param1, param2);
0490
0491   return result;
0492 }
0493
0494 HcalTimingParam HcalDbHardcode::makeTimingParam(HcalGenericDetId fId) const {
0495   int nhits = 0;
0496   float phase = 0.0;
0497   float rms = 0.0;
0498   if (fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel) {
0499     nhits = 4;
0500     phase = 4.5;
0501     rms = 6.5;
0502   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenEndcap) {
0503     nhits = 4;
0504     phase = 9.3;
0505     rms = 7.8;
0506   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenOuter) {
0507     nhits = 4;
0508     phase = 8.6;
0509     rms = 2.3;
0510   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenForward) {
0511     nhits = 4;
0512     phase = 12.4;
0513     rms = 12.29;
0514   }
0515   HcalTimingParam result(fId.rawId(), nhits, phase, rms);
0516
0517   return result;
0518 }
0519
0520 #define EMAP_NHBHECR 9
0521 #define EMAP_NHFCR 3
0522 #define EMAP_NHOCR 4
0523 #define EMAP_NFBR 8
0524 #define EMAP_NFCH 3
0525 #define EMAP_NHTRS 3
0526 #define EMAP_NHSETS 4
0527 #define EMAP_NTOPBOT 2
0528 #define EMAP_NHTRSHO 4
0529 #define EMAP_NHSETSHO 3
0530
0531 std::unique_ptr<HcalDcsMap> HcalDbHardcode::makeHardcodeDcsMap() const {
0532   HcalDcsMapAddons::Helper dcs_map_helper;
0533   dcs_map_helper.mapGeomId2DcsId(HcalDetId(HcalBarrel, -16, 1, 1),
0534                                  HcalDcsDetId(HcalDcsBarrel, -1, 1, HcalDcsDetId::HV, 2));
0535   dcs_map_helper.mapGeomId2DcsId(HcalDetId(HcalForward, -41, 3, 1),
0536                                  HcalDcsDetId(HcalDcsForward, -1, 1, HcalDcsDetId::DYN8, 1));
0537   dcs_map_helper.mapGeomId2DcsId(HcalDetId(HcalForward, -26, 25, 2),
0538                                  HcalDcsDetId(HcalDcsForward, -1, 7, HcalDcsDetId::HV, 1));
0539   dcs_map_helper.mapGeomId2DcsId(HcalDetId(HcalBarrel, -15, 68, 1),
0540                                  HcalDcsDetId(HcalDcsBarrel, -1, 18, HcalDcsDetId::HV, 3));
0541   dcs_map_helper.mapGeomId2DcsId(HcalDetId(HcalOuter, -14, 1, 4),
0542                                  HcalDcsDetId(HcalDcsOuter, -2, 2, HcalDcsDetId::HV, 4));
0543   dcs_map_helper.mapGeomId2DcsId(HcalDetId(HcalForward, 41, 71, 2),
0544                                  HcalDcsDetId(HcalDcsForward, 1, 4, HcalDcsDetId::DYN8, 3));
0545   return std::make_unique<HcalDcsMap>(dcs_map_helper);
0546 }
0547
0548 std::unique_ptr<HcalElectronicsMap> HcalDbHardcode::makeHardcodeMap(const std::vector<HcalGenericDetId>& cells) const {
0549   static const int kUTCAMask = 0x4000000;         //set bit 26 for uTCA version
0550   static const int kLinearIndexMax = 0x7FFFF;     //19 bits
0551   static const int kTriggerBitMask = 0x02000000;  //2^25
0552   uint32_t counter = 0;
0553   uint32_t counterTrig = 0;
0554   HcalElectronicsMapAddons::Helper emapHelper;
0555   for (const auto& fId : cells) {
0556     if (fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel ||
0557         fId.genericSubdet() == HcalGenericDetId::HcalGenEndcap ||
0558         fId.genericSubdet() == HcalGenericDetId::HcalGenForward ||
0559         fId.genericSubdet() == HcalGenericDetId::HcalGenOuter || fId.genericSubdet() == HcalGenericDetId::HcalGenZDC) {
0560       ++counter;
0561       assert(counter < kLinearIndexMax);
0562       uint32_t raw = counter;
0563       raw |= kUTCAMask;
0564       HcalElectronicsId elId(raw);
0565       emapHelper.mapEId2chId(elId, fId);
0566     } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenTriggerTower) {
0567       ++counterTrig;
0568       assert(counterTrig < kLinearIndexMax);
0569       uint32_t raw = counterTrig;
0570       raw |= kUTCAMask | kTriggerBitMask;
0571       HcalElectronicsId elId(raw);
0572       emapHelper.mapEId2tId(elId, fId);
0573     }
0574   }
0575   return std::make_unique<HcalElectronicsMap>(emapHelper);
0576 }
0577
0578 std::unique_ptr<HcalFrontEndMap> HcalDbHardcode::makeHardcodeFrontEndMap(
0579     const std::vector<HcalGenericDetId>& cells) const {
0580   HcalFrontEndMapAddons::Helper emapHelper;
0581   std::stringstream mystream;
0582   std::string detector[5] = {"XX", "HB", "HE", "HO", "HF"};
0583   for (const auto& fId : cells) {
0584     if (fId.isHcalDetId()) {
0585       HcalDetId id = HcalDetId(fId.rawId());
0586       HcalSubdetector subdet = id.subdet();
0587       int ieta = id.ietaAbs();
0588       int iside = id.zside();
0589       int iphi = id.iphi();
0590       std::string det, rbx;
0591       int irm(0);
0592       char tempbuff[30];
0593       char sidesign = (iside == -1) ? 'M' : 'P';
0594       if (subdet == HcalBarrel || subdet == HcalEndcap) {
0595         det = detector[subdet];
0596         irm = (iphi + 1) % 4 + 1;
0597         int iwedge(0);
0598         if (ieta >= 21 && (irm == 1 || irm == 3))
0599           iwedge = (iphi + 1 + irm + 1) / 4;
0600         else
0601           iwedge = (iphi + irm + 1) / 4;
0602         if (iwedge > 18)
0603           iwedge -= 18;
0604         sprintf(tempbuff, "%s%c%2.2i%c", det.c_str(), sidesign, iwedge, '\0');
0605         mystream << tempbuff;
0606         rbx = mystream.str();
0607         mystream.str("");
0608         emapHelper.loadObject(id, irm, rbx);
0609       } else if (subdet == HcalForward) {
0610         det = detector[subdet];
0611         int hfphi(0);
0612         if ((iside == 1 && ieta == 40) || (iside == -1 && ieta == 41)) {
0613           irm = ((iphi + 1) / 2) % 36 + 1;
0614           hfphi = ((iphi + 1) / 6) % 12 + 1;
0615         } else {
0616           irm = (iphi + 1) / 2;
0617           hfphi = (iphi - 1) / 6 + 1;
0618         }
0619         irm = (irm - 1) % 3 + 1;
0620         sprintf(tempbuff, "%s%c%2.2i%c", det.c_str(), sidesign, hfphi, '\0');
0621         mystream << tempbuff;
0622         rbx = mystream.str();
0623         mystream.str("");
0624         emapHelper.loadObject(id, irm, rbx);
0625       } else if (subdet == HcalOuter) {
0626         det = detector[subdet];
0627         int ring(0), sector(0);
0628         if (ieta <= 4)
0629           ring = 0;
0630         else if (ieta >= 5 && ieta <= 10)
0631           ring = 1;
0632         else
0633           ring = 2;
0634         for (int i = -2; i < iphi; i += 6)
0635           sector++;
0636         if (sector > 12)
0637           sector = 1;
0638         irm = ((iphi + 1) / 2) % 6 + 1;
0639         if (ring != 0 && sector % 2 != 0)
0640           sector++;
0641         if (ring == 0)
0642           sprintf(tempbuff, "%s%i%2.2d", det.c_str(), ring, sector);
0643         else
0644           sprintf(tempbuff, "%s%i%c%2.2d", det.c_str(), ring, sidesign, sector);
0645         mystream << tempbuff;
0646         rbx = mystream.str();
0647         mystream.str("");
0648         emapHelper.loadObject(id, irm, rbx);
0649       }
0650     }
0651   }
0652   return std::make_unique<HcalFrontEndMap>(emapHelper);
0653 }
0654
0655 int HcalDbHardcode::getLayersInDepth(int ieta, int depth, const HcalTopology* topo) {
0656   //check for cached value
0657   auto eta_depth_pair = std::make_pair(ieta, depth);
0658   auto nLayers = theLayersInDepths_.find(eta_depth_pair);
0659   if (nLayers != theLayersInDepths_.end()) {
0660     return nLayers->second;
0661   } else {
0662     std::vector<int> segmentation;
0663     topo->getDepthSegmentation(ieta, segmentation);
0664     //assume depth segmentation vector is sorted
0665     int nLayersInDepth = std::distance(std::lower_bound(segmentation.begin(), segmentation.end(), depth),
0666                                        std::upper_bound(segmentation.begin(), segmentation.end(), depth));
0667     theLayersInDepths_.insert(std::make_pair(eta_depth_pair, nLayersInDepth));
0668     return nLayersInDepth;
0669   }
0670 }
0671
0672 bool HcalDbHardcode::isHEPlan1(HcalGenericDetId fId) const {
0673   if (fId.isHcalDetId()) {
0674     HcalDetId hid(fId);
0675     //special mixed case for HE 2017
0676     if (hid.zside() == 1 && (hid.iphi() == 63 || hid.iphi() == 64 || hid.iphi() == 65 || hid.iphi() == 66))
0677       return true;
0678   }
0679   return false;
0680 }
0681
0682 HcalSiPMParameter HcalDbHardcode::makeHardcodeSiPMParameter(HcalGenericDetId fId,
0683                                                             const HcalTopology* topo,
0684                                                             double intlumi) {
0685   // SiPMParameter defined for each DetId the following quantities:
0686   //  SiPM type, PhotoElectronToAnalog, Dark Current, two auxiliary words
0687   //  (the second of those containing float noise correlation coefficient
0688   //  These numbers come from some measurements done with SiPMs
0689   // rule for type: cells with >4 layers use larger device (3.3mm diameter), otherwise 2.8mm
0690   HcalSiPMType theType = HcalNoSiPM;
0691   double thePe2fC = getParameters(fId).photoelectronsToAnalog();
0692   double theDC = getParameters(fId).darkCurrent(0, intlumi);
0693   double theNoiseCN = getParameters(fId).noiseCorrelation(0);
0694   if (fId.genericSubdet() == HcalGenericDetId::HcalGenBarrel) {
0695     if (useHBUpgrade_) {
0696       HcalDetId hid(fId);
0697       int nLayersInDepth = getLayersInDepth(hid.ietaAbs(), hid.depth(), topo);
0698       if (nLayersInDepth > 4) {
0699         theType = HcalHBHamamatsu2;
0700         theDC = getParameters(fId).darkCurrent(1, intlumi);
0701         theNoiseCN = getParameters(fId).noiseCorrelation(1);
0702       } else {
0703         theType = HcalHBHamamatsu1;
0704         theDC = getParameters(fId).darkCurrent(0, intlumi);
0705         theNoiseCN = getParameters(fId).noiseCorrelation(0);
0706       }
0707     } else
0708       theType = HcalHPD;
0709   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenEndcap) {
0710     if (useHEUpgrade_ || (testHEPlan1_ && isHEPlan1(fId))) {
0711       HcalDetId hid(fId);
0712       int nLayersInDepth = getLayersInDepth(hid.ietaAbs(), hid.depth(), topo);
0713       if (nLayersInDepth > 4) {
0714         theType = HcalHEHamamatsu2;
0715         theDC = getParameters(fId).darkCurrent(1, intlumi);
0716         theNoiseCN = getParameters(fId).noiseCorrelation(1);
0717       } else {
0718         theType = HcalHEHamamatsu1;
0719         theDC = getParameters(fId).darkCurrent(0, intlumi);
0720         theNoiseCN = getParameters(fId).noiseCorrelation(0);
0721       }
0722     } else
0723       theType = HcalHPD;
0724   } else if (fId.genericSubdet() == HcalGenericDetId::HcalGenOuter) {
0725     if (useHOUpgrade_)
0726       theType = HcalHOHamamatsu;
0727     else
0728       theType = HcalHPD;
0729   }
0730
0731   return HcalSiPMParameter(fId.rawId(), theType, thePe2fC, theDC, 0, (float)theNoiseCN);
0732 }
0733
0734 std::unique_ptr<HcalSiPMCharacteristics> HcalDbHardcode::makeHardcodeSiPMCharacteristics() const {
0735   // SiPMCharacteristics are constants for each type of SiPM:
0736   // Type, # of pixels, 3 parameters for non-linearity, cross talk parameter, ..
0737   // Obtained from data sheet and measurements
0738   // types (in order): HcalHOZecotek=1, HcalHOHamamatsu, HcalHEHamamatsu1, HcalHEHamamatsu2, HcalHBHamamatsu1, HcalHBHamamatsu2, HcalHPD
0739   HcalSiPMCharacteristicsAddons::Helper sipmHelper;
0740   for (unsigned ip = 0; ip < theSiPMCharacteristics_.size(); ++ip) {
0741     auto& ps = theSiPMCharacteristics_[ip];
0742     sipmHelper.loadObject(ip + 1,
0743                           ps.getParameter<int>("pixels"),
0744                           ps.getParameter<double>("nonlin1"),
0745                           ps.getParameter<double>("nonlin2"),
0746                           ps.getParameter<double>("nonlin3"),
0747                           ps.getParameter<double>("crosstalk"),
0748                           0,
0749                           0);
0750   }
0751   return std::make_unique<HcalSiPMCharacteristics>(sipmHelper);
0752 }
0753
0754 HcalTPChannelParameter HcalDbHardcode::makeHardcodeTPChannelParameter(HcalGenericDetId fId) const {
0755   // For each detId parameters for trigger primitive
0756   // mask for channel validity and self trigger information, fine grain
0757   // bit information and auxiliary words
0758   uint32_t bitInfo = ((44 << 16) | 30);
0759   return HcalTPChannelParameter(fId.rawId(), 0, bitInfo, 0, 0);
0760 }
0761
0762 void HcalDbHardcode::makeHardcodeTPParameters(HcalTPParameters& tppar) const {
0763   // Parameters for a given TP algorithm:
0764   // FineGrain Algorithm Version for HBHE, ADC threshold fof TDC mask of HF,
0765   // TDC mask for HF, Self Trigger bits, auxiliary words
0766   tppar.loadObject(0, 0, 0xFFFFFFFFFFFFFFFF, 0, 0, 0);
0767 }