From 8b3847344de5862308e21cd5c42e5ff3f1ca3f04 Mon Sep 17 00:00:00 2001 From: Matej Roguljic Date: Wed, 23 Apr 2025 15:19:32 +0200 Subject: [PATCH 1/5] change pixel cpe algo to handle broken clusters --- .../interface/SiPixelTemplate.h | 3 + .../SiPixelTransient/interface/SiPixelUtils.h | 32 +++-- .../SiPixelTransient/src/SiPixelTemplate.cc | 40 ++++-- .../SiPixelTransient/src/SiPixelUtils.cc | 53 ++++--- .../python/siPixelGoodEdgeAlgo_cff.py | 4 + .../interface/PixelCPEClusterRepair.h | 1 + .../interface/PixelCPEFastParamsHost.h | 4 +- .../interface/PixelCPEGeneric.h | 2 + .../interface/PixelCPETemplateReco.h | 1 + .../interface/SiPixelTemplateReco.h | 21 ++- .../interface/pixelCPEforDevice.h | 78 ++++++++++- .../PixelCPEFastParamsESProducerAlpaka.cc | 6 +- .../python/PixelCPEClusterRepair_cfi.py | 5 + .../python/PixelCPEESProducers_cff.py | 6 + .../python/PixelCPEGeneric_cfi.py | 4 + .../python/PixelCPETemplateReco_cfi.py | 3 + .../src/PixelCPEClusterRepair.cc | 9 +- .../src/PixelCPEFastParamsHost.cc | 7 +- .../SiPixelRecHits/src/PixelCPEGeneric.cc | 13 +- .../src/PixelCPETemplateReco.cc | 5 +- .../SiPixelRecHits/src/SiPixelTemplateReco.cc | 130 ++++++++++++++++-- 21 files changed, 364 insertions(+), 63 deletions(-) create mode 100644 Configuration/ProcessModifiers/python/siPixelGoodEdgeAlgo_cff.py diff --git a/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h b/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h index 2b9990fd83b5d..6587f2e1dbd5b 100644 --- a/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h +++ b/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h @@ -285,6 +285,9 @@ class SiPixelTemplate { // initialize the rest; static void postInit(std::vector& thePixelTemp_); + // Interpolate with y Gaussian Parameter interpolation to be used with goodEdge reconstruction algorithm + bool interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx, bool goodEdgeAlgo); + // Interpolate input alpha and beta angles to produce a working template for each individual hit. bool interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx); diff --git a/CondFormats/SiPixelTransient/interface/SiPixelUtils.h b/CondFormats/SiPixelTransient/interface/SiPixelUtils.h index 4c50e85b4d719..6c97a47195220 100644 --- a/CondFormats/SiPixelTransient/interface/SiPixelUtils.h +++ b/CondFormats/SiPixelTransient/interface/SiPixelUtils.h @@ -2,21 +2,23 @@ #define CondFormats_SiPixelTransient_SiPixelUtils_h namespace siPixelUtils { - float generic_position_formula(int size, //!< Size of this projection. - int q_f, //!< Charge in the first pixel. - int q_l, //!< Charge in the last pixel. - float upper_edge_first_pix, //!< As the name says. - float lower_edge_last_pix, //!< As the name says. - float lorentz_shift, //!< L-width - float theThickness, //detector thickness - float cot_angle, //!< cot of alpha_ or beta_ - float pitch, //!< thePitchX or thePitchY - float pitchfraction_first, //!< true if the first is big - float pitchfraction_last, //!< true if the last is big - float eff_charge_cut_low, //!< Use edge if > W_eff (in pix) &&& - float eff_charge_cut_high, //!< Use edge if < W_eff (in pix) &&& - float size_cut //!< Use edge when size == cuts - ); + float generic_position_formula( + int size, //!< Size of this projection. + int q_f, //!< Charge in the first pixel. + int q_l, //!< Charge in the last pixel. + float upper_edge_first_pix, //!< As the name says. + float lower_edge_last_pix, //!< As the name says. + float lorentz_shift, //!< L-width + float theThickness, //detector thickness + float cot_angle, //!< cot of alpha_ or beta_ + float pitch, //!< thePitchX or thePitchY + float pitchfraction_first, //!< true if the first is big + float pitchfraction_last, //!< true if the last is big + float eff_charge_cut_low, //!< Use edge if > W_eff (in pix) &&& + float eff_charge_cut_high, //!< Use edge if < W_eff (in pix) &&& + float size_cut, //!< Use edge when size == cuts + float delta_length_cut = 2., //!< if charge len - cls size > this (in pix), use one-sided reco + bool goodEdgeAlgo = false); } // namespace siPixelUtils #endif diff --git a/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc b/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc index 0e3c83505634a..b6037b8f897f6 100644 --- a/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc +++ b/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc @@ -82,6 +82,7 @@ // V10.21 - Address runtime issues in pushfile() for gcc 7.X due to using tempfile as char string + misc. cleanup [Petar] // V10.22 - Move templateStore to the heap, fix variable name in pushfile() [Petar] // V10.24 - Add sideload() + associated gymnastics [Petar and Oz] +// V10.25 - Restore y-residual Gaussian parameters [Morris] // Created by Morris Swartz on 10/27/06. // @@ -1324,15 +1325,18 @@ void SiPixelTemplate::postInit(std::vector& thePixelTemp_) //! \param locBx - (input) the sign of this quantity is used to determine whether to flip cot(alpha/beta)<0 quantities from cot(alpha/beta)>0 (FPix only) //! for Phase 1 FPix IP-related tracks, locBx/locBz > 0 for cot(alpha) > 0 and locBx/locBz < 0 for cot(alpha) < 0 //! for Phase 1 FPix IP-related tracks, locBx > 0 for cot(beta) > 0 and locBx < 0 for cot(beta) < 0 +//! \param goodEdgeAlgo - (input) Flag to turn on the y Gaussian Parameter interpolation to be used with goodEdge reconstruction algorithm // ************************************************************************************************************ -bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx) { +bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx, bool goodEdgeAlgo) { // Interpolate for a new set of track angles +#ifndef SI_PIXEL_TEMPLATE_STANDALONE //check for nan's if (!edm::isFinite(cotalpha) || !edm::isFinite(cotbeta)) { success_ = false; return success_; } +#endif // Local variables int i, j; @@ -1553,16 +1557,19 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float l } for (i = 0; i < 4; ++i) { - yavg_[i] = (1.f - yratio_) * thePixelTemp_[index_id_].enty[ilow].yavg[i] + - yratio_ * thePixelTemp_[index_id_].enty[ihigh].yavg[i]; - if (flip_y_) { - yavg_[i] = -yavg_[i]; - } yavg_[i] = (1.f - yratio_) * enty0_->yavg[i] + yratio_ * enty1_->yavg[i]; if (flip_y_) { yavg_[i] = -yavg_[i]; } yrms_[i] = (1.f - yratio_) * enty0_->yrms[i] + yratio_ * enty1_->yrms[i]; + + if (goodEdgeAlgo) { // restore y Gaussian Parameter interpolation + ygx0_[i] = (1.f - yratio_) * enty0_->ygx0[i] + yratio_ * enty1_->ygx0[i]; + if (flip_y_) { + ygx0_[i] = -ygx0_[i]; + } + ygsig_[i] = (1.f - yratio_) * enty0_->ygsig[i] + yratio_ * enty1_->ygsig[i]; + } //if(goodEdgeAlgo) chi2yavg_[i] = (1.f - yratio_) * enty0_->chi2yavg[i] + yratio_ * enty1_->chi2yavg[i]; chi2ymin_[i] = (1.f - yratio_) * enty0_->chi2ymin[i] + yratio_ * enty1_->chi2ymin[i]; chi2xavg[i] = (1.f - yratio_) * enty0_->chi2xavg[i] + yratio_ * enty1_->chi2xavg[i]; @@ -1869,6 +1876,23 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float l return success_; } // interpolate +// ************************************************************************************************************ +//! Interpolate input alpha and beta angles to produce a working template for each individual hit. +//! \param id - (input) index of the template to use +//! \param cotalpha - (input) the cotangent of the alpha track angle (see CMS IN 2004/014) +//! \param cotbeta - (input) the cotangent of the beta track angle (see CMS IN 2004/014) +//! \param locBz - (input) the sign of this quantity is used to determine whether to flip cot(beta)<0 quantities from cot(beta)>0 (FPix only) +//! for Phase 0 FPix IP-related tracks, locBz < 0 for cot(beta) > 0 and locBz > 0 for cot(beta) < 0 +//! for Phase 1 FPix IP-related tracks, see next comment +//! \param locBx - (input) the sign of this quantity is used to determine whether to flip cot(alpha/beta)<0 quantities from cot(alpha/beta)>0 (FPix only) +//! for Phase 1 FPix IP-related tracks, locBx/locBz > 0 for cot(alpha) > 0 and locBx/locBz < 0 for cot(alpha) < 0 +//! for Phase 1 FPix IP-related tracks, locBx > 0 for cot(beta) > 0 and locBx < 0 for cot(beta) < 0 +// ************************************************************************************************************ +bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx) { + // Interpolate for a new set of track angles, but without y Gaussian Parameter interpolation to be used with goodEdge reconstruction algorithm + return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx, false); +} + // ************************************************************************************************************ //! Interpolate input alpha and beta angles to produce a working template for each individual hit. //! \param id - (input) index of the template to use @@ -1888,7 +1912,7 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta) { if (cotalpha < 0.f) { locBz = -locBx; } - return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx); + return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx, false); } // ************************************************************************************************************ @@ -1903,7 +1927,7 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float l // Local variables float locBx = 1.f; - return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx); + return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx, false); } // ************************************************************************************************************************************* diff --git a/CondFormats/SiPixelTransient/src/SiPixelUtils.cc b/CondFormats/SiPixelTransient/src/SiPixelUtils.cc index 9678753e4d0b7..06da0087cbc49 100644 --- a/CondFormats/SiPixelTransient/src/SiPixelUtils.cc +++ b/CondFormats/SiPixelTransient/src/SiPixelUtils.cc @@ -13,21 +13,23 @@ namespace siPixelUtils { //! X and Y, in the interest of the simplicity of the code, all parameters //! are passed by the caller. //----------------------------------------------------------------------------- - float generic_position_formula(int size, //!< Size of this projection. - int q_f, //!< Charge in the first pixel. - int q_l, //!< Charge in the last pixel. - float upper_edge_first_pix, //!< As the name says. - float lower_edge_last_pix, //!< As the name says. - float lorentz_shift, //!< L-shift at half thickness - float theThickness, //detector thickness - float cot_angle, //!< cot of alpha_ or beta_ - float pitch, //!< thePitchX or thePitchY - float pitchfraction_first, - float pitchfraction_last, - float eff_charge_cut_low, //!< Use edge if > w_eff &&& - float eff_charge_cut_high, //!< Use edge if < w_eff &&& - float size_cut //!< Use edge when size == cuts - ) { + float generic_position_formula( + int size, //!< Size of this projection. + int q_f, //!< Charge in the first pixel. + int q_l, //!< Charge in the last pixel. + float upper_edge_first_pix, //!< As the name says. + float lower_edge_last_pix, //!< As the name says. + float lorentz_shift, //!< L-shift at half thickness + float theThickness, //detector thickness + float cot_angle, //!< cot of alpha_ or beta_ + float pitch, //!< thePitchX or thePitchY + float pitchfraction_first, + float pitchfraction_last, + float eff_charge_cut_low, //!< Use edge if > w_eff &&& + float eff_charge_cut_high, //!< Use edge if < w_eff &&& + float size_cut, //!< Use edge when size == cuts + float delta_length_cut, //!< if charge len - cls size > this (in pix), use one-sided reco + bool goodEdgeAlgo) { float geom_center = 0.5f * (upper_edge_first_pix + lower_edge_last_pix); //--- The case of only one pixel in this projection is separate. Note that @@ -56,11 +58,26 @@ namespace siPixelUtils { //--- it with an *average* effective charge width, which is the average //--- length of the edge pixels. // - // bool usedEdgeAlgo = false; if ((size >= size_cut) || ((w_eff / pitch < eff_charge_cut_low) | (w_eff / pitch > eff_charge_cut_high))) { w_eff = pitch * 0.5f * sum_of_edge; // ave. length of edge pixels (first+last) (cm) - // usedEdgeAlgo = true; - } + + if (goodEdgeAlgo) { + float delta = w_eff - 0.5 * sum_of_edge * pitch; + if (delta / pitch > delta_length_cut) { + // define the centers of the first last last pixel coordinates + float x1 = upper_edge_first_pix - 0.5 * pitchfraction_first * pitch; + float x2 = lower_edge_last_pix + 0.5 * pitchfraction_last * pitch; + // observed cluster is much shorter than expected, use one-sided reco + if (w_pred > 0.f) { + float hit_pos = x1 + 0.5 * w_pred; + return hit_pos; + } else { + float hit_pos = x2 + 0.5 * w_pred; + return hit_pos; + } + } //if (delta / pitch > delta_length_cut) + } //if(goodEdgeAlgo) + } //if(size >= size_cut) ||... //--- Finally, compute the position in this projection float q_diff = q_l - q_f; diff --git a/Configuration/ProcessModifiers/python/siPixelGoodEdgeAlgo_cff.py b/Configuration/ProcessModifiers/python/siPixelGoodEdgeAlgo_cff.py new file mode 100644 index 0000000000000..c37354bbc59b7 --- /dev/null +++ b/Configuration/ProcessModifiers/python/siPixelGoodEdgeAlgo_cff.py @@ -0,0 +1,4 @@ +import FWCore.ParameterSet.Config as cms + +# This modifier enables the good edge algorithm in pixel hit reconstruction that handles broken/truncated pixel cluster caused by radiation damage +siPixelGoodEdgeAlgo = cms.Modifier() diff --git a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEClusterRepair.h b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEClusterRepair.h index 0e44e3907b9dd..5843e75cfc78d 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEClusterRepair.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEClusterRepair.h @@ -102,6 +102,7 @@ class PixelCPEClusterRepair : public PixelCPEBase { std::vector thePixelTemp2D_; int speed_; + bool goodEdgeAlgo_; bool UseClusterSplitter_; diff --git a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEFastParamsHost.h b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEFastParamsHost.h index 7d57c46dd7a13..eba4942684bc9 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEFastParamsHost.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEFastParamsHost.h @@ -26,7 +26,8 @@ class PixelCPEFastParamsHost : public PixelCPEGenericBase { const TrackerTopology& ttopo, const SiPixelLorentzAngle* lorentzAngle, const SiPixelGenErrorDBObject* genErrorDBObject, - const SiPixelLorentzAngle* lorentzAngleWidth); + const SiPixelLorentzAngle* lorentzAngleWidth, + const bool goodEdgeAlgo); // non-copyable PixelCPEFastParamsHost(PixelCPEFastParamsHost const&) = delete; @@ -61,6 +62,7 @@ class PixelCPEFastParamsHost : public PixelCPEGenericBase { void fillParamsForDevice(); Buffer buffer_; + bool goodEdgeAlgo_; }; #endif // RecoLocalTracker_SiPixelRecHits_interface_PixelCPEFastParamsHost_h diff --git a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEGeneric.h b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEGeneric.h index bad29fb56aefc..11ebf7efd7a1c 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEGeneric.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPEGeneric.h @@ -81,12 +81,14 @@ class PixelCPEGeneric : public PixelCPEGenericBase { float the_eff_charge_cut_highY; float the_size_cutX; float the_size_cutY; + float delta_length_cut; bool inflate_errors; bool inflate_all_errors_no_trk_angle; bool DoCosmics_; bool IrradiationBiasCorrection_; + bool goodEdgeAlgo_; bool isPhase2_; bool NoTemplateErrorsWhenNoTrkAngles_; diff --git a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPETemplateReco.h b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPETemplateReco.h index dc92fae67d4c6..b4f034f9925a2 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/PixelCPETemplateReco.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/PixelCPETemplateReco.h @@ -79,6 +79,7 @@ class PixelCPETemplateReco : public PixelCPEBase { const std::vector *thePixelTemp_; int speed_; + bool goodEdgeAlgo_; bool UseClusterSplitter_; diff --git a/RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateReco.h b/RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateReco.h index 51e434e812a60..0cd9e6bbb6906 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateReco.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/SiPixelTemplateReco.h @@ -99,7 +99,26 @@ namespace SiPixelTemplateReco { std::vector >& zeropix, float& probQ, int& nypix, - int& nxpix); + int& nxpix, + bool goodEdgeAlgo); + + int PixelTempReco1D(int id, + float cotalpha, + float cotbeta, + float locBz, + float locBx, + ClusMatrix& cluster, + SiPixelTemplate& templ, + float& yrec, + float& sigmay, + float& proby, + float& xrec, + float& sigmax, + float& probx, + int& qbin, + int speed, + float& probQ, + bool goodEdgeAlgo); int PixelTempReco1D(int id, float cotalpha, diff --git a/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h b/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h index a97add7edb7b3..0d8b6d862fcfa 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h @@ -64,6 +64,8 @@ namespace pixelCPEforDevice { uint16_t maxModuleStride; uint8_t numberOfLaddersInBarrel; + + bool goodEdgeAlgo_; }; struct DetParams { @@ -122,6 +124,53 @@ namespace pixelCPEforDevice { cotbeta = gvy * gvz; } + constexpr inline float check_for_short(uint16_t upper_edge_first_pix, //!< As the name says. + uint16_t lower_edge_last_pix, //!< As the name says. + float lorentz_shift, //!< L-shift at half thickness + float theThickness, //detector thickness + float cot_angle, //!< cot of alpha_ or beta_ + float pitch, //!< thePitchX or thePitchY + bool first_is_big, //!< true if the first is big + bool last_is_big) //!< true if the last is big + { + //checks if the observed cluster is much shorter than expected + //if yes, returns the position predicted using a one-sided reconstruction + //if not, returns a dummy value -999. + float delta_length_cut = 2.; //to be fine tuned later + float w_inner = pitch * (lower_edge_last_pix - upper_edge_first_pix); + float w_pred = theThickness * cot_angle - lorentz_shift; + + float w_eff = std::abs(w_pred) - w_inner; + + float pitchfraction_first = 1.0; + float pitchfraction_last = 1.0; + + if (first_is_big) + pitchfraction_first += 1.0f; + if (last_is_big) + pitchfraction_last += 1.0f; + + float sum_of_edge = pitchfraction_first + pitchfraction_last; + float delta = w_eff - 0.5 * sum_of_edge * pitch; + + if (delta / pitch > delta_length_cut) { + // define the centers of the first and last pixel coordinates + float x1 = pitch * (upper_edge_first_pix - 0.5 * pitchfraction_first); + float x2 = pitch * (lower_edge_last_pix + 0.5 * pitchfraction_last); + // observed cluster is much shorter than expected, use one-sided reco + if (w_pred > 0.f) { + float hit_pos = x1 + 0.5 * w_pred; + return hit_pos; + } else { + float hit_pos = x2 + 0.5 * w_pred; + return hit_pos; + } + } //if(delta/pitch > delta_length_cut) + else { + return -999.; + } + } + constexpr inline float correction(int sizeM1, int q_f, //!< Charge in the first pixel. int q_l, //!< Charge in the last pixel. @@ -250,7 +299,34 @@ namespace pixelCPEforDevice { computeAnglesFromDet(detParams, xPos, yPos, cotalpha, cotbeta); auto thickness = detParams.isBarrel ? comParams.theThicknessB : comParams.theThicknessE; - + if (comParams.goodEdgeAlgo_) { + float xPos_alt = check_for_short(llxl, + urxl, + detParams.chargeWidthX, + thickness, + cotalpha, + detParams.thePitchX, + TrackerTraits::isBigPixX(cp.minRow[ic]), + TrackerTraits::isBigPixX(cp.maxRow[ic])); + + float yPos_alt = check_for_short(llyl, + uryl, + detParams.chargeWidthY, + thickness, + cotbeta, + detParams.thePitchY, + TrackerTraits::isBigPixY(cp.minCol[ic]), + TrackerTraits::isBigPixY(cp.maxCol[ic])); + + if (xPos_alt != -999.) { + // observed cluster is much shorter than expected, use one-sided reco + xPos = xPos_alt; + } + if (yPos_alt != -999.) { + // observed cluster is much shorter than expected, use one-sided reco + yPos = yPos_alt; + } + } //if(comParams.goodEdgeAlgo_) auto xCorr = correction(cp.maxRow[ic] - cp.minRow[ic], cp.q_f_X[ic], cp.q_l_X[ic], diff --git a/RecoLocalTracker/SiPixelRecHits/plugins/alpaka/PixelCPEFastParamsESProducerAlpaka.cc b/RecoLocalTracker/SiPixelRecHits/plugins/alpaka/PixelCPEFastParamsESProducerAlpaka.cc index 2872cedd14aeb..efcb88557424c 100644 --- a/RecoLocalTracker/SiPixelRecHits/plugins/alpaka/PixelCPEFastParamsESProducerAlpaka.cc +++ b/RecoLocalTracker/SiPixelRecHits/plugins/alpaka/PixelCPEFastParamsESProducerAlpaka.cc @@ -42,6 +42,7 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { edm::ParameterSet pset_; bool useErrorsFromTemplates_; + bool goodEdgeAlgo_; }; using namespace edm; @@ -52,6 +53,7 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { auto const& myname = p.getParameter("ComponentName"); auto const& magname = p.getParameter("MagneticFieldRecord"); useErrorsFromTemplates_ = p.getParameter("UseErrorsFromTemplates"); + goodEdgeAlgo_ = p.getParameter("GoodEdgeAlgo"); auto cc = setWhatProduced(this, myname); magfieldToken_ = cc.consumes(magname); @@ -84,7 +86,8 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { iRecord.get(hTTToken_), &iRecord.get(lorentzAngleToken_), genErrorDBObjectProduct, - lorentzAngleWidthProduct); + lorentzAngleWidthProduct, + goodEdgeAlgo_); } template @@ -103,6 +106,7 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { desc.add("EdgeClusterErrorY", 85.0); desc.add("UseErrorsFromTemplates", true); desc.add("TruncatePixelCharge", true); + desc.add("GoodEdgeAlgo", false); std::string name = "PixelCPEFastParams"; name += TrackerTraits::nameModifier; diff --git a/RecoLocalTracker/SiPixelRecHits/python/PixelCPEClusterRepair_cfi.py b/RecoLocalTracker/SiPixelRecHits/python/PixelCPEClusterRepair_cfi.py index 7823d45144ae2..0eb7dbbb0d371 100644 --- a/RecoLocalTracker/SiPixelRecHits/python/PixelCPEClusterRepair_cfi.py +++ b/RecoLocalTracker/SiPixelRecHits/python/PixelCPEClusterRepair_cfi.py @@ -6,3 +6,8 @@ ComponentName = "PixelCPEClusterRepairWithoutProbQ", speed = 0 ) + +# Enable the good edge algorithm in pixel hit reconstruction that handles broken/truncated pixel cluster caused by radiation damage +from Configuration.ProcessModifiers.siPixelGoodEdgeAlgo_cff import siPixelGoodEdgeAlgo +siPixelGoodEdgeAlgo.toModify(templates2, GoodEdgeAlgo = True) +siPixelGoodEdgeAlgo.toModify(templates2_speed0, GoodEdgeAlgo = True) diff --git a/RecoLocalTracker/SiPixelRecHits/python/PixelCPEESProducers_cff.py b/RecoLocalTracker/SiPixelRecHits/python/PixelCPEESProducers_cff.py index d0532d33ff260..781e92f7c2ca8 100644 --- a/RecoLocalTracker/SiPixelRecHits/python/PixelCPEESProducers_cff.py +++ b/RecoLocalTracker/SiPixelRecHits/python/PixelCPEESProducers_cff.py @@ -1,5 +1,6 @@ import FWCore.ParameterSet.Config as cms from Configuration.ProcessModifiers.alpaka_cff import alpaka +from Configuration.ProcessModifiers.siPixelGoodEdgeAlgo_cff import siPixelGoodEdgeAlgo # # Load all Pixel Cluster Position Estimator ESProducers @@ -24,3 +25,8 @@ def _addProcessCPEsAlpaka(process): modifyConfigurationForAlpakaCPEs_ = alpaka.makeProcessModifier(_addProcessCPEsAlpaka) +def _enableGoodEdgeAlgoAlpaka(process): + _addProcessCPEsAlpaka(process) + process.pixelCPEFastParamsESProducerAlpakaPhase1.GoodEdgeAlgo = True + +enableGoodEdgeAlgoForAlpakaCPEs_ = (alpaka & siPixelGoodEdgeAlgo).makeProcessModifier(_enableGoodEdgeAlgoAlpaka) diff --git a/RecoLocalTracker/SiPixelRecHits/python/PixelCPEGeneric_cfi.py b/RecoLocalTracker/SiPixelRecHits/python/PixelCPEGeneric_cfi.py index 370eabae2b06d..21bab181a2d59 100644 --- a/RecoLocalTracker/SiPixelRecHits/python/PixelCPEGeneric_cfi.py +++ b/RecoLocalTracker/SiPixelRecHits/python/PixelCPEGeneric_cfi.py @@ -8,6 +8,10 @@ from Configuration.Eras.Modifier_run3_common_cff import run3_common run3_common.toModify(PixelCPEGenericESProducer, IrradiationBiasCorrection = True) +# Enable the good edge algorithm in pixel hit reconstruction that handles broken/truncated pixel cluster caused by radiation damage +from Configuration.ProcessModifiers.siPixelGoodEdgeAlgo_cff import siPixelGoodEdgeAlgo +siPixelGoodEdgeAlgo.toModify(PixelCPEGenericESProducer, GoodEdgeAlgo = True) + # customize the Pixel CPE generic producer for phase2 from Configuration.Eras.Modifier_phase2_tracker_cff import phase2_tracker phase2_tracker.toModify(PixelCPEGenericESProducer, diff --git a/RecoLocalTracker/SiPixelRecHits/python/PixelCPETemplateReco_cfi.py b/RecoLocalTracker/SiPixelRecHits/python/PixelCPETemplateReco_cfi.py index 3c74bbe3c75ff..b2287a01d3875 100644 --- a/RecoLocalTracker/SiPixelRecHits/python/PixelCPETemplateReco_cfi.py +++ b/RecoLocalTracker/SiPixelRecHits/python/PixelCPETemplateReco_cfi.py @@ -3,3 +3,6 @@ from RecoLocalTracker.SiPixelRecHits._templates_default_cfi import _templates_default templates = _templates_default.clone() +# Enable the good edge algorithm in pixel hit reconstruction that handles broken/truncated pixel cluster caused by radiation damage +from Configuration.ProcessModifiers.siPixelGoodEdgeAlgo_cff import siPixelGoodEdgeAlgo +siPixelGoodEdgeAlgo.toModify(templates, GoodEdgeAlgo = True) diff --git a/RecoLocalTracker/SiPixelRecHits/src/PixelCPEClusterRepair.cc b/RecoLocalTracker/SiPixelRecHits/src/PixelCPEClusterRepair.cc index d2d8502daeea4..24fdb00d9edb3 100644 --- a/RecoLocalTracker/SiPixelRecHits/src/PixelCPEClusterRepair.cc +++ b/RecoLocalTracker/SiPixelRecHits/src/PixelCPEClusterRepair.cc @@ -75,7 +75,9 @@ PixelCPEClusterRepair::PixelCPEClusterRepair(edm::ParameterSet const& conf, } speed_ = conf.getParameter("speed"); + goodEdgeAlgo_ = conf.getParameter("GoodEdgeAlgo"); LogDebug("PixelCPEClusterRepair::PixelCPEClusterRepair:") << "Template speed = " << speed_ << "\n"; + LogDebug("PixelCPEClusterRepair::PixelCPEClusterRepair:") << "GoodEdgeAlgo = " << goodEdgeAlgo_ << "\n"; // this returns the magnetic field value in kgauss (1T = 10 kgauss) int theMagField = mag->nominalValue(); @@ -355,7 +357,8 @@ void PixelCPEClusterRepair::callTempReco1D(DetParam const& theDetParam, zeropix, theClusterParam.probabilityQ_, nypix, - nxpix); + nxpix, + goodEdgeAlgo_); // ****************************************************************** //--- Check exit status @@ -549,7 +552,8 @@ void PixelCPEClusterRepair::checkRecommend2D(DetParam const& theDetParam, } // The 1d pixel template SiPixelTemplate templ(*thePixelTemp_); - if (!templ.interpolate(ID, theClusterParam.cotalpha, theClusterParam.cotbeta, theDetParam.bz, theDetParam.bx)) { + if (!templ.interpolate( + ID, theClusterParam.cotalpha, theClusterParam.cotbeta, theDetParam.bz, theDetParam.bx, goodEdgeAlgo_)) { //error setting up template, return false theClusterParam.recommended2D_ = false; return; @@ -718,6 +722,7 @@ void PixelCPEClusterRepair::fillPSetDescription(edm::ParameterSetDescription& de desc.add("forwardTemplateID", 0); desc.add("directoryWithTemplates", 0); desc.add("speed", -2); + desc.add("GoodEdgeAlgo", false); desc.add("UseClusterSplitter", false); desc.add("MaxSizeMismatchInY", 0.3); desc.add("MinChargeRatio", 0.8); diff --git a/RecoLocalTracker/SiPixelRecHits/src/PixelCPEFastParamsHost.cc b/RecoLocalTracker/SiPixelRecHits/src/PixelCPEFastParamsHost.cc index b290aabc194d1..5dfd279ae89f2 100644 --- a/RecoLocalTracker/SiPixelRecHits/src/PixelCPEFastParamsHost.cc +++ b/RecoLocalTracker/SiPixelRecHits/src/PixelCPEFastParamsHost.cc @@ -20,7 +20,8 @@ PixelCPEFastParamsHost::PixelCPEFastParamsHost(edm::ParameterSet const TrackerTopology& ttopo, const SiPixelLorentzAngle* lorentzAngle, const SiPixelGenErrorDBObject* genErrorDBObject, - const SiPixelLorentzAngle* lorentzAngleWidth) + const SiPixelLorentzAngle* lorentzAngleWidth, + const bool goodEdgeAlgo) : PixelCPEGenericBase(conf, mag, geom, ttopo, lorentzAngle, genErrorDBObject, lorentzAngleWidth), buffer_(cms::alpakatools::make_host_buffer>()) { // Use errors from templates or from GenError @@ -30,7 +31,7 @@ PixelCPEFastParamsHost::PixelCPEFastParamsHost(edm::ParameterSet << "ERROR: GenErrors not filled correctly. Check the sqlite file. Using SiPixelTemplateDBObject version " << (*genErrorDBObject_).version(); } - + goodEdgeAlgo_ = goodEdgeAlgo; fillParamsForDevice(); } @@ -43,7 +44,7 @@ void PixelCPEFastParamsHost::fillParamsForDevice() { buffer_->commonParams().theThicknessB = m_DetParams.front().theThickness; buffer_->commonParams().theThicknessE = m_DetParams.back().theThickness; buffer_->commonParams().numberOfLaddersInBarrel = TrackerTraits::numberOfLaddersInBarrel; - + buffer_->commonParams().goodEdgeAlgo_ = goodEdgeAlgo_; LogDebug("PixelCPEFastParamsHost") << "thickness " << buffer_->commonParams().theThicknessB << ' ' << buffer_->commonParams().theThicknessE; diff --git a/RecoLocalTracker/SiPixelRecHits/src/PixelCPEGeneric.cc b/RecoLocalTracker/SiPixelRecHits/src/PixelCPEGeneric.cc index d3ddb480b2c69..7b17b5340505b 100644 --- a/RecoLocalTracker/SiPixelRecHits/src/PixelCPEGeneric.cc +++ b/RecoLocalTracker/SiPixelRecHits/src/PixelCPEGeneric.cc @@ -44,6 +44,7 @@ PixelCPEGeneric::PixelCPEGeneric(edm::ParameterSet const& conf, the_eff_charge_cut_highY = conf.getParameter("eff_charge_cut_highY"); the_size_cutX = conf.getParameter("size_cutX"); the_size_cutY = conf.getParameter("size_cutY"); + delta_length_cut = 2.0f; //seems to be a good value, if needed, will make this externally settable // Externally settable flags to inflate errors inflate_errors = conf.getParameter("inflate_errors"); @@ -51,6 +52,8 @@ PixelCPEGeneric::PixelCPEGeneric(edm::ParameterSet const& conf, NoTemplateErrorsWhenNoTrkAngles_ = conf.getParameter("NoTemplateErrorsWhenNoTrkAngles"); IrradiationBiasCorrection_ = conf.getParameter("IrradiationBiasCorrection"); + goodEdgeAlgo_ = conf.getParameter("GoodEdgeAlgo"); //use only one edge reco if cluster seems broken + DoCosmics_ = conf.getParameter("DoCosmics"); isPhase2_ = conf.getParameter("isPhase2"); @@ -92,6 +95,7 @@ PixelCPEGeneric::PixelCPEGeneric(edm::ParameterSet const& conf, cout << "(int)useErrorsFromTemplates_ = " << (int)useErrorsFromTemplates_ << endl; cout << "truncatePixelCharge_ = " << (int)truncatePixelCharge_ << endl; cout << "IrradiationBiasCorrection_ = " << (int)IrradiationBiasCorrection_ << endl; + cout << "goodEdgeAlgo_ = " << (int)goodEdgeAlgo_ << endl; cout << "(int)DoCosmics_ = " << (int)DoCosmics_ << endl; cout << "(int)LoadTemplatesFromDB_ = " << (int)LoadTemplatesFromDB_ << endl; #endif @@ -266,7 +270,9 @@ LocalPoint PixelCPEGeneric::localPosition(DetParam const& theDetParam, ClusterPa theDetParam.theTopol->pixelFractionInX(theClusterParam.theCluster->maxPixelRow()), the_eff_charge_cut_lowX, the_eff_charge_cut_highX, - the_size_cutX); // cut for eff charge width &&& + the_size_cutX, // cut for eff charge width &&& + delta_length_cut, + goodEdgeAlgo_); // apply the lorentz offset correction xPos = xPos + shiftX; @@ -290,7 +296,9 @@ LocalPoint PixelCPEGeneric::localPosition(DetParam const& theDetParam, ClusterPa theDetParam.theTopol->pixelFractionInY(theClusterParam.theCluster->maxPixelCol()), the_eff_charge_cut_lowY, the_eff_charge_cut_highY, - the_size_cutY); // cut for eff charge width &&& + the_size_cutY, // cut for eff charge width &&& + delta_length_cut, + goodEdgeAlgo_); // apply the lorentz offset correction yPos = yPos + shiftY; @@ -447,6 +455,7 @@ void PixelCPEGeneric::fillPSetDescription(edm::ParameterSetDescription& desc) { desc.add("UseErrorsFromTemplates", true); desc.add("TruncatePixelCharge", true); desc.add("IrradiationBiasCorrection", false); + desc.add("GoodEdgeAlgo", false); desc.add("DoCosmics", false); desc.add("isPhase2", false); desc.add("SmallPitch", false); diff --git a/RecoLocalTracker/SiPixelRecHits/src/PixelCPETemplateReco.cc b/RecoLocalTracker/SiPixelRecHits/src/PixelCPETemplateReco.cc index beb2884ecd288..de04276201101 100644 --- a/RecoLocalTracker/SiPixelRecHits/src/PixelCPETemplateReco.cc +++ b/RecoLocalTracker/SiPixelRecHits/src/PixelCPETemplateReco.cc @@ -81,6 +81,7 @@ PixelCPETemplateReco::PixelCPETemplateReco(edm::ParameterSet const& conf, << " not loaded correctly from text file. Reconstruction will fail.\n\n"; } + goodEdgeAlgo_ = conf.getParameter("GoodEdgeAlgo"); speed_ = conf.getParameter("speed"); LogDebug("PixelCPETemplateReco::PixelCPETemplateReco:") << "Template speed = " << speed_ << "\n"; @@ -249,7 +250,8 @@ LocalPoint PixelCPETemplateReco::localPosition(DetParam const& theDetParam, Clus theClusterParam.templProbX_, theClusterParam.templQbin_, speed_, - theClusterParam.templProbQ_); + theClusterParam.templProbQ_, + goodEdgeAlgo_); // ****************************************************************** @@ -529,4 +531,5 @@ void PixelCPETemplateReco::fillPSetDescription(edm::ParameterSetDescription& des desc.add("directoryWithTemplates", 0); desc.add("speed", -2); desc.add("UseClusterSplitter", false); + desc.add("GoodEdgeAlgo", false); } diff --git a/RecoLocalTracker/SiPixelRecHits/src/SiPixelTemplateReco.cc b/RecoLocalTracker/SiPixelRecHits/src/SiPixelTemplateReco.cc index 78d0a97ee20de..4eeeff19f91b6 100644 --- a/RecoLocalTracker/SiPixelRecHits/src/SiPixelTemplateReco.cc +++ b/RecoLocalTracker/SiPixelRecHits/src/SiPixelTemplateReco.cc @@ -45,6 +45,8 @@ // V10.00 - Use new template object to reco Phase 1 FPix hits // V10.01 - Fix memory overwriting bug // V10.10 - Change VVIObjF so it only reads kappa +// V10.30 - Use "good" end to center the y-projections [for high eta radiation damage], +// - use Gaussian fit paramters for y-bias and y-errors, remove chi2min_y warning // // // Created by Morris Swartz on 10/27/06. @@ -131,6 +133,7 @@ using namespace SiPixelTemplateReco; //! \param probQ - (output) the Vavilov-distribution-based cluster charge probability //! \param nypix - (output) the projected y-size of the cluster //! \param nxpix - (output) the projected x-size of the cluster +//! \param goodEdgeAlgo - (input) bool to indicate whether to use centering based on the good end of the cluster (compatible only with templates derived using the same) // ************************************************************************************************************************************* int SiPixelTemplateReco::PixelTempReco1D(int id, float cotalpha, @@ -151,7 +154,8 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, std::vector >& zeropix, float& probQ, int& nypix, - int& nxpix) + int& nxpix, + bool goodEdgeAlgo) { // Local variables @@ -186,7 +190,7 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, // check to see of the track direction is in the physical range of the loaded template if (id >= 0) { //if id < 0 bypass interpolation (used in calibration) - if (!templ.interpolate(id, cotalpha, cotbeta, locBz, locBx)) { + if (!templ.interpolate(id, cotalpha, cotbeta, locBz, locBx, goodEdgeAlgo)) { if (theVerboseLevel > 2) { LOGDEBUG("SiPixelTemplateReco") << "input cluster direction cot(alpha) = " << cotalpha << ", cot(beta) = " << cotbeta << ", local B_z = " << locBz @@ -463,8 +467,18 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, // next, center the cluster on template center if necessary - midpix = (fypix + lypix) / 2; - shifty = templ.cytemp() - midpix; + if (goodEdgeAlgo) { + if (cotbeta >= 0) { + midpix = (int)(fypix + 0.5 * cotbeta * templ.zsize() / ysize); + } else { + midpix = (int)(lypix + 0.5 * cotbeta * templ.zsize() / ysize); + } + shifty = BHY - midpix; + } //if(goodEdgeAlgo) + else { + midpix = (fypix + lypix) / 2; + shifty = templ.cytemp() - midpix; + } //else // calculate new cluster boundaries @@ -913,12 +927,19 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, // Now calculate the mean bias correction and uncertainties float qyfrac = (qfy - qly) / (qfy + qly); - bias = templ.yflcorr(binq, qyfrac) + templ.yavg(binq); + if (goodEdgeAlgo) { + bias = templ.yflcorr(binq, qyfrac) + templ.ygx0(binq); + } else { + bias = templ.yflcorr(binq, qyfrac) + templ.yavg(binq); + } // uncertainty and final correction depend upon charge bin - yrec = (0.125f * binl + BHY - 2.5f + rat * (binh - binl) * 0.125f - (float)shifty + originy) * ysize - bias; - sigmay = templ.yrms(binq); + if (goodEdgeAlgo) { + sigmay = templ.ygsig(binq); + } else { + sigmay = templ.yrms(binq); + } // Do goodness of fit test in y @@ -1202,6 +1223,89 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, return 0; } // PixelTempReco1D +// ************************************************************************************************************************************* +// Overload parameter list for compatibility with older versions +//! Reconstruct the best estimate of the hit position for pixel clusters. +//! \param id - (input) identifier of the template to use +//! \param cotalpha - (input) the cotangent of the alpha track angle (see CMS IN 2004/014) +//! \param cotbeta - (input) the cotangent of the beta track angle (see CMS IN 2004/014) +//! \param locBz - (input) the sign of this quantity is used to determine whether to flip cot(beta)<0 quantities from cot(beta)>0 (FPix only) +//! for Phase 0 FPix IP-related tracks, locBz < 0 for cot(beta) > 0 and locBz > 0 for cot(beta) < 0 +//! for Phase 1 FPix IP-related tracks, see next comment +//! \param locBx - (input) the sign of this quantity is used to determine whether to flip cot(alpha/beta)<0 quantities from cot(alpha/beta)>0 (FPix only) +//! for Phase 1 FPix IP-related tracks, locBx/locBz > 0 for cot(alpha) > 0 and locBx/locBz < 0 for cot(alpha) < 0 +//! for Phase 1 FPix IP-related tracks, locBx > 0 for cot(beta) > 0 and locBx < 0 for cot(beta) < 0//! \param cotbeta - (input) the cotangent of the beta track angle (see CMS IN 2004/014) +//! \param cluster - (input) boost multi_array container of 7x21 array of pixel signals, +//! origin of local coords (0,0) at center of pixel cluster[0][0]. +//! \param ydouble - (input) STL vector of 21 element array to flag a double-pixel +//! \param xdouble - (input) STL vector of 7 element array to flag a double-pixel +//! \param templ - (input) the template used in the reconstruction +//! \param yrec - (output) best estimate of y-coordinate of hit in microns +//! \param sigmay - (output) best estimate of uncertainty on yrec in microns +//! \param proby - (output) probability describing goodness-of-fit for y-reco +//! \param xrec - (output) best estimate of x-coordinate of hit in microns +//! \param sigmax - (output) best estimate of uncertainty on xrec in microns +//! \param probx - (output) probability describing goodness-of-fit for x-reco +//! \param qbin - (output) index (0-4) describing the charge of the cluster +//! [0: 1.5 > zeropix; + int nypix, nxpix; + + return SiPixelTemplateReco::PixelTempReco1D(id, + cotalpha, + cotbeta, + locBz, + locBx, + cluster, + templ, + yrec, + sigmay, + proby, + xrec, + sigmax, + probx, + qbin, + speed, + deadpix, + zeropix, + probQ, + nypix, + nxpix, + goodEdgeAlgo); + +} // PixelTempReco1D + // ************************************************************************************************************************************* // Overload parameter list for compatibility with older versions //! Reconstruct the best estimate of the hit position for pixel clusters. @@ -1258,6 +1362,7 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, const bool deadpix = false; std::vector > zeropix; int nypix, nxpix; + const bool goodEdgeAlgo = false; return SiPixelTemplateReco::PixelTempReco1D(id, cotalpha, @@ -1278,7 +1383,8 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, zeropix, probQ, nypix, - nxpix); + nxpix, + goodEdgeAlgo); } // PixelTempReco1D @@ -1328,6 +1434,7 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, { // Local variables const bool deadpix = false; + const bool goodEdgeAlgo = false; std::vector > zeropix; int nypix, nxpix; float locBx, locBz; @@ -1359,7 +1466,8 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, zeropix, probQ, nypix, - nxpix); + nxpix, + goodEdgeAlgo); } // PixelTempReco1D @@ -1405,6 +1513,7 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, { // Local variables const bool deadpix = false; + const bool goodEdgeAlgo = false; std::vector > zeropix; int nypix, nxpix; float locBx, locBz; @@ -1441,6 +1550,7 @@ int SiPixelTemplateReco::PixelTempReco1D(int id, zeropix, probQ, nypix, - nxpix); + nxpix, + goodEdgeAlgo); } // PixelTempReco1D From be025db71b93a7966a07aa95954d91ce539af948 Mon Sep 17 00:00:00 2001 From: Matej Roguljic Date: Thu, 8 May 2025 09:10:16 +0200 Subject: [PATCH 2/5] fix delta in goodEdgeAlgo and don't run it on x projection --- .../SiPixelTransient/src/SiPixelUtils.cc | 8 +- .../interface/pixelCPEforDevice.h | 74 +++++++++---------- 2 files changed, 41 insertions(+), 41 deletions(-) diff --git a/CondFormats/SiPixelTransient/src/SiPixelUtils.cc b/CondFormats/SiPixelTransient/src/SiPixelUtils.cc index 06da0087cbc49..ba88e8687dcea 100644 --- a/CondFormats/SiPixelTransient/src/SiPixelUtils.cc +++ b/CondFormats/SiPixelTransient/src/SiPixelUtils.cc @@ -52,6 +52,7 @@ namespace siPixelUtils { //--- The `effective' charge width -- particle's path in first and last pixels only float w_eff = std::abs(w_pred) - w_inner; + float delta = w_eff - 0.5 * sum_of_edge * pitch; //--- If the observed charge width is inconsistent with the expectations //--- based on the track, do *not* use w_pred-w_innner. Instead, replace @@ -62,16 +63,15 @@ namespace siPixelUtils { w_eff = pitch * 0.5f * sum_of_edge; // ave. length of edge pixels (first+last) (cm) if (goodEdgeAlgo) { - float delta = w_eff - 0.5 * sum_of_edge * pitch; if (delta / pitch > delta_length_cut) { - // define the centers of the first last last pixel coordinates - float x1 = upper_edge_first_pix - 0.5 * pitchfraction_first * pitch; - float x2 = lower_edge_last_pix + 0.5 * pitchfraction_last * pitch; // observed cluster is much shorter than expected, use one-sided reco if (w_pred > 0.f) { + // x1,x2 are centers of the first last last pixel coordinates + float x1 = upper_edge_first_pix - 0.5 * pitchfraction_first * pitch; float hit_pos = x1 + 0.5 * w_pred; return hit_pos; } else { + float x2 = lower_edge_last_pix + 0.5 * pitchfraction_last * pitch; float hit_pos = x2 + 0.5 * w_pred; return hit_pos; } diff --git a/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h b/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h index 0d8b6d862fcfa..8d98e0b87c697 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h @@ -124,14 +124,14 @@ namespace pixelCPEforDevice { cotbeta = gvy * gvz; } - constexpr inline float check_for_short(uint16_t upper_edge_first_pix, //!< As the name says. - uint16_t lower_edge_last_pix, //!< As the name says. - float lorentz_shift, //!< L-shift at half thickness - float theThickness, //detector thickness - float cot_angle, //!< cot of alpha_ or beta_ - float pitch, //!< thePitchX or thePitchY - bool first_is_big, //!< true if the first is big - bool last_is_big) //!< true if the last is big + constexpr inline float updatePositionForShortCluster(uint16_t upper_edge_first_pix, //!< As the name says. + uint16_t lower_edge_last_pix, //!< As the name says. + float lorentz_shift, //!< L-shift at half thickness + float theThickness, //detector thickness + float cot_angle, //!< cot of alpha_ or beta_ + float pitch, //!< thePitchX or thePitchY + bool first_is_big, //!< true if the first is big + bool last_is_big) //!< true if the last is big { //checks if the observed cluster is much shorter than expected //if yes, returns the position predicted using a one-sided reconstruction @@ -151,23 +151,23 @@ namespace pixelCPEforDevice { pitchfraction_last += 1.0f; float sum_of_edge = pitchfraction_first + pitchfraction_last; - float delta = w_eff - 0.5 * sum_of_edge * pitch; + float delta = w_eff - 0.5f * sum_of_edge * pitch; if (delta / pitch > delta_length_cut) { // define the centers of the first and last pixel coordinates - float x1 = pitch * (upper_edge_first_pix - 0.5 * pitchfraction_first); - float x2 = pitch * (lower_edge_last_pix + 0.5 * pitchfraction_last); // observed cluster is much shorter than expected, use one-sided reco if (w_pred > 0.f) { - float hit_pos = x1 + 0.5 * w_pred; + float x1 = pitch * (upper_edge_first_pix - 0.5f * pitchfraction_first); + float hit_pos = x1 + 0.5f * w_pred; return hit_pos; } else { - float hit_pos = x2 + 0.5 * w_pred; + float x2 = pitch * (lower_edge_last_pix + 0.5f * pitchfraction_last); + float hit_pos = x2 + 0.5f * w_pred; return hit_pos; } } //if(delta/pitch > delta_length_cut) else { - return -999.; + return -999.f; } } @@ -300,29 +300,29 @@ namespace pixelCPEforDevice { auto thickness = detParams.isBarrel ? comParams.theThicknessB : comParams.theThicknessE; if (comParams.goodEdgeAlgo_) { - float xPos_alt = check_for_short(llxl, - urxl, - detParams.chargeWidthX, - thickness, - cotalpha, - detParams.thePitchX, - TrackerTraits::isBigPixX(cp.minRow[ic]), - TrackerTraits::isBigPixX(cp.maxRow[ic])); - - float yPos_alt = check_for_short(llyl, - uryl, - detParams.chargeWidthY, - thickness, - cotbeta, - detParams.thePitchY, - TrackerTraits::isBigPixY(cp.minCol[ic]), - TrackerTraits::isBigPixY(cp.maxCol[ic])); - - if (xPos_alt != -999.) { - // observed cluster is much shorter than expected, use one-sided reco - xPos = xPos_alt; - } - if (yPos_alt != -999.) { + //Even unbroken clusters are short in x, so they will not be "shorter than expected" + //We save time by not checking this + // float xPos_alt = updatePositionForShortCluster(llxl, + // urxl, + // detParams.chargeWidthX, + // thickness, + // cotalpha, + // detParams.thePitchX, + // TrackerTraits::isBigPixX(cp.minRow[ic]), + // TrackerTraits::isBigPixX(cp.maxRow[ic])); + // if (xPos_alt != -999.f) { + // // observed cluster is much shorter than expected, use one-sided reco + // xPos = xPos_alt; + // } + float yPos_alt = updatePositionForShortCluster(llyl, + uryl, + detParams.chargeWidthY, + thickness, + cotbeta, + detParams.thePitchY, + TrackerTraits::isBigPixY(cp.minCol[ic]), + TrackerTraits::isBigPixY(cp.maxCol[ic])); + if (yPos_alt != -999.f) { // observed cluster is much shorter than expected, use one-sided reco yPos = yPos_alt; } From 75a8efc611e67ca566ed9109589736e1a34a7909 Mon Sep 17 00:00:00 2001 From: Matej Roguljic Date: Fri, 9 May 2025 17:49:27 +0200 Subject: [PATCH 3/5] fix goodEdgeAlgo CPE alpaka implementation --- .../interface/pixelCPEforDevice.h | 78 +++++++++---------- 1 file changed, 38 insertions(+), 40 deletions(-) diff --git a/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h b/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h index 8d98e0b87c697..5cc2cd94f1564 100644 --- a/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h +++ b/RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforDevice.h @@ -131,11 +131,13 @@ namespace pixelCPEforDevice { float cot_angle, //!< cot of alpha_ or beta_ float pitch, //!< thePitchX or thePitchY bool first_is_big, //!< true if the first is big - bool last_is_big) //!< true if the last is big - { + bool last_is_big, //!< true if the last is big + int sizeM1) { //checks if the observed cluster is much shorter than expected //if yes, returns the position predicted using a one-sided reconstruction //if not, returns a dummy value -999. + if (0 == sizeM1) // not applicable to size 1 + return -999.f; float delta_length_cut = 2.; //to be fine tuned later float w_inner = pitch * (lower_edge_last_pix - upper_edge_first_pix); float w_pred = theThickness * cot_angle - lorentz_shift; @@ -299,32 +301,28 @@ namespace pixelCPEforDevice { computeAnglesFromDet(detParams, xPos, yPos, cotalpha, cotbeta); auto thickness = detParams.isBarrel ? comParams.theThicknessB : comParams.theThicknessE; + + bool useShortYClusterCPE = false; + float yCorr = 0.f; + //Even unbroken clusters are short in x, so they will not be "shorter than expected" + //We save time by not checking for x if (comParams.goodEdgeAlgo_) { - //Even unbroken clusters are short in x, so they will not be "shorter than expected" - //We save time by not checking this - // float xPos_alt = updatePositionForShortCluster(llxl, - // urxl, - // detParams.chargeWidthX, - // thickness, - // cotalpha, - // detParams.thePitchX, - // TrackerTraits::isBigPixX(cp.minRow[ic]), - // TrackerTraits::isBigPixX(cp.maxRow[ic])); - // if (xPos_alt != -999.f) { - // // observed cluster is much shorter than expected, use one-sided reco - // xPos = xPos_alt; - // } - float yPos_alt = updatePositionForShortCluster(llyl, - uryl, - detParams.chargeWidthY, - thickness, - cotbeta, - detParams.thePitchY, - TrackerTraits::isBigPixY(cp.minCol[ic]), - TrackerTraits::isBigPixY(cp.maxCol[ic])); - if (yPos_alt != -999.f) { - // observed cluster is much shorter than expected, use one-sided reco - yPos = yPos_alt; + float yPosAlt = updatePositionForShortCluster(llyl, + uryl, + detParams.chargeWidthY, + thickness, + cotbeta, + detParams.thePitchY, + TrackerTraits::isBigPixY(cp.minCol[ic]), + TrackerTraits::isBigPixY(cp.maxCol[ic]), + cp.maxCol[ic] - cp.minCol[ic]); + if (yPosAlt != -999.f) { + // observed cluster is much shorter than expected, use one-sided reco + yPos = yPosAlt - + (0.5f * (float)detParams.nCols + TrackerTraits::bigPixYCorrection) * + detParams.thePitchY; // Adjust yPosAlt to match module centered coordinate system of main algorithm + yCorr = detParams.shiftY; + useShortYClusterCPE = true; } } //if(comParams.goodEdgeAlgo_) auto xCorr = correction(cp.maxRow[ic] - cp.minRow[ic], @@ -338,19 +336,19 @@ namespace pixelCPEforDevice { detParams.thePitchX, TrackerTraits::isBigPixX(cp.minRow[ic]), TrackerTraits::isBigPixX(cp.maxRow[ic])); - - auto yCorr = correction(cp.maxCol[ic] - cp.minCol[ic], - cp.q_f_Y[ic], - cp.q_l_Y[ic], - llyl, - uryl, - detParams.chargeWidthY, // lorentz shift in cm - thickness, - cotbeta, - detParams.thePitchY, - TrackerTraits::isBigPixY(cp.minCol[ic]), - TrackerTraits::isBigPixY(cp.maxCol[ic])); - + if (!useShortYClusterCPE) { + yCorr = correction(cp.maxCol[ic] - cp.minCol[ic], + cp.q_f_Y[ic], + cp.q_l_Y[ic], + llyl, + uryl, + detParams.chargeWidthY, // lorentz shift in cm + thickness, + cotbeta, + detParams.thePitchY, + TrackerTraits::isBigPixY(cp.minCol[ic]), + TrackerTraits::isBigPixY(cp.maxCol[ic])); + } //if(!useShortYClusterCPE) cp.xpos[ic] = xPos + xCorr; cp.ypos[ic] = yPos + yCorr; } From 9d233f49a7dd1847a8e9aeb58d77b6a851057a59 Mon Sep 17 00:00:00 2001 From: Matej Roguljic Date: Mon, 12 May 2025 08:09:27 +0200 Subject: [PATCH 4/5] add goodEdgeAlgo wfs and hlt customisation --- Configuration/PyReleaseValidation/README.md | 2 +- .../python/upgradeWorkflowComponents.py | 36 ++++++++++++++++++- .../Configuration/python/customizeHLT.py | 15 ++++++++ 3 files changed, 51 insertions(+), 2 deletions(-) create mode 100644 RecoLocalTracker/Configuration/python/customizeHLT.py diff --git a/Configuration/PyReleaseValidation/README.md b/Configuration/PyReleaseValidation/README.md index 7135a5592a1d2..2e2441e115c3e 100644 --- a/Configuration/PyReleaseValidation/README.md +++ b/Configuration/PyReleaseValidation/README.md @@ -115,4 +115,4 @@ The offsets currently in use are: * 0.112: Activate OuterTracker inefficiency (PS-p: bias rails inefficiency; PS-s and SS: 1% bad strips) * 0.113: Activate OuterTracker inefficiency (PS-p: bias rails inefficiency; PS-s and SS: 5% bad strips) * 0.114: Activate OuterTracker inefficiency (PS-p: bias rails inefficiency; PS-s and SS: 10% bad strips) -* 0.141: Activate emulation of the signal shape of the InnerTracker FE chip (CROC) +* 0.141: Activate emulation of the signal shape of the InnerTracker FE chip (CROC) \ No newline at end of file diff --git a/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py b/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py index e16e9aedc0fd0..4ee11bf7ba5e6 100644 --- a/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py +++ b/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py @@ -640,6 +640,40 @@ def condition(self, fragment, stepList, key, hasHarvest): offset = 0.18, ) +# pixel GoodEdgeAlgo CPE workflows +class UpgradeWorkflow_siPixelGoodEdgeAlgo(UpgradeWorkflow): + def setup_(self, step, stepName, stepDict, k, properties): + if 'HLTOnly'==step: + stepDict[stepName][k] = merge([{'--customise' : 'RecoLocalTracker/Configuration/customizeHLT.customiseFor47966'}, stepDict[step][k]]) + if 'Reco' in step: + stepDict[stepName][k] = merge([{'--procModifiers': 'siPixelGoodEdgeAlgo'}, stepDict[step][k]]) + def condition(self, fragment, stepList, key, hasHarvest): + result = (fragment=="QCD_Pt_1800_2400_14" or fragment=="TTbar_14TeV" ) and any(y in key for y in ['2025']) + return result +upgradeWFs['siPixelGoodEdgeAlgo'] = UpgradeWorkflow_siPixelGoodEdgeAlgo( + steps = [ + 'HLTOnly', + 'DigiTrigger', + 'Reco', + 'RecoFakeHLT', + 'RecoGlobal', + 'RecoGlobalFakeHLT', + 'RecoNano', + 'RecoNanoFakeHLT', + ], + PU = [ + 'HLTOnly', + 'DigiTrigger', + 'Reco', + 'RecoFakeHLT', + 'RecoGlobal', + 'RecoGlobalFakeHLT', + 'RecoNano', + 'RecoNanoFakeHLT', + ], + suffix = '_siPixelGoodEdgeAlgo', + offset = 0.181, +) #Workflow to enable displacedRegionalStep tracking iteration class UpgradeWorkflow_displacedRegional(UpgradeWorkflowTracking): @@ -3705,4 +3739,4 @@ def __init__(self, howMuch, dataset): ('Hydjet_Quenched_MinBias_5362GeV_cfi', UpgradeFragment(U2000by1,'HydjetQMinBias_5362GeV')), ('Hydjet_Quenched_MinBias_5519GeV_cfi', UpgradeFragment(U2000by1,'HydjetQMinBias_5519GeV')), ('SingleMuPt15Eta0_0p4_cfi', UpgradeFragment(Kby(9,100),'SingleMuPt15Eta0p_0p4')), -]) +]) \ No newline at end of file diff --git a/RecoLocalTracker/Configuration/python/customizeHLT.py b/RecoLocalTracker/Configuration/python/customizeHLT.py new file mode 100644 index 0000000000000..df447b10e15f0 --- /dev/null +++ b/RecoLocalTracker/Configuration/python/customizeHLT.py @@ -0,0 +1,15 @@ +import FWCore.ParameterSet.Config as cms + +# helper functions +from HLTrigger.Configuration.common import * + + +# Enable good edge CPE algorithm in pixel local reconstruction +def customiseFor47966(process): + for prod in esproducers_by_type(process, 'PixelCPEGenericESProducer', 'PixelCPEFastParamsESProducerAlpakaPhase1@alpaka'): + if not hasattr(prod, 'GoodEdgeAlgo'): + prod.GoodEdgeAlgo = cms.bool(True) + else: + prod.GoodEdgeAlgo = True + + return process \ No newline at end of file From 61beee252fd9a4492af0012afe5c9d76c9677c05 Mon Sep 17 00:00:00 2001 From: Matej Roguljic Date: Wed, 21 May 2025 13:39:39 +0200 Subject: [PATCH 5/5] Apply suggestions from code review pixel cpe goodEdgeAlgo: simplified generic implementation and resolved wf collision Co-authored-by: Dinko F. --- .../interface/SiPixelTemplate.h | 4 +--- .../SiPixelTransient/src/SiPixelTemplate.cc | 21 ++----------------- Configuration/PyReleaseValidation/README.md | 3 ++- .../python/upgradeWorkflowComponents.py | 4 ++-- .../Configuration/python/customizeHLT.py | 3 ++- 5 files changed, 9 insertions(+), 26 deletions(-) diff --git a/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h b/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h index 6587f2e1dbd5b..7408ec5c8698a 100644 --- a/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h +++ b/CondFormats/SiPixelTransient/interface/SiPixelTemplate.h @@ -286,10 +286,8 @@ class SiPixelTemplate { static void postInit(std::vector& thePixelTemp_); // Interpolate with y Gaussian Parameter interpolation to be used with goodEdge reconstruction algorithm - bool interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx, bool goodEdgeAlgo); - // Interpolate input alpha and beta angles to produce a working template for each individual hit. - bool interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx); + bool interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx, bool goodEdgeAlgo = false); // Interpolate input alpha and beta angles to produce a working template for each individual hit. bool interpolate(int id, float cotalpha, float cotbeta, float locBz); diff --git a/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc b/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc index b6037b8f897f6..c91bfe4a0f25d 100644 --- a/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc +++ b/CondFormats/SiPixelTransient/src/SiPixelTemplate.cc @@ -1876,23 +1876,6 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float l return success_; } // interpolate -// ************************************************************************************************************ -//! Interpolate input alpha and beta angles to produce a working template for each individual hit. -//! \param id - (input) index of the template to use -//! \param cotalpha - (input) the cotangent of the alpha track angle (see CMS IN 2004/014) -//! \param cotbeta - (input) the cotangent of the beta track angle (see CMS IN 2004/014) -//! \param locBz - (input) the sign of this quantity is used to determine whether to flip cot(beta)<0 quantities from cot(beta)>0 (FPix only) -//! for Phase 0 FPix IP-related tracks, locBz < 0 for cot(beta) > 0 and locBz > 0 for cot(beta) < 0 -//! for Phase 1 FPix IP-related tracks, see next comment -//! \param locBx - (input) the sign of this quantity is used to determine whether to flip cot(alpha/beta)<0 quantities from cot(alpha/beta)>0 (FPix only) -//! for Phase 1 FPix IP-related tracks, locBx/locBz > 0 for cot(alpha) > 0 and locBx/locBz < 0 for cot(alpha) < 0 -//! for Phase 1 FPix IP-related tracks, locBx > 0 for cot(beta) > 0 and locBx < 0 for cot(beta) < 0 -// ************************************************************************************************************ -bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx) { - // Interpolate for a new set of track angles, but without y Gaussian Parameter interpolation to be used with goodEdge reconstruction algorithm - return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx, false); -} - // ************************************************************************************************************ //! Interpolate input alpha and beta angles to produce a working template for each individual hit. //! \param id - (input) index of the template to use @@ -1912,7 +1895,7 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta) { if (cotalpha < 0.f) { locBz = -locBx; } - return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx, false); + return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx); } // ************************************************************************************************************ @@ -1927,7 +1910,7 @@ bool SiPixelTemplate::interpolate(int id, float cotalpha, float cotbeta, float l // Local variables float locBx = 1.f; - return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx, false); + return SiPixelTemplate::interpolate(id, cotalpha, cotbeta, locBz, locBx); } // ************************************************************************************************************************************* diff --git a/Configuration/PyReleaseValidation/README.md b/Configuration/PyReleaseValidation/README.md index 2e2441e115c3e..d96d0a7e80f57 100644 --- a/Configuration/PyReleaseValidation/README.md +++ b/Configuration/PyReleaseValidation/README.md @@ -115,4 +115,5 @@ The offsets currently in use are: * 0.112: Activate OuterTracker inefficiency (PS-p: bias rails inefficiency; PS-s and SS: 1% bad strips) * 0.113: Activate OuterTracker inefficiency (PS-p: bias rails inefficiency; PS-s and SS: 5% bad strips) * 0.114: Activate OuterTracker inefficiency (PS-p: bias rails inefficiency; PS-s and SS: 10% bad strips) -* 0.141: Activate emulation of the signal shape of the InnerTracker FE chip (CROC) \ No newline at end of file +* 0.141: Activate emulation of the signal shape of the InnerTracker FE chip (CROC) +* 0.186: Run-3 goodEdgeAlgo CPE \ No newline at end of file diff --git a/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py b/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py index 4ee11bf7ba5e6..14e171ca390a4 100644 --- a/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py +++ b/Configuration/PyReleaseValidation/python/upgradeWorkflowComponents.py @@ -672,7 +672,7 @@ def condition(self, fragment, stepList, key, hasHarvest): 'RecoNanoFakeHLT', ], suffix = '_siPixelGoodEdgeAlgo', - offset = 0.181, + offset = 0.186, ) #Workflow to enable displacedRegionalStep tracking iteration @@ -3739,4 +3739,4 @@ def __init__(self, howMuch, dataset): ('Hydjet_Quenched_MinBias_5362GeV_cfi', UpgradeFragment(U2000by1,'HydjetQMinBias_5362GeV')), ('Hydjet_Quenched_MinBias_5519GeV_cfi', UpgradeFragment(U2000by1,'HydjetQMinBias_5519GeV')), ('SingleMuPt15Eta0_0p4_cfi', UpgradeFragment(Kby(9,100),'SingleMuPt15Eta0p_0p4')), -]) \ No newline at end of file +]) diff --git a/RecoLocalTracker/Configuration/python/customizeHLT.py b/RecoLocalTracker/Configuration/python/customizeHLT.py index df447b10e15f0..238cf6e6d1070 100644 --- a/RecoLocalTracker/Configuration/python/customizeHLT.py +++ b/RecoLocalTracker/Configuration/python/customizeHLT.py @@ -12,4 +12,5 @@ def customiseFor47966(process): else: prod.GoodEdgeAlgo = True - return process \ No newline at end of file + return process + \ No newline at end of file