Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-06-25 05:16:13 

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 
0007 //#define DEBUG
0008 
0009 // MessageLogger
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 
0012 // Magnetic field
0013 #include "MagneticField/Engine/interface/MagneticField.h"
0014 
0015 // The template header files
0016 #include "RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateReco.h"
0017 
0018 // Commented for now (3/10/17) until we figure out how to resuscitate 2D template splitter
0019 /// #include "RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateSplit.h"
0020 
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022 
0023 #include <vector>
0024 #include "boost/multi_array.hpp"
0025 
0026 #include <iostream>
0027 
0028 using namespace SiPixelTemplateReco;
0029 //using namespace SiPixelTemplateSplit;
0030 using namespace std;
0031 
0032 namespace {
0033   constexpr float micronsToCm = 1.0e-4;
0034   constexpr int cluster_matrix_size_x = 13;
0035   constexpr int cluster_matrix_size_y = 21;
0036 }  // namespace
0037 
0038 //-----------------------------------------------------------------------------
0039 //  Constructor.
0040 //
0041 //-----------------------------------------------------------------------------
0042 PixelCPETemplateReco::PixelCPETemplateReco(edm::ParameterSet const& conf,
0043                                            const MagneticField* mag,
0044                                            const TrackerGeometry& geom,
0045                                            const TrackerTopology& ttopo,
0046                                            const SiPixelLorentzAngle* lorentzAngle,
0047                                            const std::vector<SiPixelTemplateStore>* templateStore,
0048                                            const SiPixelTemplateDBObject* templateDBobject)
0049     : PixelCPEBase(conf, mag, geom, ttopo, lorentzAngle, nullptr, templateDBobject, nullptr, 1) {
0050   //cout << endl;
0051   //cout << "Constructing PixelCPETemplateReco::PixelCPETemplateReco(...)................................................." << endl;
0052   //cout << endl;
0053 
0054   // Configurable parameters
0055   //DoCosmics_ = conf.getParameter<bool>("DoCosmics"); // Not used in templates
0056   //LoadTemplatesFromDB_ = conf.getParameter<bool>("LoadTemplatesFromDB"); // Moved to Base
0057 
0058   //cout << " PixelCPETemplateReco : (int)LoadTemplatesFromDB_ = " << (int)LoadTemplatesFromDB_ << endl;
0059   //cout << "field_magnitude = " << field_magnitude << endl;
0060 
0061   // configuration parameter to decide between DB or text file template access
0062 
0063   if (LoadTemplatesFromDB_) {
0064     //cout << "PixelCPETemplateReco: Loading templates from database (DB) --------- " << endl;
0065     thePixelTemp_ = templateStore;
0066   } else {
0067     //cout << "PixelCPETemplateReco : Loading templates for barrel and forward from ASCII files ----------" << endl;
0068     barrelTemplateID_ = conf.getParameter<int>("barrelTemplateID");
0069     forwardTemplateID_ = conf.getParameter<int>("forwardTemplateID");
0070     templateDir_ = conf.getParameter<int>("directoryWithTemplates");
0071 
0072     thePixelTemp_ = &thePixelTempCache_;
0073     if (!SiPixelTemplate::pushfile(barrelTemplateID_, thePixelTempCache_, templateDir_))
0074       throw cms::Exception("PixelCPETemplateReco")
0075           << "\nERROR: Template ID " << barrelTemplateID_
0076           << " not loaded correctly from text file. Reconstruction will fail.\n\n";
0077 
0078     if (!SiPixelTemplate::pushfile(forwardTemplateID_, thePixelTempCache_, templateDir_))
0079       throw cms::Exception("PixelCPETemplateReco")
0080           << "\nERROR: Template ID " << forwardTemplateID_
0081           << " not loaded correctly from text file. Reconstruction will fail.\n\n";
0082   }
0083 
0084   goodEdgeAlgo_ = conf.getParameter<bool>("GoodEdgeAlgo");
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.theTopol->isItBigPixelInX(irow + row_offset);
0210 
0211   // y directions (longer), columns
0212   for (int icol = 0; icol < mcol; ++icol)
0213     ydouble[icol] = theDetParam.theTopol->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                                          goodEdgeAlgo_);
0255 
0256   // ******************************************************************
0257 
0258   // Check exit status
0259   if UNLIKELY (theClusterParam.ierr != 0) {
0260     LogDebug("PixelCPETemplateReco::localPosition")
0261         << "reconstruction failed with error " << theClusterParam.ierr << "\n";
0262 
0263     // Template reco has failed, compute position estimates based on cluster center of gravity + Lorentz drift
0264     // Future improvement would be to call generic reco instead
0265 
0266     // ggiurgiu@jhu.edu, 21/09/2010 : trk angles needed to correct for bows/kinks
0267     if (theClusterParam.with_track_angle) {
0268       theClusterParam.templXrec_ =
0269           theDetParam.theTopol->localX(theClusterParam.theCluster->x(), theClusterParam.loc_trk_pred) + lorentzshiftX;
0270       theClusterParam.templYrec_ =
0271           theDetParam.theTopol->localY(theClusterParam.theCluster->y(), theClusterParam.loc_trk_pred) + lorentzshiftY;
0272     } else {
0273       edm::LogError("PixelCPETemplateReco") << "@SUB = PixelCPETemplateReco::localPosition"
0274                                             << "Should never be here. PixelCPETemplateReco should always be called "
0275                                                "with track angles. This is a bad error !!! ";
0276 
0277       theClusterParam.templXrec_ = theDetParam.theTopol->localX(theClusterParam.theCluster->x()) + lorentzshiftX;
0278       theClusterParam.templYrec_ = theDetParam.theTopol->localY(theClusterParam.theCluster->y()) + lorentzshiftY;
0279     }
0280   } else if UNLIKELY (UseClusterSplitter_ && theClusterParam.templQbin_ == 0) {
0281     edm::LogError("PixelCPETemplateReco") << " PixelCPETemplateReco: Qbin = 0 but using cluster splitter, we should "
0282                                              "never be here !!!!!!!!!!!!!!!!!!!!!! \n"
0283                                           << "(int)UseClusterSplitter_ = " << (int)UseClusterSplitter_ << endl;
0284 
0285     //ierr =
0286     //PixelTempSplit( ID, fpix, cotalpha_, cotbeta_,
0287     //      clust_array_2d, ydouble, xdouble,
0288     //      templ,
0289     //      templYrec1_, templYrec2_, templSigmaY_, templProbY_,
0290     //      templXrec1_, templXrec2_, templSigmaX_, templProbX_,
0291     //      templQbin_ );
0292 
0293     // Commented for now (3/10/17) until we figure out how to resuscitate 2D template splitter
0294     ///      std::vector< SiPixelTemplateStore2D > thePixelTemp2D_;
0295     ///SiPixelTemplate2D::pushfile(ID, thePixelTemp2D_);
0296     ///SiPixelTemplate2D templ2D_(thePixelTemp2D_);
0297 
0298     theClusterParam.ierr = -123;
0299     /*
0300        float dchisq;
0301        float templProbQ_;
0302        SiPixelTemplateSplit::PixelTempSplit( ID, theClusterParam.cotalpha, theClusterParam.cotbeta,
0303        clust_array_2d,
0304        ydouble, xdouble,
0305        templ,
0306        templYrec1_, templYrec2_, theClusterParam.templSigmaY_, theClusterParam.templProbY_,
0307        templXrec1_, templXrec2_, theClusterParam.templSigmaX_, theClusterParam.templProbX_,
0308        theClusterParam.templQbin_,
0309        templProbQ_,
0310        true,
0311        dchisq,
0312        templ2D_ );
0313        
0314        */
0315     if (theClusterParam.ierr != 0) {
0316       // Template reco has failed, compute position estimates based on cluster center of gravity + Lorentz drift
0317       // Future improvement would be to call generic reco instead
0318 
0319       // ggiurgiu@jhu.edu, 12/09/2010 : trk angles needed to correct for bows/kinks
0320       if (theClusterParam.with_track_angle) {
0321         theClusterParam.templXrec_ =
0322             theDetParam.theTopol->localX(theClusterParam.theCluster->x(), theClusterParam.loc_trk_pred) + lorentzshiftX;
0323         theClusterParam.templYrec_ =
0324             theDetParam.theTopol->localY(theClusterParam.theCluster->y(), theClusterParam.loc_trk_pred) + lorentzshiftY;
0325       } else {
0326         edm::LogError("PixelCPETemplateReco") << "@SUB = PixelCPETemplateReco::localPosition"
0327                                               << "Should never be here. PixelCPETemplateReco should always be called "
0328                                                  "with track angles. This is a bad error !!! ";
0329         theClusterParam.templXrec_ = theDetParam.theTopol->localX(theClusterParam.theCluster->x()) + lorentzshiftX;
0330         theClusterParam.templYrec_ = theDetParam.theTopol->localY(theClusterParam.theCluster->y()) + lorentzshiftY;
0331       }
0332     } else {
0333       // go from micrometer to centimeter
0334       templXrec1_ *= micronsToCm;
0335       templYrec1_ *= micronsToCm;
0336       templXrec2_ *= micronsToCm;
0337       templYrec2_ *= micronsToCm;
0338 
0339       // go back to the module coordinate system
0340       templXrec1_ += lp.x();
0341       templYrec1_ += lp.y();
0342       templXrec2_ += lp.x();
0343       templYrec2_ += lp.y();
0344 
0345       // calculate distance from each hit to the track and choose the hit closest to the track
0346       float distX1 = std::abs(templXrec1_ - theClusterParam.trk_lp_x);
0347       float distX2 = std::abs(templXrec2_ - theClusterParam.trk_lp_x);
0348       float distY1 = std::abs(templYrec1_ - theClusterParam.trk_lp_y);
0349       float distY2 = std::abs(templYrec2_ - theClusterParam.trk_lp_y);
0350       theClusterParam.templXrec_ = (distX1 < distX2 ? templXrec1_ : templXrec2_);
0351       theClusterParam.templYrec_ = (distY1 < distY2 ? templYrec1_ : templYrec2_);
0352     }
0353   }  // else if ( UseClusterSplitter_ && templQbin_ == 0 )
0354 
0355   else  // apparenly this is the good one!
0356   {
0357     // go from micrometer to centimeter
0358     theClusterParam.templXrec_ *= micronsToCm;
0359     theClusterParam.templYrec_ *= micronsToCm;
0360 
0361     // go back to the module coordinate system
0362     theClusterParam.templXrec_ += lp.x();
0363     theClusterParam.templYrec_ += lp.y();
0364 
0365     // Compute the Alignment Group Corrections [template ID should already be selected from call to reco procedure]
0366     if (doLorentzFromAlignment_) {
0367       // Do only if the lotentzshift has meaningfull numbers
0368       if (theDetParam.lorentzShiftInCmX != 0.0 || theDetParam.lorentzShiftInCmY != 0.0) {
0369         // the LA width/shift returned by templates use (+)
0370         // the LA width/shift produced by PixelCPEBase for positive LA is (-)
0371         // correct this by iserting (-)
0372         //float temp1 = -micronsToCm*templ.lorxwidth();  // old
0373         //float temp2 = -micronsToCm*templ.lorywidth();  // does not incl 1/2
0374         float templateLorbiasCmX = -micronsToCm * templ.lorxbias();  // new
0375         float templateLorbiasCmY = -micronsToCm * templ.lorybias();  //incl. 1/2
0376         // now, correctly, we can use the difference of shifts
0377         //theClusterParam.templXrec_ += 0.5*(theDetParam.lorentzShiftInCmX - templateLorbiasCmX);
0378         //theClusterParam.templYrec_ += 0.5*(theDetParam.lorentzShiftInCmY - templateLorbiasCmY);
0379         theClusterParam.templXrec_ += (0.5 * (theDetParam.lorentzShiftInCmX) - templateLorbiasCmX);
0380         theClusterParam.templYrec_ += (0.5 * (theDetParam.lorentzShiftInCmY) - templateLorbiasCmY);
0381         //cout << "Templates: la lorentz offset = "
0382         //   <<(0.5*(theDetParam.lorentzShiftInCmX)-templateLorbiasCmX)
0383         //   <<" "<<templateLorbiasCmX<<" "<<templateLorbiasCmY
0384         //   <<" "<<temp1<<" "<<temp2
0385         //   <<" "<<theDetParam.lorentzShiftInCmX
0386         //   <<" "<<theDetParam.lorentzShiftInCmY
0387         //   << endl; //dk
0388       }  //else {cout<<" LA is 0, disable offset corrections "<<endl;} //dk
0389     }  //else {cout<<" Do not do LA offset correction "<<endl;} //dk
0390   }
0391 
0392   // Save probabilities and qBin in the quantities given to us by the base class
0393   // (for which there are also inline getters).  &&& templProbX_ etc. should be retired...
0394   theClusterParam.probabilityX_ = theClusterParam.templProbX_;
0395   theClusterParam.probabilityY_ = theClusterParam.templProbY_;
0396   theClusterParam.probabilityQ_ = theClusterParam.templProbQ_;
0397   theClusterParam.qBin_ = theClusterParam.templQbin_;
0398 
0399   if (theClusterParam.ierr == 0)  // always true here
0400     theClusterParam.hasFilledProb_ = true;
0401 
0402   return LocalPoint(theClusterParam.templXrec_, theClusterParam.templYrec_);
0403 }
0404 
0405 //------------------------------------------------------------------
0406 //  localError() relies on localPosition() being called FIRST!!!
0407 //------------------------------------------------------------------
0408 LocalError PixelCPETemplateReco::localError(DetParam const& theDetParam, ClusterParam& theClusterParamBase) const {
0409   ClusterParamTemplate& theClusterParam = static_cast<ClusterParamTemplate&>(theClusterParamBase);
0410 
0411   //cout << endl;
0412   //cout << "Set PixelCPETemplate errors .............................................." << endl;
0413 
0414   //cout << "CPETemplate : " << endl;
0415 
0416   //--- Default is the maximum error used for edge clusters.
0417   //--- (never used, in fact: let comment it out, shut up the complains of the static analyzer, and save a few CPU cycles)
0418   //   const float sig12 = 1./sqrt(12.0);
0419   //   float xerr = theDetParam.thePitchX *sig12;
0420   //   float yerr = theDetParam.thePitchY *sig12;
0421   float xerr, yerr;
0422 
0423   // Check if the errors were already set at the clusters splitting level
0424   if (theClusterParam.theCluster->getSplitClusterErrorX() > 0.0f &&
0425       theClusterParam.theCluster->getSplitClusterErrorX() < clusterSplitMaxError_ &&
0426       theClusterParam.theCluster->getSplitClusterErrorY() > 0.0f &&
0427       theClusterParam.theCluster->getSplitClusterErrorY() < clusterSplitMaxError_) {
0428     xerr = theClusterParam.theCluster->getSplitClusterErrorX() * micronsToCm;
0429     yerr = theClusterParam.theCluster->getSplitClusterErrorY() * micronsToCm;
0430 
0431     //cout << "Errors set at cluster splitting level : " << endl;
0432     //cout << "xerr = " << xerr << endl;
0433     //cout << "yerr = " << yerr << endl;
0434   } else {
0435     // If errors are not split at the cluster splitting level, set the errors here
0436 
0437     //cout  << "Errors are not split at the cluster splitting level, set the errors here : " << endl;
0438 
0439     int maxPixelCol = theClusterParam.theCluster->maxPixelCol();
0440     int maxPixelRow = theClusterParam.theCluster->maxPixelRow();
0441     int minPixelCol = theClusterParam.theCluster->minPixelCol();
0442     int minPixelRow = theClusterParam.theCluster->minPixelRow();
0443 
0444     //--- Are we near either of the edges?
0445     bool edgex =
0446         (theDetParam.theTopol->isItEdgePixelInX(minPixelRow) || theDetParam.theTopol->isItEdgePixelInX(maxPixelRow));
0447     bool edgey =
0448         (theDetParam.theTopol->isItEdgePixelInY(minPixelCol) || theDetParam.theTopol->isItEdgePixelInY(maxPixelCol));
0449 
0450     if (theClusterParam.ierr != 0) {
0451       // If reconstruction fails the hit position is calculated from cluster center of gravity
0452       // corrected in x by average Lorentz drift. Assign huge errors.
0453       //xerr = 10.0 * (float)theClusterParam.theCluster->sizeX() * xerr;
0454       //yerr = 10.0 * (float)theClusterParam.theCluster->sizeX() * yerr;
0455 
0456       if (!GeomDetEnumerators::isTrackerPixel(theDetParam.thePart))
0457         throw cms::Exception("PixelCPETemplateReco::localPosition :") << "A non-pixel detector type in here?";
0458 
0459       // Assign better errors based on the residuals for failed template cases
0460       if (GeomDetEnumerators::isBarrel(theDetParam.thePart)) {
0461         xerr = 55.0f * micronsToCm;
0462         yerr = 36.0f * micronsToCm;
0463       } else {
0464         xerr = 42.0f * micronsToCm;
0465         yerr = 39.0f * micronsToCm;
0466       }
0467 
0468       //cout << "xerr = " << xerr << endl;
0469       //cout << "yerr = " << yerr << endl;
0470 
0471       //return LocalError(xerr*xerr, 0, yerr*yerr);
0472     } else if (edgex || edgey) {
0473       // for edge pixels assign errors according to observed residual RMS
0474       if (edgex && !edgey) {
0475         xerr = xEdgeXError_ * micronsToCm;
0476         yerr = xEdgeYError_ * micronsToCm;
0477       } else if (!edgex && edgey) {
0478         xerr = yEdgeXError_ * micronsToCm;
0479         yerr = yEdgeYError_ * micronsToCm;
0480       } else if (edgex && edgey) {
0481         xerr = bothEdgeXError_ * micronsToCm;
0482         yerr = bothEdgeYError_ * micronsToCm;
0483       } else {
0484         throw cms::Exception(" PixelCPETemplateReco::localError: Something wrong with pixel edge flag !!!");
0485       }
0486 
0487       //cout << "xerr = " << xerr << endl;
0488       //cout << "yerr = " << yerr << endl;
0489     } else {
0490       // &&& need a class const
0491       //const float micronsToCm = 1.0e-4;
0492 
0493       xerr = theClusterParam.templSigmaX_ * micronsToCm;
0494       yerr = theClusterParam.templSigmaY_ * micronsToCm;
0495 
0496       //cout << "xerr = " << xerr << endl;
0497       //cout << "yerr = " << yerr << endl;
0498 
0499       // &&& should also check ierr (saved as class variable) and return
0500       // &&& nonsense (another class static) if the template fit failed.
0501     }
0502 
0503     if (theVerboseLevel > 9) {
0504       LogDebug("PixelCPETemplateReco") << " Sizex = " << theClusterParam.theCluster->sizeX()
0505                                        << " Sizey = " << theClusterParam.theCluster->sizeY() << " Edgex = " << edgex
0506                                        << " Edgey = " << edgey << " ErrX  = " << xerr << " ErrY  = " << yerr;
0507     }
0508 
0509   }  // else
0510 
0511   if (!(xerr > 0.0f))
0512     throw cms::Exception("PixelCPETemplateReco::localError")
0513         << "\nERROR: Negative pixel error xerr = " << xerr << "\n\n";
0514 
0515   if (!(yerr > 0.0f))
0516     throw cms::Exception("PixelCPETemplateReco::localError")
0517         << "\nERROR: Negative pixel error yerr = " << yerr << "\n\n";
0518 
0519   //cout << "Final errors set to: " << endl;
0520   //cout << "xerr = " << xerr << endl;
0521   //cout << "yerr = " << yerr << endl;
0522   //cout << "Out of PixelCPETemplateREco..........................................................................." << endl;
0523   //cout << endl;
0524 
0525   return LocalError(xerr * xerr, 0, yerr * yerr);
0526 }
0527 
0528 void PixelCPETemplateReco::fillPSetDescription(edm::ParameterSetDescription& desc) {
0529   desc.add<int>("barrelTemplateID", 0);
0530   desc.add<int>("forwardTemplateID", 0);
0531   desc.add<int>("directoryWithTemplates", 0);
0532   desc.add<int>("speed", -2);
0533   desc.add<bool>("UseClusterSplitter", false);
0534   desc.add<bool>("GoodEdgeAlgo", false);
0535 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!