GlobalCaloTrigger/test/produceTrivialCalibrationLut.cc

0001 #include "L1Trigger/GlobalCaloTrigger/test/produceTrivialCalibrationLut.h"
0002
0003 #include "CondFormats/L1TObjects/interface/L1GctJetFinderParams.h"
0004
0005 #include "CondFormats/L1TObjects/interface/L1CaloEtScale.h"
0006
0007 #include <math.h>
0008 #include <cassert>
0009
0010 produceTrivialCalibrationLut::produceTrivialCalibrationLut()
0011     : m_htScaleLSB(0.5),
0012       m_jetCalibFunc(L1GctJetFinderParams::NUMBER_ETA_VALUES),
0013       m_tauCalibFunc(L1GctJetFinderParams::N_CENTRAL_ETA_VALUES),
0014       m_jetEtScaleInputLsb(0.5),
0015       m_corrFunType(1),
0016       m_jfPars(new L1GctJetFinderParams()) {
0017   const double jetThresholds[64] = {0.,   10.,  12.,  14.,  15.,  18.,  20.,  22.,  24.,  25.,  28.,  30.,  32.,
0018                                     35.,  37.,  40.,  45.,  50.,  55.,  60.,  65.,  70.,  75.,  80.,  85.,  90.,
0019                                     100., 110., 120., 125., 130., 140., 150., 160., 170., 175., 180., 190., 200.,
0020                                     215., 225., 235., 250., 275., 300., 325., 350., 375., 400., 425., 450., 475.,
0021                                     500., 525., 550., 575., 600., 625., 650., 675., 700., 725., 750., 775.};
0022   for (unsigned i = 0; i < 64; i++) {
0023     m_jetEtThresholds.push_back(jetThresholds[i]);
0024   }
0025
0026   setupJfPars();
0027 }
0028
0029 produceTrivialCalibrationLut::~produceTrivialCalibrationLut() {}
0030
0031 void produceTrivialCalibrationLut::setPowerSeriesCorrectionType() {
0032   m_corrFunType = 1;
0033   for (unsigned i = 0; i < m_jetCalibFunc.size(); i++) {
0034     m_jetCalibFunc.at(i).clear();
0035   }
0036   for (unsigned i = 0; i < m_tauCalibFunc.size(); i++) {
0037     m_tauCalibFunc.at(i).clear();
0038   }
0039   setupJfPars();
0040 }
0041
0042 void produceTrivialCalibrationLut::setOrcaStyleCorrectionType() {
0043   const double y[11] = {47.4, 49.4, 47.1, 49.3, 48.2, 42.0, 33.8, 17.1, 13.3, 12.4, 9.3};
0044   const double A[11] = {-20.7, -22.5, -22.2, -22.9, -24.5, -23.9, -22.1, -6.6, -4.5, -3.8, 1.3};
0045   const double B[11] = {0.7922, 0.7867, 0.7645, 0.7331, 0.7706, 0.7945, 0.8202, 0.6958, 0.7071, 0.6558, 0.2719};
0046   const double C[11] = {
0047       9.53E-5, 9.60E-5, 12.09E-5, 12.21E-5, 12.80E-5, 14.58E-5, 14.03E-5, 6.88E-5, 7.26E-5, 48.90E-5, 341.8E-5};
0048
0049   m_corrFunType = 2;
0050   for (unsigned i = 0; i < m_jetCalibFunc.size(); i++) {
0051     m_jetCalibFunc.at(i).clear();
0052     m_jetCalibFunc.at(i).push_back(y[i]);
0053     m_jetCalibFunc.at(i).push_back(A[i]);
0054     m_jetCalibFunc.at(i).push_back(B[i]);
0055     m_jetCalibFunc.at(i).push_back(C[i]);
0056   }
0057   for (unsigned i = 0; i < m_tauCalibFunc.size(); i++) {
0058     m_tauCalibFunc.at(i).clear();
0059     m_tauCalibFunc.at(i).push_back(y[i]);
0060     m_tauCalibFunc.at(i).push_back(A[i]);
0061     m_tauCalibFunc.at(i).push_back(B[i]);
0062     m_tauCalibFunc.at(i).push_back(C[i]);
0063   }
0064   setOrcaStyleParams();
0065   setupJfPars();
0066 }
0067
0068 const produceTrivialCalibrationLut::lutPtrVector produceTrivialCalibrationLut::produce() const {
0069   L1CaloEtScale* jetScale = new L1CaloEtScale(m_jetEtScaleInputLsb, m_jetEtThresholds);
0070
0071   lutPtrVector lutVector;
0072   lutPtr nextLut(new L1GctJetEtCalibrationLut());
0073
0074   for (unsigned ieta = 0; ieta < L1GctJetFinderParams::NUMBER_ETA_VALUES; ieta++) {
0075     nextLut->setEtaBin(ieta);
0076     nextLut->setFunction(m_jfPars);
0077     nextLut->setOutputEtScale(jetScale);
0078     lutVector.push_back(nextLut);
0079     nextLut.reset(new L1GctJetEtCalibrationLut());
0080   }
0081
0082   return lutVector;
0083 }
0084
0085 void produceTrivialCalibrationLut::setupJfPars() {
0086   m_jfPars->setRegionEtLsb(m_jetEtScaleInputLsb);
0087
0088   m_jfPars->setJetEtCalibrationParams(m_corrFunType, m_jetCalibFunc, m_tauCalibFunc);
0089
0090   m_jfPars->setHtSumParams(m_htScaleLSB, 0.0, 0.0);
0091 }
0092
0093 //--------------------------------------------------------------------------
0094 //
0095 // For ORCA-style calibration, we extend the calibration function downwards
0096 // in energy in an automated way here.
0097 void produceTrivialCalibrationLut::setOrcaStyleParams() {
0098   for (unsigned i = 0; i < m_jetCalibFunc.size(); ++i) {
0099     setOrcaStyleParamsForBin(m_jetCalibFunc.at(i));
0100   }
0101   for (unsigned i = 0; i < m_tauCalibFunc.size(); ++i) {
0102     setOrcaStyleParamsForBin(m_tauCalibFunc.at(i));
0103   }
0104 }
0105
0106 // The ORCA-style calibration function is a series of inverted quadratic functions
0107 // (ie x = A + B.y + C.y^2)
0108 //
0109 // This assumes that just one high-energy set of parameters is supplied, together with
0110 // an energy threshold value (in terms of y, the "true" jet Et). It calculates a set of
0111 // parameters to be applied to lower-Et jets subject to the following constraints:
0112 //  (i) The calibration function is continuous at the high threshold value;
0113 //  (ii) Its slope is also continuous at the high threshold value;
0114 //  (iii) At the low (zero-suppression) threshold, the calibration function returns y=x.
0115 //
0116 void produceTrivialCalibrationLut::setOrcaStyleParamsForBin(std::vector<double>& paramsForBin) {
0117   assert(paramsForBin.size() == 4);
0118
0119   double x2 = paramsForBin.at(0);
0120   double A = paramsForBin.at(1);
0121   double B = paramsForBin.at(2);
0122   double C = paramsForBin.at(3);
0123
0124   double y2 = 2 * (x2 - A) / (B + sqrt(B * B - 4 * (A - x2) * C));
0125
0126   double y1 = 5.0;  // This was zero suppresion threshold parameter, but such a parameter is no longer defined
0127   double g = (y1 - (A + y1 * (B + y1 * C))) / (pow((y2 - y1), 2));
0128   A = A + g * y2 * y2;
0129   B = B - 2.0 * g * y2;
0130   C = C + g;
0131
0132   paramsForBin.push_back(y1);
0133   paramsForBin.push_back(A);
0134   paramsForBin.push_back(B);
0135   paramsForBin.push_back(C);
0136 }
0137 //
0138 //--------------------------------------------------------------------------