Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:00:40

0001 // Include our own header first
0002 #include "RecoLocalTracker/SiPixelRecHits/interface/PixelCPETemplateReco.h"
0003 
0004 // Geometry services
0005 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0006 #include "Geometry/TrackerGeometryBuilder/interface/RectangularPixelTopology.h"
0007 
0008 //#define DEBUG
0009 
0010 // MessageLogger
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 
0013 // Magnetic field
0014 #include "MagneticField/Engine/interface/MagneticField.h"
0015 
0016 // The template header files
0017 #include "RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateReco.h"
0018 
0019 // Commented for now (3/10/17) until we figure out how to resuscitate 2D template splitter
0020 /// #include "RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateSplit.h"
0021 
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023 
0024 #include <vector>
0025 #include "boost/multi_array.hpp"
0026 
0027 #include <iostream>
0028 
0029 using namespace SiPixelTemplateReco;
0030 //using namespace SiPixelTemplateSplit;
0031 using namespace std;
0032 
0033 namespace {
0034   constexpr float micronsToCm = 1.0e-4;
0035   constexpr int cluster_matrix_size_x = 13;
0036   constexpr int cluster_matrix_size_y = 21;
0037 }  // namespace
0038 
0039 //-----------------------------------------------------------------------------
0040 //  Constructor.
0041 //
0042 //-----------------------------------------------------------------------------
0043 PixelCPETemplateReco::PixelCPETemplateReco(edm::ParameterSet const& conf,
0044                                            const MagneticField* mag,
0045                                            const TrackerGeometry& geom,
0046                                            const TrackerTopology& ttopo,
0047                                            const SiPixelLorentzAngle* lorentzAngle,
0048                                            const std::vector<SiPixelTemplateStore>* templateStore,
0049                                            const SiPixelTemplateDBObject* templateDBobject)
0050     : PixelCPEBase(conf, mag, geom, ttopo, lorentzAngle, nullptr, templateDBobject, nullptr, 1) {
0051   //cout << endl;
0052   //cout << "Constructing PixelCPETemplateReco::PixelCPETemplateReco(...)................................................." << endl;
0053   //cout << endl;
0054 
0055   // Configurable parameters
0056   //DoCosmics_ = conf.getParameter<bool>("DoCosmics"); // Not used in templates
0057   //LoadTemplatesFromDB_ = conf.getParameter<bool>("LoadTemplatesFromDB"); // Moved to Base
0058 
0059   //cout << " PixelCPETemplateReco : (int)LoadTemplatesFromDB_ = " << (int)LoadTemplatesFromDB_ << endl;
0060   //cout << "field_magnitude = " << field_magnitude << endl;
0061 
0062   // configuration parameter to decide between DB or text file template access
0063 
0064   if (LoadTemplatesFromDB_) {
0065     //cout << "PixelCPETemplateReco: Loading templates from database (DB) --------- " << endl;
0066     thePixelTemp_ = templateStore;
0067   } else {
0068     //cout << "PixelCPETemplateReco : Loading templates for barrel and forward from ASCII files ----------" << endl;
0069     barrelTemplateID_ = conf.getParameter<int>("barrelTemplateID");
0070     forwardTemplateID_ = conf.getParameter<int>("forwardTemplateID");
0071     templateDir_ = conf.getParameter<int>("directoryWithTemplates");
0072 
0073     thePixelTemp_ = &thePixelTempCache_;
0074     if (!SiPixelTemplate::pushfile(barrelTemplateID_, thePixelTempCache_, templateDir_))
0075       throw cms::Exception("PixelCPETemplateReco")
0076           << "\nERROR: Template ID " << barrelTemplateID_
0077           << " not loaded correctly from text file. Reconstruction will fail.\n\n";
0078 
0079     if (!SiPixelTemplate::pushfile(forwardTemplateID_, thePixelTempCache_, templateDir_))
0080       throw cms::Exception("PixelCPETemplateReco")
0081           << "\nERROR: Template ID " << forwardTemplateID_
0082           << " not loaded correctly from text file. Reconstruction will fail.\n\n";
0083   }
0084 
0085   speed_ = conf.getParameter<int>("speed");
0086   LogDebug("PixelCPETemplateReco::PixelCPETemplateReco:") << "Template speed = " << speed_ << "\n";
0087 
0088   UseClusterSplitter_ = conf.getParameter<bool>("UseClusterSplitter");
0089 }
0090 
0091 //-----------------------------------------------------------------------------
0092 //  Clean up.
0093 //-----------------------------------------------------------------------------
0094 PixelCPETemplateReco::~PixelCPETemplateReco() {}
0095 
0096 std::unique_ptr<PixelCPEBase::ClusterParam> PixelCPETemplateReco::createClusterParam(const SiPixelCluster& cl) const {
0097   return std::make_unique<ClusterParamTemplate>(cl);
0098 }
0099 
0100 //------------------------------------------------------------------
0101 //  Public methods mandated by the base class.
0102 //------------------------------------------------------------------
0103 
0104 //------------------------------------------------------------------
0105 //  The main call to the template code.
0106 //------------------------------------------------------------------
0107 LocalPoint PixelCPETemplateReco::localPosition(DetParam const& theDetParam, ClusterParam& theClusterParamBase) const {
0108   ClusterParamTemplate& theClusterParam = static_cast<ClusterParamTemplate&>(theClusterParamBase);
0109 
0110   if (!GeomDetEnumerators::isTrackerPixel(theDetParam.thePart))
0111     throw cms::Exception("PixelCPETemplateReco::localPosition :") << "A non-pixel detector type in here?";
0112   //  barrel(false) or forward(true)
0113   const bool fpix = GeomDetEnumerators::isEndcap(theDetParam.thePart);
0114 
0115   int ID = -9999;
0116   if (LoadTemplatesFromDB_) {
0117     int ID0 = templateDBobject_->getTemplateID(theDetParam.theDet->geographicalId());  // just to comapre
0118     ID = theDetParam.detTemplateId;
0119     if (ID0 != ID)
0120       edm::LogError("PixelCPETemplateReco") << " different id" << ID << " " << ID0 << endl;
0121   } else {  // from asci file
0122     if (!fpix)
0123       ID = barrelTemplateID_;  // barrel
0124     else
0125       ID = forwardTemplateID_;  // forward
0126   }
0127   //cout << "PixelCPETemplateReco : ID = " << ID << endl;
0128 
0129   SiPixelTemplate templ(*thePixelTemp_);
0130 
0131   // Preparing to retrieve ADC counts from the SiPixeltheClusterParam.theCluster->  In the cluster,
0132   // we have the following:
0133   //   int minPixelRow(); // Minimum pixel index in the x direction (low edge).
0134   //   int maxPixelRow(); // Maximum pixel index in the x direction (top edge).
0135   //   int minPixelCol(); // Minimum pixel index in the y direction (left edge).
0136   //   int maxPixelCol(); // Maximum pixel index in the y direction (right edge).
0137   // So the pixels from minPixelRow() will go into clust_array_2d[0][*],
0138   // and the pixels from minPixelCol() will go into clust_array_2d[*][0].
0139 
0140   int row_offset = theClusterParam.theCluster->minPixelRow();
0141   int col_offset = theClusterParam.theCluster->minPixelCol();
0142 
0143   // Store the coordinates of the center of the (0,0) pixel of the array that
0144   // gets passed to PixelTempReco1D
0145   // Will add these values to the output of  PixelTempReco1D
0146   float tmp_x = float(row_offset) + 0.5f;
0147   float tmp_y = float(col_offset) + 0.5f;
0148 
0149   // Store these offsets (to be added later) in a LocalPoint after tranforming
0150   // them from measurement units (pixel units) to local coordinates (cm)
0151   //
0152   //
0153 
0154   // In case of template reco failure, these are the lorentz drift corrections
0155   // to be applied
0156   float lorentzshiftX = 0.5f * theDetParam.lorentzShiftInCmX;
0157   float lorentzshiftY = 0.5f * theDetParam.lorentzShiftInCmY;
0158 
0159   // ggiurgiu@jhu.edu 12/09/2010 : update call with trk angles needed for bow/kink corrections
0160   LocalPoint lp;
0161 
0162   if (theClusterParam.with_track_angle)
0163     lp = theDetParam.theTopol->localPosition(MeasurementPoint(tmp_x, tmp_y), theClusterParam.loc_trk_pred);
0164   else {
0165     edm::LogError("PixelCPETemplateReco") << "@SUB = PixelCPETemplateReco::localPosition"
0166                                           << "Should never be here. PixelCPETemplateReco should always be called with "
0167                                              "track angles. This is a bad error !!! ";
0168 
0169     lp = theDetParam.theTopol->localPosition(MeasurementPoint(tmp_x, tmp_y));
0170   }
0171 
0172   // first compute matrix size
0173   int mrow = 0, mcol = 0;
0174   for (int i = 0; i != theClusterParam.theCluster->size(); ++i) {
0175     auto pix = theClusterParam.theCluster->pixel(i);
0176     int irow = int(pix.x);
0177     int icol = int(pix.y);
0178     mrow = std::max(mrow, irow);
0179     mcol = std::max(mcol, icol);
0180   }
0181   mrow -= row_offset;
0182   mrow += 1;
0183   mrow = std::min(mrow, cluster_matrix_size_x);
0184   mcol -= col_offset;
0185   mcol += 1;
0186   mcol = std::min(mcol, cluster_matrix_size_y);
0187   assert(mrow > 0);
0188   assert(mcol > 0);
0189 
0190   float clustMatrix[mrow][mcol];
0191   memset(clustMatrix, 0, sizeof(float) * mrow * mcol);
0192 
0193   // Copy clust's pixels (calibrated in electrons) into clusMatrix;
0194   for (int i = 0; i != theClusterParam.theCluster->size(); ++i) {
0195     auto pix = theClusterParam.theCluster->pixel(i);
0196     int irow = int(pix.x) - row_offset;
0197     int icol = int(pix.y) - col_offset;
0198 
0199     // Gavril : what do we do here if the row/column is larger than cluster_matrix_size_x/cluster_matrix_size_y  ?
0200     // Ignore them for the moment...
0201     if ((irow < mrow) & (icol < mcol))
0202       clustMatrix[irow][icol] = float(pix.adc);
0203   }
0204 
0205   // Make and fill the bool arrays flagging double pixels
0206   bool xdouble[mrow], ydouble[mcol];
0207   // x directions (shorter), rows
0208   for (int irow = 0; irow < mrow; ++irow)
0209     xdouble[irow] = theDetParam.theRecTopol->isItBigPixelInX(irow + row_offset);
0210 
0211   // y directions (longer), columns
0212   for (int icol = 0; icol < mcol; ++icol)
0213     ydouble[icol] = theDetParam.theRecTopol->isItBigPixelInY(icol + col_offset);
0214 
0215   SiPixelTemplateReco::ClusMatrix clusterPayload{&clustMatrix[0][0], xdouble, ydouble, mrow, mcol};
0216 
0217   // Output:
0218   float nonsense = -99999.9f;  // nonsense init value
0219   theClusterParam.templXrec_ = theClusterParam.templYrec_ = theClusterParam.templSigmaX_ =
0220       theClusterParam.templSigmaY_ = nonsense;
0221   // If the template recontruction fails, we want to return 1.0 for now
0222   theClusterParam.templProbY_ = theClusterParam.templProbX_ = theClusterParam.templProbQ_ = 1.0f;
0223   theClusterParam.templQbin_ = 0;
0224   // We have a boolean denoting whether the reco failed or not
0225   theClusterParam.hasFilledProb_ = false;
0226 
0227   float templYrec1_ = nonsense;
0228   float templXrec1_ = nonsense;
0229   float templYrec2_ = nonsense;
0230   float templXrec2_ = nonsense;
0231 
0232   // ******************************************************************
0233   // Do it! Use cotalpha_ and cotbeta_ calculated in PixelCPEBase
0234 
0235   float locBz = theDetParam.bz;
0236   float locBx = theDetParam.bx;
0237 
0238   theClusterParam.ierr = PixelTempReco1D(ID,
0239                                          theClusterParam.cotalpha,
0240                                          theClusterParam.cotbeta,
0241                                          locBz,
0242                                          locBx,
0243                                          clusterPayload,
0244                                          templ,
0245                                          theClusterParam.templYrec_,
0246                                          theClusterParam.templSigmaY_,
0247                                          theClusterParam.templProbY_,
0248                                          theClusterParam.templXrec_,
0249                                          theClusterParam.templSigmaX_,
0250                                          theClusterParam.templProbX_,
0251                                          theClusterParam.templQbin_,
0252                                          speed_,
0253                                          theClusterParam.templProbQ_);
0254 
0255   // ******************************************************************
0256 
0257   // Check exit status
0258   if UNLIKELY (theClusterParam.ierr != 0) {
0259     LogDebug("PixelCPETemplateReco::localPosition")
0260         << "reconstruction failed with error " << theClusterParam.ierr << "\n";
0261 
0262     // Template reco has failed, compute position estimates based on cluster center of gravity + Lorentz drift
0263     // Future improvement would be to call generic reco instead
0264 
0265     // ggiurgiu@jhu.edu, 21/09/2010 : trk angles needed to correct for bows/kinks
0266     if (theClusterParam.with_track_angle) {
0267       theClusterParam.templXrec_ =
0268           theDetParam.theTopol->localX(theClusterParam.theCluster->x(), theClusterParam.loc_trk_pred) + lorentzshiftX;
0269       theClusterParam.templYrec_ =
0270           theDetParam.theTopol->localY(theClusterParam.theCluster->y(), theClusterParam.loc_trk_pred) + lorentzshiftY;
0271     } else {
0272       edm::LogError("PixelCPETemplateReco") << "@SUB = PixelCPETemplateReco::localPosition"
0273                                             << "Should never be here. PixelCPETemplateReco should always be called "
0274                                                "with track angles. This is a bad error !!! ";
0275 
0276       theClusterParam.templXrec_ = theDetParam.theTopol->localX(theClusterParam.theCluster->x()) + lorentzshiftX;
0277       theClusterParam.templYrec_ = theDetParam.theTopol->localY(theClusterParam.theCluster->y()) + lorentzshiftY;
0278     }
0279   } else if UNLIKELY (UseClusterSplitter_ && theClusterParam.templQbin_ == 0) {
0280     edm::LogError("PixelCPETemplateReco") << " PixelCPETemplateReco: Qbin = 0 but using cluster splitter, we should "
0281                                              "never be here !!!!!!!!!!!!!!!!!!!!!! \n"
0282                                           << "(int)UseClusterSplitter_ = " << (int)UseClusterSplitter_ << endl;
0283 
0284     //ierr =
0285     //PixelTempSplit( ID, fpix, cotalpha_, cotbeta_,
0286     //      clust_array_2d, ydouble, xdouble,
0287     //      templ,
0288     //      templYrec1_, templYrec2_, templSigmaY_, templProbY_,
0289     //      templXrec1_, templXrec2_, templSigmaX_, templProbX_,
0290     //      templQbin_ );
0291 
0292     // Commented for now (3/10/17) until we figure out how to resuscitate 2D template splitter
0293     ///      std::vector< SiPixelTemplateStore2D > thePixelTemp2D_;
0294     ///SiPixelTemplate2D::pushfile(ID, thePixelTemp2D_);
0295     ///SiPixelTemplate2D templ2D_(thePixelTemp2D_);
0296 
0297     theClusterParam.ierr = -123;
0298     /*
0299        float dchisq;
0300        float templProbQ_;
0301        SiPixelTemplateSplit::PixelTempSplit( ID, theClusterParam.cotalpha, theClusterParam.cotbeta,
0302        clust_array_2d,
0303        ydouble, xdouble,
0304        templ,
0305        templYrec1_, templYrec2_, theClusterParam.templSigmaY_, theClusterParam.templProbY_,
0306        templXrec1_, templXrec2_, theClusterParam.templSigmaX_, theClusterParam.templProbX_,
0307        theClusterParam.templQbin_,
0308        templProbQ_,
0309        true,
0310        dchisq,
0311        templ2D_ );
0312        
0313        */
0314     if (theClusterParam.ierr != 0) {
0315       // Template reco has failed, compute position estimates based on cluster center of gravity + Lorentz drift
0316       // Future improvement would be to call generic reco instead
0317 
0318       // ggiurgiu@jhu.edu, 12/09/2010 : trk angles needed to correct for bows/kinks
0319       if (theClusterParam.with_track_angle) {
0320         theClusterParam.templXrec_ =
0321             theDetParam.theTopol->localX(theClusterParam.theCluster->x(), theClusterParam.loc_trk_pred) + lorentzshiftX;
0322         theClusterParam.templYrec_ =
0323             theDetParam.theTopol->localY(theClusterParam.theCluster->y(), theClusterParam.loc_trk_pred) + lorentzshiftY;
0324       } else {
0325         edm::LogError("PixelCPETemplateReco") << "@SUB = PixelCPETemplateReco::localPosition"
0326                                               << "Should never be here. PixelCPETemplateReco should always be called "
0327                                                  "with track angles. This is a bad error !!! ";
0328         theClusterParam.templXrec_ = theDetParam.theTopol->localX(theClusterParam.theCluster->x()) + lorentzshiftX;
0329         theClusterParam.templYrec_ = theDetParam.theTopol->localY(theClusterParam.theCluster->y()) + lorentzshiftY;
0330       }
0331     } else {
0332       // go from micrometer to centimeter
0333       templXrec1_ *= micronsToCm;
0334       templYrec1_ *= micronsToCm;
0335       templXrec2_ *= micronsToCm;
0336       templYrec2_ *= micronsToCm;
0337 
0338       // go back to the module coordinate system
0339       templXrec1_ += lp.x();
0340       templYrec1_ += lp.y();
0341       templXrec2_ += lp.x();
0342       templYrec2_ += lp.y();
0343 
0344       // calculate distance from each hit to the track and choose the hit closest to the track
0345       float distX1 = std::abs(templXrec1_ - theClusterParam.trk_lp_x);
0346       float distX2 = std::abs(templXrec2_ - theClusterParam.trk_lp_x);
0347       float distY1 = std::abs(templYrec1_ - theClusterParam.trk_lp_y);
0348       float distY2 = std::abs(templYrec2_ - theClusterParam.trk_lp_y);
0349       theClusterParam.templXrec_ = (distX1 < distX2 ? templXrec1_ : templXrec2_);
0350       theClusterParam.templYrec_ = (distY1 < distY2 ? templYrec1_ : templYrec2_);
0351     }
0352   }  // else if ( UseClusterSplitter_ && templQbin_ == 0 )
0353 
0354   else  // apparenly this is the good one!
0355   {
0356     // go from micrometer to centimeter
0357     theClusterParam.templXrec_ *= micronsToCm;
0358     theClusterParam.templYrec_ *= micronsToCm;
0359 
0360     // go back to the module coordinate system
0361     theClusterParam.templXrec_ += lp.x();
0362     theClusterParam.templYrec_ += lp.y();
0363 
0364     // Compute the Alignment Group Corrections [template ID should already be selected from call to reco procedure]
0365     if (doLorentzFromAlignment_) {
0366       // Do only if the lotentzshift has meaningfull numbers
0367       if (theDetParam.lorentzShiftInCmX != 0.0 || theDetParam.lorentzShiftInCmY != 0.0) {
0368         // the LA width/shift returned by templates use (+)
0369         // the LA width/shift produced by PixelCPEBase for positive LA is (-)
0370         // correct this by iserting (-)
0371         //float temp1 = -micronsToCm*templ.lorxwidth();  // old
0372         //float temp2 = -micronsToCm*templ.lorywidth();  // does not incl 1/2
0373         float templateLorbiasCmX = -micronsToCm * templ.lorxbias();  // new
0374         float templateLorbiasCmY = -micronsToCm * templ.lorybias();  //incl. 1/2
0375         // now, correctly, we can use the difference of shifts
0376         //theClusterParam.templXrec_ += 0.5*(theDetParam.lorentzShiftInCmX - templateLorbiasCmX);
0377         //theClusterParam.templYrec_ += 0.5*(theDetParam.lorentzShiftInCmY - templateLorbiasCmY);
0378         theClusterParam.templXrec_ += (0.5 * (theDetParam.lorentzShiftInCmX) - templateLorbiasCmX);
0379         theClusterParam.templYrec_ += (0.5 * (theDetParam.lorentzShiftInCmY) - templateLorbiasCmY);
0380         //cout << "Templates: la lorentz offset = "
0381         //   <<(0.5*(theDetParam.lorentzShiftInCmX)-templateLorbiasCmX)
0382         //   <<" "<<templateLorbiasCmX<<" "<<templateLorbiasCmY
0383         //   <<" "<<temp1<<" "<<temp2
0384         //   <<" "<<theDetParam.lorentzShiftInCmX
0385         //   <<" "<<theDetParam.lorentzShiftInCmY
0386         //   << endl; //dk
0387       }  //else {cout<<" LA is 0, disable offset corrections "<<endl;} //dk
0388     }    //else {cout<<" Do not do LA offset correction "<<endl;} //dk
0389   }
0390 
0391   // Save probabilities and qBin in the quantities given to us by the base class
0392   // (for which there are also inline getters).  &&& templProbX_ etc. should be retired...
0393   theClusterParam.probabilityX_ = theClusterParam.templProbX_;
0394   theClusterParam.probabilityY_ = theClusterParam.templProbY_;
0395   theClusterParam.probabilityQ_ = theClusterParam.templProbQ_;
0396   theClusterParam.qBin_ = theClusterParam.templQbin_;
0397 
0398   if (theClusterParam.ierr == 0)  // always true here
0399     theClusterParam.hasFilledProb_ = true;
0400 
0401   return LocalPoint(theClusterParam.templXrec_, theClusterParam.templYrec_);
0402 }
0403 
0404 //------------------------------------------------------------------
0405 //  localError() relies on localPosition() being called FIRST!!!
0406 //------------------------------------------------------------------
0407 LocalError PixelCPETemplateReco::localError(DetParam const& theDetParam, ClusterParam& theClusterParamBase) const {
0408   ClusterParamTemplate& theClusterParam = static_cast<ClusterParamTemplate&>(theClusterParamBase);
0409 
0410   //cout << endl;
0411   //cout << "Set PixelCPETemplate errors .............................................." << endl;
0412 
0413   //cout << "CPETemplate : " << endl;
0414 
0415   //--- Default is the maximum error used for edge clusters.
0416   //--- (never used, in fact: let comment it out, shut up the complains of the static analyzer, and save a few CPU cycles)
0417   //   const float sig12 = 1./sqrt(12.0);
0418   //   float xerr = theDetParam.thePitchX *sig12;
0419   //   float yerr = theDetParam.thePitchY *sig12;
0420   float xerr, yerr;
0421 
0422   // Check if the errors were already set at the clusters splitting level
0423   if (theClusterParam.theCluster->getSplitClusterErrorX() > 0.0f &&
0424       theClusterParam.theCluster->getSplitClusterErrorX() < clusterSplitMaxError_ &&
0425       theClusterParam.theCluster->getSplitClusterErrorY() > 0.0f &&
0426       theClusterParam.theCluster->getSplitClusterErrorY() < clusterSplitMaxError_) {
0427     xerr = theClusterParam.theCluster->getSplitClusterErrorX() * micronsToCm;
0428     yerr = theClusterParam.theCluster->getSplitClusterErrorY() * micronsToCm;
0429 
0430     //cout << "Errors set at cluster splitting level : " << endl;
0431     //cout << "xerr = " << xerr << endl;
0432     //cout << "yerr = " << yerr << endl;
0433   } else {
0434     // If errors are not split at the cluster splitting level, set the errors here
0435 
0436     //cout  << "Errors are not split at the cluster splitting level, set the errors here : " << endl;
0437 
0438     int maxPixelCol = theClusterParam.theCluster->maxPixelCol();
0439     int maxPixelRow = theClusterParam.theCluster->maxPixelRow();
0440     int minPixelCol = theClusterParam.theCluster->minPixelCol();
0441     int minPixelRow = theClusterParam.theCluster->minPixelRow();
0442 
0443     //--- Are we near either of the edges?
0444     bool edgex = (theDetParam.theRecTopol->isItEdgePixelInX(minPixelRow) ||
0445                   theDetParam.theRecTopol->isItEdgePixelInX(maxPixelRow));
0446     bool edgey = (theDetParam.theRecTopol->isItEdgePixelInY(minPixelCol) ||
0447                   theDetParam.theRecTopol->isItEdgePixelInY(maxPixelCol));
0448 
0449     if (theClusterParam.ierr != 0) {
0450       // If reconstruction fails the hit position is calculated from cluster center of gravity
0451       // corrected in x by average Lorentz drift. Assign huge errors.
0452       //xerr = 10.0 * (float)theClusterParam.theCluster->sizeX() * xerr;
0453       //yerr = 10.0 * (float)theClusterParam.theCluster->sizeX() * yerr;
0454 
0455       if (!GeomDetEnumerators::isTrackerPixel(theDetParam.thePart))
0456         throw cms::Exception("PixelCPETemplateReco::localPosition :") << "A non-pixel detector type in here?";
0457 
0458       // Assign better errors based on the residuals for failed template cases
0459       if (GeomDetEnumerators::isBarrel(theDetParam.thePart)) {
0460         xerr = 55.0f * micronsToCm;
0461         yerr = 36.0f * micronsToCm;
0462       } else {
0463         xerr = 42.0f * micronsToCm;
0464         yerr = 39.0f * micronsToCm;
0465       }
0466 
0467       //cout << "xerr = " << xerr << endl;
0468       //cout << "yerr = " << yerr << endl;
0469 
0470       //return LocalError(xerr*xerr, 0, yerr*yerr);
0471     } else if (edgex || edgey) {
0472       // for edge pixels assign errors according to observed residual RMS
0473       if (edgex && !edgey) {
0474         xerr = xEdgeXError_ * micronsToCm;
0475         yerr = xEdgeYError_ * micronsToCm;
0476       } else if (!edgex && edgey) {
0477         xerr = yEdgeXError_ * micronsToCm;
0478         yerr = yEdgeYError_ * micronsToCm;
0479       } else if (edgex && edgey) {
0480         xerr = bothEdgeXError_ * micronsToCm;
0481         yerr = bothEdgeYError_ * micronsToCm;
0482       } else {
0483         throw cms::Exception(" PixelCPETemplateReco::localError: Something wrong with pixel edge flag !!!");
0484       }
0485 
0486       //cout << "xerr = " << xerr << endl;
0487       //cout << "yerr = " << yerr << endl;
0488     } else {
0489       // &&& need a class const
0490       //const float micronsToCm = 1.0e-4;
0491 
0492       xerr = theClusterParam.templSigmaX_ * micronsToCm;
0493       yerr = theClusterParam.templSigmaY_ * micronsToCm;
0494 
0495       //cout << "xerr = " << xerr << endl;
0496       //cout << "yerr = " << yerr << endl;
0497 
0498       // &&& should also check ierr (saved as class variable) and return
0499       // &&& nonsense (another class static) if the template fit failed.
0500     }
0501 
0502     if (theVerboseLevel > 9) {
0503       LogDebug("PixelCPETemplateReco") << " Sizex = " << theClusterParam.theCluster->sizeX()
0504                                        << " Sizey = " << theClusterParam.theCluster->sizeY() << " Edgex = " << edgex
0505                                        << " Edgey = " << edgey << " ErrX  = " << xerr << " ErrY  = " << yerr;
0506     }
0507 
0508   }  // else
0509 
0510   if (!(xerr > 0.0f))
0511     throw cms::Exception("PixelCPETemplateReco::localError")
0512         << "\nERROR: Negative pixel error xerr = " << xerr << "\n\n";
0513 
0514   if (!(yerr > 0.0f))
0515     throw cms::Exception("PixelCPETemplateReco::localError")
0516         << "\nERROR: Negative pixel error yerr = " << yerr << "\n\n";
0517 
0518   //cout << "Final errors set to: " << endl;
0519   //cout << "xerr = " << xerr << endl;
0520   //cout << "yerr = " << yerr << endl;
0521   //cout << "Out of PixelCPETemplateREco..........................................................................." << endl;
0522   //cout << endl;
0523 
0524   return LocalError(xerr * xerr, 0, yerr * yerr);
0525 }
0526 
0527 void PixelCPETemplateReco::fillPSetDescription(edm::ParameterSetDescription& desc) {
0528   desc.add<int>("barrelTemplateID", 0);
0529   desc.add<int>("forwardTemplateID", 0);
0530   desc.add<int>("directoryWithTemplates", 0);
0531   desc.add<int>("speed", -2);
0532   desc.add<bool>("UseClusterSplitter", false);
0533 }