Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-02-12 09:07:52

0001 #include "RecoLocalTracker/SiPixelRecHits/interface/PixelCPEGeneric.h"
0002 
0003 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0004 #include "Geometry/TrackerGeometryBuilder/interface/RectangularPixelTopology.h"
0005 #include "DataFormats/DetId/interface/DetId.h"
0006 
0007 // Pixel templates contain the rec hit error parameterizaiton
0008 #include "CondFormats/SiPixelTransient/interface/SiPixelTemplate.h"
0009 
0010 // The generic formula
0011 #include "CondFormats/SiPixelTransient/interface/SiPixelUtils.h"
0012 
0013 // Services
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "MagneticField/Engine/interface/MagneticField.h"
0016 
0017 #include "boost/multi_array.hpp"
0018 
0019 #include <iostream>
0020 using namespace std;
0021 
0022 namespace {
0023   constexpr float micronsToCm = 1.0e-4;
0024 }  // namespace
0025 
0026 //-----------------------------------------------------------------------------
0027 //!  The constructor.
0028 //-----------------------------------------------------------------------------
0029 PixelCPEGeneric::PixelCPEGeneric(edm::ParameterSet const& conf,
0030                                  const MagneticField* mag,
0031                                  const TrackerGeometry& geom,
0032                                  const TrackerTopology& ttopo,
0033                                  const SiPixelLorentzAngle* lorentzAngle,
0034                                  const SiPixelGenErrorDBObject* genErrorDBObject,
0035                                  const SiPixelLorentzAngle* lorentzAngleWidth = nullptr)
0036     : PixelCPEGenericBase(conf, mag, geom, ttopo, lorentzAngle, genErrorDBObject, lorentzAngleWidth) {
0037   if (theVerboseLevel > 0)
0038     LogDebug("PixelCPEGeneric") << " constructing a generic algorithm for ideal pixel detector.\n"
0039                                 << " CPEGeneric:: VerboseLevel = " << theVerboseLevel;
0040 
0041   // Externally settable cuts
0042   the_eff_charge_cut_lowX = conf.getParameter<double>("eff_charge_cut_lowX");
0043   the_eff_charge_cut_lowY = conf.getParameter<double>("eff_charge_cut_lowY");
0044   the_eff_charge_cut_highX = conf.getParameter<double>("eff_charge_cut_highX");
0045   the_eff_charge_cut_highY = conf.getParameter<double>("eff_charge_cut_highY");
0046   the_size_cutX = conf.getParameter<double>("size_cutX");
0047   the_size_cutY = conf.getParameter<double>("size_cutY");
0048 
0049   // Externally settable flags to inflate errors
0050   inflate_errors = conf.getParameter<bool>("inflate_errors");
0051   inflate_all_errors_no_trk_angle = conf.getParameter<bool>("inflate_all_errors_no_trk_angle");
0052 
0053   NoTemplateErrorsWhenNoTrkAngles_ = conf.getParameter<bool>("NoTemplateErrorsWhenNoTrkAngles");
0054   IrradiationBiasCorrection_ = conf.getParameter<bool>("IrradiationBiasCorrection");
0055   DoCosmics_ = conf.getParameter<bool>("DoCosmics");
0056 
0057   // Upgrade means phase 2
0058   isPhase2_ = conf.getParameter<bool>("Upgrade");
0059 
0060   // For cosmics force the use of simple errors
0061   if ((DoCosmics_))
0062     useErrorsFromTemplates_ = false;
0063 
0064   if (!useErrorsFromTemplates_ && (truncatePixelCharge_ || IrradiationBiasCorrection_ || LoadTemplatesFromDB_)) {
0065     throw cms::Exception("PixelCPEGeneric::PixelCPEGeneric: ")
0066         << "\nERROR: useErrorsFromTemplates_ is set to False in PixelCPEGeneric_cfi.py. "
0067         << " In this case it does not make sense to set any of the following to True: "
0068         << " truncatePixelCharge_, IrradiationBiasCorrection_, DoCosmics_, LoadTemplatesFromDB_ !!!"
0069         << "\n\n";
0070   }
0071 
0072   // Use errors from templates or from GenError
0073   if (useErrorsFromTemplates_) {
0074     if (LoadTemplatesFromDB_) {  // From DB
0075       if (!SiPixelGenError::pushfile(*genErrorDBObject_, thePixelGenError_))
0076         throw cms::Exception("InvalidCalibrationLoaded")
0077             << "ERROR: GenErrors not filled correctly. Check the sqlite file. Using SiPixelTemplateDBObject version "
0078             << (*genErrorDBObject_).version();
0079       LogDebug("PixelCPEGeneric") << "Loaded genErrorDBObject v" << (*genErrorDBObject_).version();
0080     } else {  // From file
0081       if (!SiPixelGenError::pushfile(-999, thePixelGenError_))
0082         throw cms::Exception("InvalidCalibrationLoaded")
0083             << "ERROR: GenErrors not loaded correctly from text file. Reconstruction will fail.";
0084     }  // if load from DB
0085 
0086   } else {
0087 #ifdef EDM_ML_DEBUG
0088     cout << " Use simple parametrised errors " << endl;
0089 #endif
0090   }  // if ( useErrorsFromTemplates_ )
0091 
0092 #ifdef EDM_ML_DEBUG
0093   cout << "From PixelCPEGeneric::PixelCPEGeneric(...)" << endl;
0094   cout << "(int)useErrorsFromTemplates_ = " << (int)useErrorsFromTemplates_ << endl;
0095   cout << "truncatePixelCharge_         = " << (int)truncatePixelCharge_ << endl;
0096   cout << "IrradiationBiasCorrection_   = " << (int)IrradiationBiasCorrection_ << endl;
0097   cout << "(int)DoCosmics_              = " << (int)DoCosmics_ << endl;
0098   cout << "(int)LoadTemplatesFromDB_    = " << (int)LoadTemplatesFromDB_ << endl;
0099 #endif
0100 }
0101 
0102 //-----------------------------------------------------------------------------
0103 //! Hit position in the local frame (in cm).  Unlike other CPE's, this
0104 //! one converts everything from the measurement frame (in channel numbers)
0105 //! into the local frame (in centimeters).
0106 //-----------------------------------------------------------------------------
0107 LocalPoint PixelCPEGeneric::localPosition(DetParam const& theDetParam, ClusterParam& theClusterParamBase) const {
0108   ClusterParamGeneric& theClusterParam = static_cast<ClusterParamGeneric&>(theClusterParamBase);
0109 
0110   //cout<<" in PixelCPEGeneric:localPosition - "<<endl; //dk
0111 
0112   float chargeWidthX = (theDetParam.lorentzShiftInCmX * theDetParam.widthLAFractionX);
0113   float chargeWidthY = (theDetParam.lorentzShiftInCmY * theDetParam.widthLAFractionY);
0114   float shiftX = 0.5f * theDetParam.lorentzShiftInCmX;
0115   float shiftY = 0.5f * theDetParam.lorentzShiftInCmY;
0116 
0117   //cout<<" main la width "<<chargeWidthX<<" "<<chargeWidthY<<endl;
0118 
0119   if (useErrorsFromTemplates_) {
0120     float qclus = theClusterParam.theCluster->charge();
0121     float locBz = theDetParam.bz;
0122     float locBx = theDetParam.bx;
0123     //cout << "PixelCPEGeneric::localPosition(...) : locBz = " << locBz << endl;
0124 
0125     theClusterParam.pixmx = -999;     // max pixel charge for truncation of 2-D cluster
0126     theClusterParam.sigmay = -999.9;  // CPE Generic y-error for multi-pixel cluster
0127     theClusterParam.deltay = -999.9;  // CPE Generic y-bias for multi-pixel cluster
0128     theClusterParam.sigmax = -999.9;  // CPE Generic x-error for multi-pixel cluster
0129     theClusterParam.deltax = -999.9;  // CPE Generic x-bias for multi-pixel cluster
0130     theClusterParam.sy1 = -999.9;     // CPE Generic y-error for single single-pixel
0131     theClusterParam.dy1 = -999.9;     // CPE Generic y-bias for single single-pixel cluster
0132     theClusterParam.sy2 = -999.9;     // CPE Generic y-error for single double-pixel cluster
0133     theClusterParam.dy2 = -999.9;     // CPE Generic y-bias for single double-pixel cluster
0134     theClusterParam.sx1 = -999.9;     // CPE Generic x-error for single single-pixel cluster
0135     theClusterParam.dx1 = -999.9;     // CPE Generic x-bias for single single-pixel cluster
0136     theClusterParam.sx2 = -999.9;     // CPE Generic x-error for single double-pixel cluster
0137     theClusterParam.dx2 = -999.9;     // CPE Generic x-bias for single double-pixel cluster
0138 
0139     SiPixelGenError gtempl(thePixelGenError_);
0140     int gtemplID_ = theDetParam.detTemplateId;
0141 
0142     //int gtemplID0 = genErrorDBObject_->getGenErrorID(theDetParam.theDet->geographicalId().rawId());
0143     //if(gtemplID0!=gtemplID_) cout<<" different id "<< gtemplID_<<" "<<gtemplID0<<endl;
0144 
0145     theClusterParam.qBin_ = gtempl.qbin(gtemplID_,
0146                                         theClusterParam.cotalpha,
0147                                         theClusterParam.cotbeta,
0148                                         locBz,
0149                                         locBx,
0150                                         qclus,
0151                                         IrradiationBiasCorrection_,
0152                                         theClusterParam.pixmx,
0153                                         theClusterParam.sigmay,
0154                                         theClusterParam.deltay,
0155                                         theClusterParam.sigmax,
0156                                         theClusterParam.deltax,
0157                                         theClusterParam.sy1,
0158                                         theClusterParam.dy1,
0159                                         theClusterParam.sy2,
0160                                         theClusterParam.dy2,
0161                                         theClusterParam.sx1,
0162                                         theClusterParam.dx1,
0163                                         theClusterParam.sx2,
0164                                         theClusterParam.dx2);
0165 
0166     // now use the charge widths stored in the new generic template headers (change to the
0167     // incorrect sign convention of the base class)
0168     bool useLAWidthFromGenError = false;
0169     if (useLAWidthFromGenError) {
0170       chargeWidthX = (-micronsToCm * gtempl.lorxwidth());
0171       chargeWidthY = (-micronsToCm * gtempl.lorywidth());
0172       LogDebug("PixelCPE localPosition():") << "redefine la width (gen-error)" << chargeWidthX << chargeWidthY;
0173     }
0174     LogDebug("PixelCPE localPosition():") << "GenError:" << gtemplID_;
0175 
0176     // These numbers come in microns from the qbin(...) call. Transform them to cm.
0177     theClusterParam.deltax = theClusterParam.deltax * micronsToCm;
0178     theClusterParam.dx1 = theClusterParam.dx1 * micronsToCm;
0179     theClusterParam.dx2 = theClusterParam.dx2 * micronsToCm;
0180 
0181     theClusterParam.deltay = theClusterParam.deltay * micronsToCm;
0182     theClusterParam.dy1 = theClusterParam.dy1 * micronsToCm;
0183     theClusterParam.dy2 = theClusterParam.dy2 * micronsToCm;
0184 
0185     theClusterParam.sigmax = theClusterParam.sigmax * micronsToCm;
0186     theClusterParam.sx1 = theClusterParam.sx1 * micronsToCm;
0187     theClusterParam.sx2 = theClusterParam.sx2 * micronsToCm;
0188 
0189     theClusterParam.sigmay = theClusterParam.sigmay * micronsToCm;
0190     theClusterParam.sy1 = theClusterParam.sy1 * micronsToCm;
0191     theClusterParam.sy2 = theClusterParam.sy2 * micronsToCm;
0192 
0193   }  // if ( useErrorsFromTemplates_ )
0194   else {
0195     theClusterParam.qBin_ = 0;
0196   }
0197 
0198   int q_f_X;  //!< Q of the first  pixel  in X
0199   int q_l_X;  //!< Q of the last   pixel  in X
0200   int q_f_Y;  //!< Q of the first  pixel  in Y
0201   int q_l_Y;  //!< Q of the last   pixel  in Y
0202   collect_edge_charges(theClusterParam, q_f_X, q_l_X, q_f_Y, q_l_Y, useErrorsFromTemplates_ && truncatePixelCharge_);
0203 
0204   //--- Find the inner widths along X and Y in one shot.  We
0205   //--- compute the upper right corner of the inner pixels
0206   //--- (== lower left corner of upper right pixel) and
0207   //--- the lower left corner of the inner pixels
0208   //--- (== upper right corner of lower left pixel), and then
0209   //--- subtract these two points in the formula.
0210 
0211   //--- Upper Right corner of Lower Left pixel -- in measurement frame
0212   MeasurementPoint meas_URcorn_LLpix(theClusterParam.theCluster->minPixelRow() + 1.0,
0213                                      theClusterParam.theCluster->minPixelCol() + 1.0);
0214 
0215   //--- Lower Left corner of Upper Right pixel -- in measurement frame
0216   MeasurementPoint meas_LLcorn_URpix(theClusterParam.theCluster->maxPixelRow(),
0217                                      theClusterParam.theCluster->maxPixelCol());
0218 
0219   //--- These two now converted into the local
0220   LocalPoint local_URcorn_LLpix;
0221   LocalPoint local_LLcorn_URpix;
0222 
0223   // PixelCPEGeneric can be used with or without track angles
0224   // If PixelCPEGeneric is called with track angles, use them to correct for bows/kinks:
0225   if (theClusterParam.with_track_angle) {
0226     local_URcorn_LLpix = theDetParam.theTopol->localPosition(meas_URcorn_LLpix, theClusterParam.loc_trk_pred);
0227     local_LLcorn_URpix = theDetParam.theTopol->localPosition(meas_LLcorn_URpix, theClusterParam.loc_trk_pred);
0228   } else {
0229     local_URcorn_LLpix = theDetParam.theTopol->localPosition(meas_URcorn_LLpix);
0230     local_LLcorn_URpix = theDetParam.theTopol->localPosition(meas_LLcorn_URpix);
0231   }
0232 
0233 #ifdef EDM_ML_DEBUG
0234   if (theVerboseLevel > 20) {
0235     cout << "\n\t >>> theClusterParam.theCluster->x = " << theClusterParam.theCluster->x()
0236          << "\n\t >>> theClusterParam.theCluster->y = " << theClusterParam.theCluster->y()
0237          << "\n\t >>> cluster: minRow = " << theClusterParam.theCluster->minPixelRow()
0238          << "  minCol = " << theClusterParam.theCluster->minPixelCol()
0239          << "\n\t >>> cluster: maxRow = " << theClusterParam.theCluster->maxPixelRow()
0240          << "  maxCol = " << theClusterParam.theCluster->maxPixelCol()
0241          << "\n\t >>> meas: inner lower left  = " << meas_URcorn_LLpix.x() << "," << meas_URcorn_LLpix.y()
0242          << "\n\t >>> meas: inner upper right = " << meas_LLcorn_URpix.x() << "," << meas_LLcorn_URpix.y() << endl;
0243   }
0244 #endif
0245 
0246   //--- &&& Note that the cuts below should not be hardcoded (like in Orca and
0247   //--- &&& CPEFromDetPosition/PixelCPEInitial), but rather be
0248   //--- &&& externally settable (but tracked) parameters.
0249 
0250   //--- Position, including the half lorentz shift
0251 
0252 #ifdef EDM_ML_DEBUG
0253   if (theVerboseLevel > 20)
0254     cout << "\t >>> Generic:: processing X" << endl;
0255 #endif
0256 
0257   float xPos = SiPixelUtils::generic_position_formula(
0258       theClusterParam.theCluster->sizeX(),
0259       q_f_X,
0260       q_l_X,
0261       local_URcorn_LLpix.x(),
0262       local_LLcorn_URpix.x(),
0263       chargeWidthX,  // lorentz shift in cm
0264       theDetParam.theThickness,
0265       theClusterParam.cotalpha,
0266       theDetParam.thePitchX,
0267       theDetParam.theRecTopol->isItBigPixelInX(theClusterParam.theCluster->minPixelRow()),
0268       theDetParam.theRecTopol->isItBigPixelInX(theClusterParam.theCluster->maxPixelRow()),
0269       the_eff_charge_cut_lowX,
0270       the_eff_charge_cut_highX,
0271       the_size_cutX);  // cut for eff charge width &&&
0272 
0273   // apply the lorentz offset correction
0274   xPos = xPos + shiftX;
0275 
0276 #ifdef EDM_ML_DEBUG
0277   if (theVerboseLevel > 20)
0278     cout << "\t >>> Generic:: processing Y" << endl;
0279 #endif
0280 
0281   float yPos = SiPixelUtils::generic_position_formula(
0282       theClusterParam.theCluster->sizeY(),
0283       q_f_Y,
0284       q_l_Y,
0285       local_URcorn_LLpix.y(),
0286       local_LLcorn_URpix.y(),
0287       chargeWidthY,  // lorentz shift in cm
0288       theDetParam.theThickness,
0289       theClusterParam.cotbeta,
0290       theDetParam.thePitchY,
0291       theDetParam.theRecTopol->isItBigPixelInY(theClusterParam.theCluster->minPixelCol()),
0292       theDetParam.theRecTopol->isItBigPixelInY(theClusterParam.theCluster->maxPixelCol()),
0293       the_eff_charge_cut_lowY,
0294       the_eff_charge_cut_highY,
0295       the_size_cutY);  // cut for eff charge width &&&
0296 
0297   // apply the lorentz offset correction
0298   yPos = yPos + shiftY;
0299 
0300   // Apply irradiation corrections
0301   if (IrradiationBiasCorrection_) {
0302     if (theClusterParam.theCluster->sizeX() == 1) {  // size=1
0303       // ggiurgiu@jhu.edu, 02/03/09 : for size = 1, the Lorentz shift is already accounted by the irradiation correction
0304       //float tmp1 =  (0.5 * theDetParam.lorentzShiftInCmX);
0305       //cout << "Apply correction correction_dx1 = " << theClusterParam.dx1 << " to xPos = " << xPos;
0306       xPos = xPos - (0.5f * theDetParam.lorentzShiftInCmX);
0307       // Find if pixel is double (big).
0308       bool bigInX = theDetParam.theRecTopol->isItBigPixelInX(theClusterParam.theCluster->maxPixelRow());
0309       if (!bigInX)
0310         xPos -= theClusterParam.dx1;
0311       else
0312         xPos -= theClusterParam.dx2;
0313       //cout<<" to "<<xPos<<" "<<(tmp1+theClusterParam.dx1)<<endl;
0314     } else {  // size>1
0315       //cout << "Apply correction correction_deltax = " << theClusterParam.deltax << " to xPos = " << xPos;
0316       xPos -= theClusterParam.deltax;
0317       //cout<<" to "<<xPos<<endl;
0318     }
0319 
0320     if (theClusterParam.theCluster->sizeY() == 1) {
0321       // ggiurgiu@jhu.edu, 02/03/09 : for size = 1, the Lorentz shift is already accounted by the irradiation correction
0322       yPos = yPos - (0.5f * theDetParam.lorentzShiftInCmY);
0323 
0324       // Find if pixel is double (big).
0325       bool bigInY = theDetParam.theRecTopol->isItBigPixelInY(theClusterParam.theCluster->maxPixelCol());
0326       if (!bigInY)
0327         yPos -= theClusterParam.dy1;
0328       else
0329         yPos -= theClusterParam.dy2;
0330 
0331     } else {
0332       //cout << "Apply correction correction_deltay = " << theClusterParam.deltay << " to yPos = " << yPos << endl;
0333       yPos -= theClusterParam.deltay;
0334     }
0335 
0336   }  // if ( IrradiationBiasCorrection_ )
0337 
0338   //cout<<" in PixelCPEGeneric:localPosition - pos = "<<xPos<<" "<<yPos<<endl; //dk
0339 
0340   //--- Now put the two together
0341   LocalPoint pos_in_local(xPos, yPos);
0342   return pos_in_local;
0343 }
0344 
0345 //==============  INFLATED ERROR AND ERRORS FROM DB BELOW  ================
0346 
0347 //-------------------------------------------------------------------------
0348 //  Hit error in the local frame
0349 //-------------------------------------------------------------------------
0350 LocalError PixelCPEGeneric::localError(DetParam const& theDetParam, ClusterParam& theClusterParamBase) const {
0351   ClusterParamGeneric& theClusterParam = static_cast<ClusterParamGeneric&>(theClusterParamBase);
0352 
0353   // local variables
0354   float xerr, yerr;
0355   bool edgex, edgey, bigInX, bigInY;
0356   int maxPixelCol, maxPixelRow, minPixelCol, minPixelRow;
0357   uint sizex, sizey;
0358 
0359   initializeLocalErrorVariables(xerr,
0360                                 yerr,
0361                                 edgex,
0362                                 edgey,
0363                                 bigInX,
0364                                 bigInY,
0365                                 maxPixelCol,
0366                                 maxPixelRow,
0367                                 minPixelCol,
0368                                 minPixelRow,
0369                                 sizex,
0370                                 sizey,
0371                                 theDetParam,
0372                                 theClusterParam);
0373 
0374   bool useTempErrors =
0375       useErrorsFromTemplates_ && (!NoTemplateErrorsWhenNoTrkAngles_ || theClusterParam.with_track_angle);
0376 
0377   if (int(sizex) != (maxPixelRow - minPixelRow + 1))
0378     LogDebug("PixelCPEGeneric") << " wrong x";
0379   if (int(sizey) != (maxPixelCol - minPixelCol + 1))
0380     LogDebug("PixelCPEGeneric") << " wrong y";
0381 
0382   LogDebug("PixelCPEGeneric") << " edge clus " << xerr << " " << yerr;  //dk
0383   if (bigInX || bigInY)
0384     LogDebug("PixelCPEGeneric") << " big " << bigInX << " " << bigInY;
0385   if (edgex || edgey)
0386     LogDebug("PixelCPEGeneric") << " edge " << edgex << " " << edgey;
0387   LogDebug("PixelCPEGeneric") << " before if " << useErrorsFromTemplates_ << " " << theClusterParam.qBin_;
0388   if (theClusterParam.qBin_ == 0)
0389     LogDebug("PixelCPEGeneric") << " qbin 0! " << edgex << " " << edgey << " " << bigInX << " " << bigInY << " "
0390                                 << sizex << " " << sizey;
0391 
0392   // from PixelCPEGenericBase
0393   setXYErrors(xerr, yerr, edgex, edgey, sizex, sizey, bigInX, bigInY, useTempErrors, theDetParam, theClusterParam);
0394 
0395   if (!useTempErrors) {
0396     LogDebug("PixelCPEGeneric") << "Track angles are not known.\n"
0397                                 << "Default angle estimation which assumes track from PV (0,0,0) does not work.";
0398   }
0399 
0400   if (!useTempErrors && inflate_errors) {
0401     int n_bigx = 0;
0402     int n_bigy = 0;
0403 
0404     for (int irow = 0; irow < 7; ++irow) {
0405       if (theDetParam.theRecTopol->isItBigPixelInX(irow + minPixelRow))
0406         ++n_bigx;
0407     }
0408 
0409     for (int icol = 0; icol < 21; ++icol) {
0410       if (theDetParam.theRecTopol->isItBigPixelInY(icol + minPixelCol))
0411         ++n_bigy;
0412     }
0413 
0414     xerr = (float)(sizex + n_bigx) * theDetParam.thePitchX / std::sqrt(12.0f);
0415     yerr = (float)(sizey + n_bigy) * theDetParam.thePitchY / std::sqrt(12.0f);
0416   }
0417 
0418 #ifdef EDM_ML_DEBUG
0419   if (!(xerr > 0.0))
0420     throw cms::Exception("PixelCPEGeneric::localError") << "\nERROR: Negative pixel error xerr = " << xerr << "\n\n";
0421 
0422   if (!(yerr > 0.0))
0423     throw cms::Exception("PixelCPEGeneric::localError") << "\nERROR: Negative pixel error yerr = " << yerr << "\n\n";
0424 #endif
0425 
0426   LogDebug("PixelCPEGeneric") << " errors  " << xerr << " " << yerr;  //dk
0427   if (theClusterParam.qBin_ == 0)
0428     LogDebug("PixelCPEGeneric") << " qbin 0 " << xerr << " " << yerr;
0429 
0430   auto xerr_sq = xerr * xerr;
0431   auto yerr_sq = yerr * yerr;
0432 
0433   return LocalError(xerr_sq, 0, yerr_sq);
0434 }
0435 
0436 void PixelCPEGeneric::fillPSetDescription(edm::ParameterSetDescription& desc) {
0437   PixelCPEGenericBase::fillPSetDescription(desc);
0438   desc.add<double>("eff_charge_cut_highX", 1.0);
0439   desc.add<double>("eff_charge_cut_highY", 1.0);
0440   desc.add<double>("eff_charge_cut_lowX", 0.0);
0441   desc.add<double>("eff_charge_cut_lowY", 0.0);
0442   desc.add<double>("size_cutX", 3.0);
0443   desc.add<double>("size_cutY", 3.0);
0444   desc.add<double>("EdgeClusterErrorX", 50.0);
0445   desc.add<double>("EdgeClusterErrorY", 85.0);
0446   desc.add<bool>("inflate_errors", false);
0447   desc.add<bool>("inflate_all_errors_no_trk_angle", false);
0448   desc.add<bool>("NoTemplateErrorsWhenNoTrkAngles", false);
0449   desc.add<bool>("UseErrorsFromTemplates", true);
0450   desc.add<bool>("TruncatePixelCharge", true);
0451   desc.add<bool>("IrradiationBiasCorrection", false);
0452   desc.add<bool>("DoCosmics", false);
0453   desc.add<bool>("Upgrade", false);
0454   desc.add<bool>("SmallPitch", false);
0455 }