Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-03-19 06:07:10

0001 //
0002 //  SiPixelTemplate2D.h (v2.65)
0003 //
0004 //  Full 2-D templates for cluster splitting, simulated cluster reweighting, and improved cluster probability
0005 //
0006 // Created by Morris Swartz on 12/01/09.
0007 // V1.01 - fix qavg_ filling
0008 // V1.02 - Add locBz to test if FPix use is out of range
0009 // V1.03 - Fix edge checking on final template to increase template size and to properly truncate cluster
0010 // v2.00 - Major changes to accommodate 2D reconstruction
0011 // v2.10 - Change chi2 and error scaling information to work with partially reconstructed clusters
0012 // v2.20 - Add cluster charge probability information, side loading for template generation
0013 // v2.21 - Double derivative interval [improves fit convergence]
0014 // v2.25 - Resize template store to accommodate FPix Templates
0015 // v2.30 - Fix bug found by P. Shuetze that compromises sqlite file loading
0016 // v2.35 - Add directory path selection to the ascii pushfile method
0017 // v2.50 - Change template storage to dynamically allocated 2D arrays of SiPixelTemplateEntry2D structs
0018 // v2.51 - Ensure that the derivative arrays are correctly zeroed between calls
0019 // v2.52 - Improve cosmetics for increased style points from judges
0020 // v2.60 - Fix FPix multiframe lookup problem [takes +-cotalpha and +-cotbeta]
0021 // v2.61a - Code 2.60 fix correctly
0022 // v2.65 - change double pixel flags to work with new shifted reco code definition
0023 //
0024 
0025 // Build the template storage structure from several pieces
0026 
0027 #ifndef SiPixelTemplate2D_h
0028 #define SiPixelTemplate2D_h 1
0029 
0030 #include <vector>
0031 #include <cassert>
0032 #include "boost/multi_array.hpp"
0033 
0034 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
0035 #include "CondFormats/SiPixelObjects/interface/SiPixel2DTemplateDBObject.h"
0036 #include "FWCore/Utilities/interface/Exception.h"
0037 #include "CondFormats/SiPixelTransient/interface/SiPixelTemplateDefs.h"
0038 #else
0039 #include "SiPixelTemplateDefs.h"
0040 #endif
0041 
0042 struct SiPixelTemplateEntry2D {  //!< Basic template entry corresponding to a single set of track angles
0043   int runnum;                    //!< number of pixelav run used to generate this entry
0044   float cotalpha;                //!< cot(alpha) is proportional to cluster length in x and is basis of interpolation
0045   float cotbeta;                 //!< cot(beta) is proportional to cluster length in y and is basis of interpolation
0046   float costrk[3];               //!< direction cosines of tracks used to generate this entry
0047   float qavg;                    //!< average cluster charge for this set of track angles
0048   float pixmax;                  //!< maximum charge for individual pixels in cluster
0049   float sxymax;                  //!< average pixel signal for use of the error parameterization
0050   int iymin;                     //!< the minimum nonzero pixel yindex in template (saves time during interpolation)
0051   int iymax;                     //!< the maximum nonzero pixel yindex in template (saves time during interpolation)
0052   int jxmin;                     //!< the minimum nonzero pixel xindex in template (saves time during interpolation)
0053   int jxmax;                     //!< the maximum nonzero pixel xindex in template (saves time during interpolation)
0054   float xypar[2][5];             //!< pixel uncertainty parameterization
0055   float lanpar[2][5];            //!< pixel landau distribution parameters
0056   short int xytemp[7][7][T2YSIZE][T2XSIZE];  //!< templates for y-reconstruction (binned over 1 central pixel)
0057   float chi2ppix;                            //!< average chi^2 per pixel
0058   float chi2scale;                           //!< scale factor for the chi2 distribution
0059   float chi2avgone;                          //!< average y chi^2 for 1 pixel clusters
0060   float chi2minone;                          //!< minimum of y chi^2 for 1 pixel clusters
0061   float clsleny;                             //!< cluster y-length in pixels at signal height symax/2
0062   float clslenx;                             //!< cluster x-length in pixels at signal height sxmax/2
0063   float mpvvav;      //!< most probable charge in Vavilov distribution (not actually for larger kappa)
0064   float sigmavav;    //!< "sigma" scale fctor for Vavilov distribution
0065   float kappavav;    //!< kappa parameter for Vavilov distribution
0066   float scalexavg;   //!< average x-error scale factor
0067   float scaleyavg;   //!< average y-error scale factor
0068   float delyavg;     //!< average length difference between template and cluster
0069   float delysig;     //!< rms of length difference between template and cluster
0070   float scalex[4];   //!< x-error scale factor in 4 charge bins
0071   float scaley[4];   //!< y-error scale factor in 4 charge bins
0072   float offsetx[4];  //!< x-offset in 4 charge bins
0073   float offsety[4];  //!< y-offset in 4 charge bins
0074   float spare[3];
0075 };
0076 
0077 struct SiPixelTemplateHeader2D {  //!< template header structure
0078   int ID;                         //!< template ID number
0079   int NTy;                        //!< number of Template y entries
0080   int NTyx;                       //!< number of Template y-slices of x entries
0081   int NTxx;                       //!< number of Template x-entries in each slice
0082   int Dtype;                      //!< detector type (0=BPix, 1=FPix)
0083   float qscale;                   //!< Charge scaling to match cmssw and pixelav
0084   float lorywidth;                //!< estimate of y-lorentz width for optimal resolution
0085   float lorxwidth;                //!< estimate of x-lorentz width for optimal resolution
0086   float lorybias;                 //!< estimate of y-lorentz bias
0087   float lorxbias;                 //!< estimate of x-lorentz bias
0088   float Vbias;                    //!< detector bias potential in Volts
0089   float temperature;              //!< detector temperature in deg K
0090   float fluence;                  //!< radiation fluence in n_eq/cm^2
0091   float s50;                      //!< 1/2 of the multihit dcol threshold in electrons
0092   float ss50;                     //!< 1/2 of the single hit dcol threshold in electrons
0093   char title[80];                 //!< template title
0094   int templ_version;              //!< Version number of the template to ensure code compatibility
0095   float Bfield;                   //!< Bfield in Tesla
0096   float fbin[3];                  //!< The QBin definitions in Q_clus/Q_avg
0097   float xsize;                    //!< pixel size (for future use in upgraded geometry)
0098   float ysize;                    //!< pixel size (for future use in upgraded geometry)
0099   float zsize;                    //!< pixel size (for future use in upgraded geometry)
0100 };
0101 
0102 struct SiPixelTemplateStore2D {  //!< template storage structure
0103   SiPixelTemplateHeader2D head;  //!< Header information
0104 
0105   //!< use 2d entry to store BPix and FPix entries [dynamically allocated]
0106   std::vector<std::vector<SiPixelTemplateEntry2D>> entry;
0107 };
0108 
0109 // ******************************************************************************************
0110 //! \class SiPixelTemplate2D
0111 //!
0112 //!  A template management class.  SiPixelTemplate contains thePixelTemp
0113 //!  (a std::vector  of SiPixelTemplateStore, each of which is a collection of many
0114 //!  SiPixelTemplateEntries).  Each SiPixelTemplateStore corresponds to a given detector
0115 //!  condition, and is valid for a range of runs.  We allow more than one Store since the
0116 //!  may change over time.
0117 //!
0118 //!  This class reads templates from files via pushfile() method.
0119 //!
0120 //!  The main functionality of SiPixelTemplate is xytemp(), which produces a template
0121 //!  on the fly, given a specific track's alpha and beta.  The results are kept in data
0122 //!  members and accessed via inline getters.
0123 //!
0124 //!  The resulting template is then used by PixelTempReco2D() (a global function) which
0125 //!  get the reference for SiPixelTemplate & templ and uses the current template to
0126 //!  reconstruct the SiPixelRecHit.
0127 // ******************************************************************************************
0128 class SiPixelTemplate2D {
0129 public:
0130   SiPixelTemplate2D(const std::vector<SiPixelTemplateStore2D>& thePixelTemp) : thePixelTemp_(thePixelTemp) {
0131     id_current_ = -1;
0132     index_id_ = -1;
0133     cota_current_ = 0.;
0134     cotb_current_ = 0.;
0135   }  //!< Default constructor
0136 
0137   // load the private store with info from the
0138   // file with the index (int) filenum ${dir}template_summary_zp${filenum}.out
0139 #ifdef SI_PIXEL_TEMPLATE_STANDALONE
0140   static bool pushfile(int filenum, std::vector<SiPixelTemplateStore2D>& pixelTemp, std::string dir = "");
0141 
0142   // For calibrations only: load precalculated values -- no interpolation.
0143   void sideload(SiPixelTemplateEntry2D* entry,
0144                 int iDtype,
0145                 float locBx,
0146                 float locBz,
0147                 float lorwdy,
0148                 float lorwdx,
0149                 float q50,
0150                 float fbin[3],
0151                 float xsize,
0152                 float ysize,
0153                 float zsize);
0154 
0155 #else
0156   static bool pushfile(int filenum,
0157                        std::vector<SiPixelTemplateStore2D>& pixelTemp,
0158                        std::string dir = "CalibTracker/SiPixelESProducers/data/");
0159 
0160   // Load from the DB (the default in CMSSW):
0161   static bool pushfile(const SiPixel2DTemplateDBObject& dbobject, std::vector<SiPixelTemplateStore2D>& pixelTemp);
0162 
0163 #endif
0164 
0165   //  Initialize things before interpolating
0166   bool getid(int id);
0167 
0168   bool interpolate(int id, float cotalpha, float cotbeta, float locBz, float locBx);
0169 
0170   // Interpolate input alpha and beta angles to produce a working template for each individual hit.
0171 
0172   // Works with Phase 0+1
0173   bool xytemp(float xhit,
0174               float yhit,
0175               bool ydouble[BYM2],
0176               bool xdouble[BXM2],
0177               float template2d[BXM2][BYM2],
0178               bool dervatives,
0179               float dpdx2d[2][BXM2][BYM2],
0180               float& QTemplate);
0181 
0182   // Overload for backward compatibility
0183 
0184   bool xytemp(float xhit, float yhit, bool ydouble[BYM2], bool xdouble[BXM2], float template2d[BXM2][BYM2]);
0185 
0186   // Overload for backward compatibility with re-weighting code
0187 
0188   bool xytemp(int id,
0189               float cotalpha,
0190               float cotbeta,
0191               float xhit,
0192               float yhit,
0193               std::vector<bool>& ydouble,
0194               std::vector<bool>& xdouble,
0195               float template2d[BXM2][BYM2]);
0196 
0197   void xysigma2(float qpixel, int index, float& xysig2);
0198 
0199   // Get the interpolated Landau distribution parameters
0200 
0201   void landau_par(float lanpar[2][5]);
0202 
0203   float qavg() { return qavg_; }      //!< average cluster charge for this set of track angles
0204   float pixmax() { return pixmax_; }  //!< maximum pixel charge
0205   float qscale() { return qscale_; }  //!< charge scaling factor
0206   float s50() { return s50_; }        //!< 1/2 of the pixel threshold signal in adc units
0207   float sxymax() { return sxymax_; }  //!< max pixel signal for pixel error calculation
0208   float scalex(int i) {
0209     if (checkIllegalIndex("scalex", 3, i)) {
0210       return scalex_[i];
0211     } else {
0212       return 0.f;
0213     }
0214   }  //!< x-error scale factor in 4 charge bins
0215   float scaley(int i) {
0216     if (checkIllegalIndex("scaley", 3, i)) {
0217       return scaley_[i];
0218     } else {
0219       return 0.f;
0220     }
0221   }  //!< y-error scale factor in 4 charge bins
0222   float offsetx(int i) {
0223     if (checkIllegalIndex("offsetx", 3, i)) {
0224       return offsetx_[i];
0225     } else {
0226       return 0.f;
0227     }
0228   }  //!< x-offset in 4 charge bins
0229   float offsety(int i) {
0230     if (checkIllegalIndex("offsety", 3, i)) {
0231       return offsety_[i];
0232     } else {
0233       return 0.f;
0234     }
0235   }  //!< y-offset in 4 charge bins
0236   float fbin(int i) {
0237     if (checkIllegalIndex("fbin", 2, i)) {
0238       return fbin_[i];
0239     } else {
0240       return 0.f;
0241     }
0242   }                                           //!< Return lower bound of Qbin definition
0243   float sizex() { return clslenx_; }          //!< return x size of template cluster
0244   float sizey() { return clsleny_; }          //!< return y size of template cluster
0245   float chi2ppix() { return chi2ppix_; }      //!< average chi^2 per struck pixel
0246   float chi2scale() { return chi2scale_; }    //!< scale factor for chi^2 distribution
0247   float chi2avgone() { return chi2avgone_; }  //!< average y chi^2 for 1 pixel clusters
0248   float chi2minone() { return chi2minone_; }  //!< minimum of y chi^2 for 1 pixel clusters
0249   float mpvvav() { return mpvvav_; }          //!< most probable Q in Vavilov distribution
0250   float sigmavav() { return sigmavav_; }      //!< scale factor in Vavilov distribution
0251   float kappavav() { return kappavav_; }      //!< kappa parameter in Vavilov distribution
0252   float lorydrift() { return lorydrift_; }    //!< signed lorentz y-width (microns)
0253   float lorxdrift() { return lorxdrift_; }    //!< signed lorentz x-width (microns)
0254   float clsleny() { return clsleny_; }        //!< cluster y-size
0255   float clslenx() { return clslenx_; }        //!< cluster x-size
0256   float scaleyavg() { return scaleyavg_; }    //!< y-reco error scaling factor
0257   float scalexavg() { return scalexavg_; }    //!< x-reco error scaling factor
0258   float delyavg() {
0259     return delyavg_;
0260   }  //!< average difference between clsleny_ and cluster length [with threshold effects]
0261   float delysig() { return delysig_; }  //!< rms difference between clsleny_ and cluster length [with threshold effects]
0262   float xsize() { return xsize_; }      //!< pixel x-size (microns)
0263   float ysize() { return ysize_; }      //!< pixel y-size (microns)
0264   float zsize() { return zsize_; }      //!< pixel z-size or thickness (microns)
0265   int storesize() {
0266     return (int)thePixelTemp_.size();
0267   }  //!< return the size of the template store (the number of stored IDs
0268 
0269 private:
0270   bool checkIllegalIndex(const std::string whichMethod, int indMax, int i) {
0271 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
0272     if (i < 0 || i > indMax) {
0273       throw cms::Exception("DataCorrupt")
0274           << "SiPixelTemplate2D::" << whichMethod << " called with illegal index = " << i << std::endl;
0275     }
0276 #else
0277     assert(i >= 0 && i < indMax + 1);
0278 
0279 #endif
0280     return true;
0281   }
0282 
0283   // Keep current template interpolaion parameters
0284 
0285   int id_current_;      //!< current id
0286   int index_id_;        //!< current index
0287   float cota_current_;  //!< current cot alpha
0288   float cotb_current_;  //!< current cot beta
0289   int Nyx_;             //!< number of cot(beta)-entries (columns) in template
0290   int Nxx_;             //!< number of cot(alpha)-entries (rows) in template
0291   int Dtype_;           //!< flags BPix (=0) or FPix (=1)
0292   float cotbeta0_;      //!< minimum cot(beta) covered
0293   float cotbeta1_;      //!< maximum cot(beta) covered
0294   float deltacotb_;     //!< cot(beta) bin size
0295   float cotalpha0_;     //!< minimum cot(alpha) covered
0296   float cotalpha1_;     //!< maximum cot(alpha) covered
0297   float deltacota_;     //!< cot(alpha) bin size
0298   int iy0_;             //!< index of nearest cot(beta) bin
0299   int iy1_;             //!< index of next-nearest cot(beta) bin
0300   float adcotb_;        //!< fractional pixel distance of cot(beta) from iy0_
0301   int jx0_;             //!< index of nearest cot(alpha) bin
0302   int jx1_;             //!< index of next-nearest cot(alpha) bin
0303   float adcota_;        //!< fractional pixel distance of cot(alpha) from jx0_
0304   int imin_;            //!< min y index of templated cluster
0305   int imax_;            //!< max y index of templated cluster
0306   int jmin_;            //!< min x index of templated cluster
0307   int jmax_;            //!< max x index of templated cluster
0308   bool flip_y_;         //!< flip y sign-sensitive quantities
0309   bool flip_x_;         //!< flip x sign-sensitive quantities
0310   bool success_;        //!< true if cotalpha, cotbeta are inside of the acceptance (dynamically loaded)
0311 
0312   // Keep results of last interpolation to return through member functions
0313 
0314   float qavg_;                //!< average cluster charge for this set of track angles
0315   float pixmax_;              //!< maximum pixel charge
0316   float qscale_;              //!< charge scaling factor
0317   float s50_;                 //!< 1/2 of the pixel threshold signal in adc units
0318   float sxymax_;              //!< average pixel signal for y-projection of cluster
0319   float xytemp_[BXM2][BYM2];  //!< template for xy-reconstruction
0320   float xypary0x0_[2][5];     //!< Polynomial error parameterization at ix0,iy0
0321   float xypary1x0_[2][5];     //!< Polynomial error parameterization at ix0,iy1
0322   float xypary0x1_[2][5];     //!< Polynomial error parameterization at ix1,iy0
0323   float lanpar_[2][5];        //!< Interpolated Landau parameters
0324   float chi2ppix_;            //!< average chi^2 per struck pixel
0325   float chi2scale_;           //!< scale factor for chi2 distribution
0326   float chi2avgone_;          //!< average chi^2 for 1 pixel clusters
0327   float chi2minone_;          //!< minimum of chi^2 for 1 pixel clusters
0328   float clsleny_;             //!< projected y-length of cluster
0329   float clslenx_;             //!< projected x-length of cluster
0330   float scalexavg_;           //!< average x-error scale factor
0331   float scaleyavg_;           //!< average y-error scale factor
0332   float delyavg_;             //!< average difference between clsleny_ and cluster length [with threshold effects]
0333   float delysig_;             //!< rms of difference between clsleny_ and cluster length [with threshold effects]
0334   float scalex_[4];           //!< x-error scale factor in charge bins
0335   float scaley_[4];           //!< y-error scale factor in charge bins
0336   float offsetx_[4];          //!< x-offset in charge bins
0337   float offsety_[4];          //!< y-offset in charge bins
0338   float mpvvav_;              //!< most probable Q in Vavilov distribution
0339   float sigmavav_;            //!< scale factor in Vavilov distribution
0340   float kappavav_;            //!< kappa parameter in Vavilov distribution
0341   float lorywidth_;           //!< Lorentz y-width (sign corrected for fpix frame)
0342   float lorxwidth_;           //!< Lorentz x-width
0343   float lorydrift_;           //!< Lorentz y-drift
0344   float lorxdrift_;           //!< Lorentz x-drift
0345   float xsize_;               //!< Pixel x-size
0346   float ysize_;               //!< Pixel y-size
0347   float zsize_;               //!< Pixel z-size (thickness)
0348   float fbin_[3];             //!< The QBin definitions in Q_clus/Q_avg
0349   const SiPixelTemplateEntry2D* entry00_;  // Pointer to presently interpolated point [iy,ix]
0350   const SiPixelTemplateEntry2D* entry10_;  // Pointer to presently interpolated point [iy+1,ix]
0351   const SiPixelTemplateEntry2D* entry01_;  // Pointer to presently interpolated point [iy,ix+1]
0352 
0353   // The actual template store is a std::vector container
0354   const std::vector<SiPixelTemplateStore2D>& thePixelTemp_;
0355 };
0356 
0357 #endif