CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Typedefs | Functions
SiPixelTemplateReco Namespace Reference

Typedefs

typedef boost::multi_array
< float, 2 > 		array_2d
 
typedef boost::multi_array
< float, 3 > 		array_3d
 

Functions

int PixelTempReco2D (int id, float cotalpha, float cotbeta, float locBz, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ)
 
int PixelTempReco2D (int id, float cotalpha, float cotbeta, float locBz, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, float &probQ)
 
int PixelTempReco2D (int id, float cotalpha, float cotbeta, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, float &probQ)
 
int PixelTempReco2D (int id, float cotalpha, float cotbeta, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed)
 
int PixelTempSplit (int id, bool fpix, float cotalpha, float cotbeta, array_2d cluster, std::vector< bool > ydouble, std::vector< bool > xdouble, SiPixelTemplate &templ, float &yrec1, float &yrec2, float &sigmay, float &proby, float &xrec1, float &xrec2, float &sigmax, float &probx, int &qbin, bool deadpix, std::vector< std::pair< int, int > > zeropix)
 
int PixelTempSplit (int id, bool fpix, float cotalpha, float cotbeta, array_2d cluster, std::vector< bool > ydouble, std::vector< bool > xdouble, SiPixelTemplate &templ, float &yrec1, float &yrec2, float &sigmay, float &proby, float &xrec1, float &xrec2, float &sigmax, float &probx, int &qbin)
 

Typedef Documentation

typedef boost::multi_array< float, 2 > SiPixelTemplateReco::array_2d

Definition at line 69 of file SiPixelTemplateReco.h.

typedef boost::multi_array<float, 3> SiPixelTemplateReco::array_3d

Definition at line 45 of file SiPixelTemplateSplit.h.

Function Documentation

int SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
float  locBz,
array_2d &  clust,
std::vector< bool > &  ydouble,
std::vector< bool > &  xdouble,
SiPixelTemplate &  templ,
float &  yrec,
float &  sigmay,
float &  proby,
float &  xrec,
float &  sigmax,
float &  probx,
int &  qbin,
int  speed,
bool  deadpix,
std::vector< std::pair< int, int > > &  zeropix,
float &  probQ 
)

Reconstruct the best estimate of the hit position for pixel clusters.

Parameters
id- (input) identifier of the template to use
cotalpha- (input) the cotangent of the alpha track angle (see CMS IN 2004/014)
cotbeta- (input) the cotangent of the beta track angle (see CMS IN 2004/014)
locBz- (input) the sign of the local B_z field for FPix (usually B_z<0 when cot(beta)>0 and B_z>0 when cot(beta)<0
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]. Set dead pixels to small non-zero values (0.1 e).
ydouble- (input) STL vector of 21 element array to flag a double-pixel
xdouble- (input) STL vector of 7 element array to flag a double-pixel
templ- (input) the template used in the reconstruction
yrec- (output) best estimate of y-coordinate of hit in microns
sigmay- (output) best estimate of uncertainty on yrec in microns
proby- (output) probability describing goodness-of-fit for y-reco
xrec- (output) best estimate of x-coordinate of hit in microns
sigmax- (output) best estimate of uncertainty on xrec in microns
probx- (output) probability describing goodness-of-fit for x-reco
qbin- (output) index (0-4) describing the charge of the cluster qbin = 0 Q/Q_avg > 1.5 [few % of all hits] 1 1.5 > Q/Q_avg > 1.0 [~30% of all hits] 2 1.0 > Q/Q_avg > 0.85 [~30% of all hits] 3 0.85 > Q/Q_avg > min1 [~30% of all hits] 4 min1 > Q/Q_avg > min2 [~0.1% of all hits] 5 min2 > Q/Q_avg [~0.1% of all hits]
speed- (input) switch (-2->5) trading speed vs robustness -2 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4), calculates Q probability w/ VVIObj (better but slower) -1 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4), calculates Q probability w/ TMath::VavilovI (poorer but faster) 0 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4) 1 faster, searches reduced 25 bin range (no big pix) + 33 bins (big pix at ends) at full density 2 faster yet, searches same range as 1 but at 1/2 density 3 fastest, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster) 4 fastest w/ Q prob, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster), calculates Q probability w/ VVIObj (better but slower) 5 fastest w/ Q prob, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster), calculates Q probability w/ TMath::VavilovI (poorer but faster)
deadpix- (input) bool to indicate that there are dead pixels to be included in the analysis
zeropix- (input) vector of index pairs pointing to the dead pixels
probQ- (output) the Vavilov-distribution-based cluster charge probability

Definition at line 120 of file SiPixelTemplateReco.cc.

References BHX, BHY, BXM1, BXM2, BXSIZE, BYM1, BYM2, BYM3, BYSIZE, SiPixelTemplate::chi2xavg(), SiPixelTemplate::chi2xavgone(), SiPixelTemplate::chi2xmin(), SiPixelTemplate::chi2xminone(), SiPixelTemplate::chi2yavg(), SiPixelTemplate::chi2yavgone(), SiPixelTemplate::chi2ymin(), SiPixelTemplate::chi2yminone(), delta, SiPixelTemplate::dxone(), SiPixelTemplate::dxtwo(), SiPixelTemplate::dyone(), SiPixelTemplate::dytwo(), ENDL, edm::hlt::Exception, f, VVIObj::fcn(), Gamma, i, SiPixelTemplate::interpolate(), j, gen::k, LOGDEBUG, LOGERROR, maxpix, SiPixelTemplate::pixmax(), SiPixelTemplate::qavg(), SiPixelTemplate::qmin(), SiPixelTemplate::qscale(), SiPixelTemplate::s50(), ExpressReco_HICollisions_FallBack::sigma, python.multivaluedict::sort(), SiPixelTemplate::sxmax(), SiPixelTemplate::sxone(), SiPixelTemplate::sxtwo(), SiPixelTemplate::symax(), SiPixelTemplate::syone(), SiPixelTemplate::sytwo(), theVerboseLevel, TXSIZE, TYSIZE, SiPixelTemplate::vavilov_pars(), SiPixelTemplate::xavg(), SiPixelTemplate::xflcorr(), SiPixelTemplate::xrms(), SiPixelTemplate::xsigma2(), SiPixelTemplate::xsize(), SiPixelTemplate::xtemp(), SiPixelTemplate::yavg(), SiPixelTemplate::yflcorr(), SiPixelTemplate::yrms(), SiPixelTemplate::ysigma2(), SiPixelTemplate::ysize(), and SiPixelTemplate::ytemp().

Referenced by PixelCPETemplateReco::localPosition(), and PixelTempReco2D().

126 {
127  // Local variables
128  int i, j, k, minbin, binl, binh, binq, midpix, fypix, nypix, lypix, logypx;
129  int fxpix, nxpix, lxpix, logxpx, shifty, shiftx, nyzero[TYSIZE];
130  int nclusx, nclusy;
131  int deltaj, jmin, jmax, fxbin, lxbin, fybin, lybin, djy, djx;
132  int fypix2D, lypix2D, fxpix2D, lxpix2D;
133  float sythr, sxthr, rnorm, delta, sigma, sigavg, pseudopix, qscale, q50;
134  float ss2, ssa, sa2, ssba, saba, sba2, rat, fq, qtotal, qpixel;
135  float originx, originy, qfy, qly, qfx, qlx, bias, maxpix, minmax;
136  double chi2x, meanx, chi2y, meany, chi2ymin, chi2xmin, chi21max;
137  double hchi2, hndof, prvav, mpv, sigmaQ, kappa, xvav, beta2;
138  float ytemp[41][BYSIZE], xtemp[41][BXSIZE], ysum[BYSIZE], xsum[BXSIZE], ysort[BYSIZE], xsort[BXSIZE];
139  float chi2ybin[41], chi2xbin[41], ysig2[BYSIZE], xsig2[BXSIZE];
140  float yw2[BYSIZE], xw2[BXSIZE], ysw[BYSIZE], xsw[BXSIZE];
141  bool yd[BYSIZE], xd[BXSIZE], anyyd, anyxd, calc_probQ, use_VVIObj;
142  float ysize, xsize;
143  const float probmin={1.110223e-16};
144  const float probQmin={1.e-5};
145 
146 // The minimum chi2 for a valid one pixel cluster = pseudopixel contribution only
147 
148  const double mean1pix={0.100}, chi21min={0.160};
149 
150  // First, interpolate the template needed to analyze this cluster
151  // check to see of the track direction is in the physical range of the loaded template
152 
153  if(!templ.interpolate(id, cotalpha, cotbeta, locBz)) {
154  if (theVerboseLevel > 2) {LOGDEBUG("SiPixelTemplateReco") << "input cluster direction cot(alpha) = " << cotalpha << ", cot(beta) = " << cotbeta<< ", local B_z = " << locBz << ", template ID = " << id << ", no reconstruction performed" << ENDL;}
155  return 20;
156  }
157 
158  // Make a local copy of the cluster container so that we can muck with it
159 
160  array_2d cluster = clust;
161 
162  // Check to see if Q probability is selected
163 
164  calc_probQ = false;
165  use_VVIObj = false;
166  if(speed < 0) {
167  calc_probQ = true;
168  if(speed < -1) use_VVIObj = true;
169  speed = 0;
170  }
171 
172  if(speed > 3) {
173  calc_probQ = true;
174  if(speed < 5) use_VVIObj = true;
175  speed = 3;
176  }
177 
178  // Get pixel dimensions from the template (to allow multiple detectors in the future)
179 
180  xsize = templ.xsize();
181  ysize = templ.ysize();
182 
183  // Define size of pseudopixel
184 
185  q50 = templ.s50();
186  pseudopix = 0.2*q50;
187 
188  // Get charge scaling factor
189 
190  qscale = templ.qscale();
191 
192  // Check that the cluster container is (up to) a 7x21 matrix and matches the dimensions of the double pixel flags
193 
194  if(cluster.num_dimensions() != 2) {
195  LOGERROR("SiPixelTemplateReco") << "input cluster container (BOOST Multiarray) has wrong number of dimensions" << ENDL;
196  return 3;
197  }
198  nclusx = (int)cluster.shape()[0];
199  nclusy = (int)cluster.shape()[1];
200  if(nclusx != (int)xdouble.size()) {
201  LOGERROR("SiPixelTemplateReco") << "input cluster container x-size is not equal to double pixel flag container size" << ENDL;
202  return 4;
203  }
204  if(nclusy != (int)ydouble.size()) {
205  LOGERROR("SiPixelTemplateReco") << "input cluster container y-size is not equal to double pixel flag container size" << ENDL;
206  return 5;
207  }
208 
209  // enforce maximum size
210 
211  if(nclusx > TXSIZE) {nclusx = TXSIZE;}
212  if(nclusy > TYSIZE) {nclusy = TYSIZE;}
213 
214  // First, rescale all pixel charges
215 
216  for(j=0; j<nclusx; ++j)
217  for(i=0; i<nclusy; ++i)
218  if(cluster[j][i] > 0) {cluster[j][i] *= qscale;}
219 
220 
221  // Next, sum the total charge and "decapitate" big pixels
222 
223  qtotal = 0.;
224  minmax = templ.pixmax();
225  for(i=0; i<nclusy; ++i) {
226  maxpix = minmax;
227  if(ydouble[i]) {maxpix *=2.;}
228  for(j=0; j<nclusx; ++j) {
229  qtotal += cluster[j][i];
230  if(cluster[j][i] > maxpix) {cluster[j][i] = maxpix;}
231  }
232  }
233 
234  // Do the cluster repair here
235 
236  if(deadpix) {
237  fypix = BYM3; lypix = -1;
238  for(i=0; i<nclusy; ++i) {
239  ysum[i] = 0.; nyzero[i] = 0;
240  // Do preliminary cluster projection in y
241  for(j=0; j<nclusx; ++j) {
242  ysum[i] += cluster[j][i];
243  }
244  if(ysum[i] > 0.) {
245  // identify ends of cluster to determine what the missing charge should be
246  if(i < fypix) {fypix = i;}
247  if(i > lypix) {lypix = i;}
248  }
249  }
250 
251  // Now loop over dead pixel list and "fix" everything
252 
253  //First see if the cluster ends are redefined and that we have only one dead pixel per column
254 
255  std::vector<std::pair<int, int> >::const_iterator zeroIter = zeropix.begin(), zeroEnd = zeropix.end();
256  for ( ; zeroIter != zeroEnd; ++zeroIter ) {
257  i = zeroIter->second;
258  if(i<0 || i>TYSIZE-1) {LOGERROR("SiPixelTemplateReco") << "dead pixel column y-index " << i << ", no reconstruction performed" << ENDL;
259  return 11;}
260 
261  // count the number of dead pixels in each column
262  ++nyzero[i];
263  // allow them to redefine the cluster ends
264  if(i < fypix) {fypix = i;}
265  if(i > lypix) {lypix = i;}
266  }
267 
268  nypix = lypix-fypix+1;
269 
270  // Now adjust the charge in the dead pixels to sum to 0.5*truncation value in the end columns and the truncation value in the interior columns
271 
272  for (zeroIter = zeropix.begin(); zeroIter != zeroEnd; ++zeroIter ) {
273  i = zeroIter->second; j = zeroIter->first;
274  if(j<0 || j>TXSIZE-1) {LOGERROR("SiPixelTemplateReco") << "dead pixel column x-index " << j << ", no reconstruction performed" << ENDL;
275  return 12;}
276  if((i == fypix || i == lypix) && nypix > 1) {maxpix = templ.symax()/2.;} else {maxpix = templ.symax();}
277  if(ydouble[i]) {maxpix *=2.;}
278  if(nyzero[i] > 0 && nyzero[i] < 3) {qpixel = (maxpix - ysum[i])/(float)nyzero[i];} else {qpixel = 1.;}
279  if(qpixel < 1.) {qpixel = 1.;}
280  cluster[j][i] = qpixel;
281  // Adjust the total cluster charge to reflect the charge of the "repaired" cluster
282  qtotal += qpixel;
283  }
284  // End of cluster repair section
285  }
286 
287  // Next, make y-projection of the cluster and copy the double pixel flags into a 25 element container
288 
289  for(i=0; i<BYSIZE; ++i) { ysum[i] = 0.; yd[i] = false;}
290  k=0;
291  anyyd = false;
292  for(i=0; i<nclusy; ++i) {
293  for(j=0; j<nclusx; ++j) {
294  ysum[k] += cluster[j][i];
295  }
296 
297  // If this is a double pixel, put 1/2 of the charge in 2 consective single pixels
298 
299  if(ydouble[i]) {
300  ysum[k] /= 2.;
301  ysum[k+1] = ysum[k];
302  yd[k] = true;
303  yd[k+1] = false;
304  k=k+2;
305  anyyd = true;
306  } else {
307  yd[k] = false;
308  ++k;
309  }
310  if(k > BYM1) {break;}
311  }
312 
313  // Next, make x-projection of the cluster and copy the double pixel flags into an 11 element container
314 
315  for(i=0; i<BXSIZE; ++i) { xsum[i] = 0.; xd[i] = false;}
316  k=0;
317  anyxd = false;
318  for(j=0; j<nclusx; ++j) {
319  for(i=0; i<nclusy; ++i) {
320  xsum[k] += cluster[j][i];
321  }
322 
323  // If this is a double pixel, put 1/2 of the charge in 2 consective single pixels
324 
325  if(xdouble[j]) {
326  xsum[k] /= 2.;
327  xsum[k+1] = xsum[k];
328  xd[k]=true;
329  xd[k+1]=false;
330  k=k+2;
331  anyxd = true;
332  } else {
333  xd[k]=false;
334  ++k;
335  }
336  if(k > BXM1) {break;}
337  }
338 
339  // next, identify the y-cluster ends, count total pixels, nypix, and logical pixels, logypx
340 
341  fypix=-1;
342  nypix=0;
343  lypix=0;
344  logypx=0;
345  for(i=0; i<BYSIZE; ++i) {
346  if(ysum[i] > 0.) {
347  if(fypix == -1) {fypix = i;}
348  if(!yd[i]) {
349  ysort[logypx] = ysum[i];
350  ++logypx;
351  }
352  ++nypix;
353  lypix = i;
354  }
355  }
356 
357  // dlengthy = (float)nypix - templ.clsleny();
358 
359  // Make sure cluster is continuous
360 
361  if((lypix-fypix+1) != nypix || nypix == 0) {
362  LOGDEBUG("SiPixelTemplateReco") << "y-length of pixel cluster doesn't agree with number of pixels above threshold" << ENDL;
363  if (theVerboseLevel > 2) {
364  LOGDEBUG("SiPixelTemplateReco") << "ysum[] = ";
365  for(i=0; i<BYSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << ysum[i] << ", ";}
366  LOGDEBUG("SiPixelTemplateReco") << ysum[BYSIZE-1] << ENDL;
367  }
368 
369  return 1;
370  }
371 
372  // If cluster is longer than max template size, technique fails
373 
374  if(nypix > TYSIZE) {
375  LOGDEBUG("SiPixelTemplateReco") << "y-length of pixel cluster is larger than maximum template size" << ENDL;
376  if (theVerboseLevel > 2) {
377  LOGDEBUG("SiPixelTemplateReco") << "ysum[] = ";
378  for(i=0; i<BYSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << ysum[i] << ", ";}
379  LOGDEBUG("SiPixelTemplateReco") << ysum[BYSIZE-1] << ENDL;
380  }
381 
382  return 6;
383  }
384 
385  // Remember these numbers for later
386 
387  fypix2D = fypix;
388  lypix2D = lypix;
389 
390  // next, center the cluster on pixel 12 if necessary
391 
392  midpix = (fypix+lypix)/2;
393  shifty = BHY - midpix;
394  if(shifty > 0) {
395  for(i=lypix; i>=fypix; --i) {
396  ysum[i+shifty] = ysum[i];
397  ysum[i] = 0.;
398  yd[i+shifty] = yd[i];
399  yd[i] = false;
400  }
401  } else if (shifty < 0) {
402  for(i=fypix; i<=lypix; ++i) {
403  ysum[i+shifty] = ysum[i];
404  ysum[i] = 0.;
405  yd[i+shifty] = yd[i];
406  yd[i] = false;
407  }
408  }
409  lypix +=shifty;
410  fypix +=shifty;
411 
412  // If the cluster boundaries are OK, add pesudopixels, otherwise quit
413 
414  if(fypix > 1 && fypix < BYM2) {
415  ysum[fypix-1] = pseudopix;
416  ysum[fypix-2] = pseudopix;
417  } else {return 8;}
418  if(lypix > 1 && lypix < BYM2) {
419  ysum[lypix+1] = pseudopix;
420  ysum[lypix+2] = pseudopix;
421  } else {return 8;}
422 
423  // finally, determine if pixel[0] is a double pixel and make an origin correction if it is
424 
425  if(ydouble[0]) {
426  originy = -0.5;
427  } else {
428  originy = 0.;
429  }
430 
431  // next, identify the x-cluster ends, count total pixels, nxpix, and logical pixels, logxpx
432 
433  fxpix=-1;
434  nxpix=0;
435  lxpix=0;
436  logxpx=0;
437  for(i=0; i<BXSIZE; ++i) {
438  if(xsum[i] > 0.) {
439  if(fxpix == -1) {fxpix = i;}
440  if(!xd[i]) {
441  xsort[logxpx] = xsum[i];
442  ++logxpx;
443  }
444  ++nxpix;
445  lxpix = i;
446  }
447  }
448 
449  // dlengthx = (float)nxpix - templ.clslenx();
450 
451  // Make sure cluster is continuous
452 
453  if((lxpix-fxpix+1) != nxpix) {
454 
455  LOGDEBUG("SiPixelTemplateReco") << "x-length of pixel cluster doesn't agree with number of pixels above threshold" << ENDL;
456  if (theVerboseLevel > 2) {
457  LOGDEBUG("SiPixelTemplateReco") << "xsum[] = ";
458  for(i=0; i<BXSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << xsum[i] << ", ";}
459  LOGDEBUG("SiPixelTemplateReco") << ysum[BXSIZE-1] << ENDL;
460  }
461 
462  return 2;
463  }
464 
465  // If cluster is longer than max template size, technique fails
466 
467  if(nxpix > TXSIZE) {
468 
469  LOGDEBUG("SiPixelTemplateReco") << "x-length of pixel cluster is larger than maximum template size" << ENDL;
470  if (theVerboseLevel > 2) {
471  LOGDEBUG("SiPixelTemplateReco") << "xsum[] = ";
472  for(i=0; i<BXSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << xsum[i] << ", ";}
473  LOGDEBUG("SiPixelTemplateReco") << ysum[BXSIZE-1] << ENDL;
474  }
475 
476  return 7;
477  }
478 
479  // Remember these numbers for later
480 
481  fxpix2D = fxpix;
482  lxpix2D = lxpix;
483 
484  // next, center the cluster on pixel 5 if necessary
485 
486  midpix = (fxpix+lxpix)/2;
487  shiftx = BHX - midpix;
488  if(shiftx > 0) {
489  for(i=lxpix; i>=fxpix; --i) {
490  xsum[i+shiftx] = xsum[i];
491  xsum[i] = 0.;
492  xd[i+shiftx] = xd[i];
493  xd[i] = false;
494  }
495  } else if (shiftx < 0) {
496  for(i=fxpix; i<=lxpix; ++i) {
497  xsum[i+shiftx] = xsum[i];
498  xsum[i] = 0.;
499  xd[i+shiftx] = xd[i];
500  xd[i] = false;
501  }
502  }
503  lxpix +=shiftx;
504  fxpix +=shiftx;
505 
506  // If the cluster boundaries are OK, add pesudopixels, otherwise quit
507 
508  if(fxpix > 1 && fxpix < BXM2) {
509  xsum[fxpix-1] = pseudopix;
510  xsum[fxpix-2] = pseudopix;
511  } else {return 9;}
512  if(lxpix > 1 && lxpix < BXM2) {
513  xsum[lxpix+1] = pseudopix;
514  xsum[lxpix+2] = pseudopix;
515  } else {return 9;}
516 
517  // finally, determine if pixel[0] is a double pixel and make an origin correction if it is
518 
519  if(xdouble[0]) {
520  originx = -0.5;
521  } else {
522  originx = 0.;
523  }
524 
525  // uncertainty and final corrections depend upon total charge bin
526 
527  fq = qtotal/templ.qavg();
528  if(fq > 1.5) {
529  binq=0;
530  } else {
531  if(fq > 1.0) {
532  binq=1;
533  } else {
534  if(fq > 0.85) {
535  binq=2;
536  } else {
537  binq=3;
538  }
539  }
540  }
541 
542  // Return the charge bin via the parameter list unless the charge is too small (then flag it)
543 
544  qbin = binq;
545  if(!deadpix && qtotal < 0.95*templ.qmin()) {qbin = 5;} else {
546  if(!deadpix && qtotal < 0.95*templ.qmin(1)) {qbin = 4;}
547  }
548  if (theVerboseLevel > 9) {
549  LOGDEBUG("SiPixelTemplateReco") <<
550  "ID = " << id <<
551  " cot(alpha) = " << cotalpha << " cot(beta) = " << cotbeta <<
552  " nclusx = " << nclusx << " nclusy = " << nclusy << ENDL;
553  }
554 
555 
556  // Next, copy the y- and x-templates to local arrays
557 
558  // First, decide on chi^2 min search parameters
559 
560 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
561  if(speed < 0 || speed > 3) {
562  throw cms::Exception("DataCorrupt") << "SiPixelTemplateReco::PixelTempReco2D called with illegal speed = " << speed << std::endl;
563  }
564 #else
565  assert(speed >= 0 && speed < 4);
566 #endif
567  fybin = 2; lybin = 38; fxbin = 2; lxbin = 38; djy = 1; djx = 1;
568  if(speed > 0) {
569  fybin = 8; lybin = 32;
570  if(yd[fypix]) {fybin = 4; lybin = 36;}
571  if(lypix > fypix) {
572  if(yd[lypix-1]) {fybin = 4; lybin = 36;}
573  }
574  fxbin = 8; lxbin = 32;
575  if(xd[fxpix]) {fxbin = 4; lxbin = 36;}
576  if(lxpix > fxpix) {
577  if(xd[lxpix-1]) {fxbin = 4; lxbin = 36;}
578  }
579  }
580 
581  if(speed > 1) {
582  djy = 2; djx = 2;
583  if(speed > 2) {
584  if(!anyyd) {djy = 4;}
585  if(!anyxd) {djx = 4;}
586  }
587  }
588 
589  if (theVerboseLevel > 9) {
590  LOGDEBUG("SiPixelTemplateReco") <<
591  "fypix " << fypix << " lypix = " << lypix <<
592  " fybin = " << fybin << " lybin = " << lybin <<
593  " djy = " << djy << " logypx = " << logypx << ENDL;
594  LOGDEBUG("SiPixelTemplateReco") <<
595  "fxpix " << fxpix << " lxpix = " << lxpix <<
596  " fxbin = " << fxbin << " lxbin = " << lxbin <<
597  " djx = " << djx << " logxpx = " << logxpx << ENDL;
598  }
599 
600  // Now do the copies
601 
602  templ.ytemp(fybin, lybin, ytemp);
603 
604  templ.xtemp(fxbin, lxbin, xtemp);
605 
606  // Do the y-reconstruction first
607 
608  // Apply the first-pass template algorithm to all clusters
609 
610  // Modify the template if double pixels are present
611 
612  if(nypix > logypx) {
613  i=fypix;
614  while(i < lypix) {
615  if(yd[i] && !yd[i+1]) {
616  for(j=fybin; j<=lybin; ++j) {
617 
618  // Sum the adjacent cells and put the average signal in both
619 
620  sigavg = (ytemp[j][i] + ytemp[j][i+1])/2.;
621  ytemp[j][i] = sigavg;
622  ytemp[j][i+1] = sigavg;
623  }
624  i += 2;
625  } else {
626  ++i;
627  }
628  }
629  }
630 
631  // Define the maximum signal to allow before de-weighting a pixel
632 
633  sythr = 1.1*(templ.symax());
634 
635  // Make sure that there will be at least two pixels that are not de-weighted
636 
637  std::sort(&ysort[0], &ysort[logypx]);
638  if(logypx == 1) {sythr = 1.01f*ysort[0];} else {
639  if (ysort[1] > sythr) { sythr = 1.01f*ysort[1]; }
640  }
641 
642  // Evaluate pixel-by-pixel uncertainties (weights) for the templ analysis
643 
644  // for(i=0; i<BYSIZE; ++i) { ysig2[i] = 0.;}
645  templ.ysigma2(fypix, lypix, sythr, ysum, ysig2);
646 
647  // Find the template bin that minimizes the Chi^2
648 
649  chi2ymin = 1.e15;
650  for(i=fybin; i<=lybin; ++i) {
651  chi2ybin[i] = -1.e15f;
652  }
653  ss2 = 0.f;
654  for(i=fypix-2; i<=lypix+2; ++i) {
655  yw2[i] = 1.f/ysig2[i];
656  ysw[i] = ysum[i]*yw2[i];
657  ss2 += ysum[i]*ysw[i];
658  }
659  minbin = -1;
660  deltaj = djy;
661  jmin = fybin;
662  jmax = lybin;
663  while(deltaj > 0) {
664  for(j=jmin; j<=jmax; j+=deltaj) {
665  if(chi2ybin[j] < -100.f) {
666  ssa = 0.f;
667  sa2 = 0.f;
668  for(i=fypix-2; i<=lypix+2; ++i) {
669  ssa += ysw[i]*ytemp[j][i];
670  sa2 += ytemp[j][i]*ytemp[j][i]*yw2[i];
671  }
672  rat=ssa/ss2;
673  if(rat <= 0.f) {LOGERROR("SiPixelTemplateReco") << "illegal chi2ymin normalization (1) = " << rat << ENDL; rat = 1.;}
674  chi2ybin[j]=ss2-2.f*ssa/rat+sa2/(rat*rat);
675  }
676  if(chi2ybin[j] < chi2ymin) {
677  chi2ymin = chi2ybin[j];
678  minbin = j;
679  }
680  }
681  deltaj /= 2;
682  if(minbin > fybin) {jmin = minbin - deltaj;} else {jmin = fybin;}
683  if(minbin < lybin) {jmax = minbin + deltaj;} else {jmax = lybin;}
684  }
685 
686  if (theVerboseLevel > 9) {
687  LOGDEBUG("SiPixelTemplateReco") <<
688  "minbin " << minbin << " chi2ymin = " << chi2ymin << ENDL;
689  }
690 
691  // Do not apply final template pass to 1-pixel clusters (use calibrated offset)
692 
693  if(logypx == 1) {
694 
695  if(nypix ==1) {
696  delta = templ.dyone();
697  sigma = templ.syone();
698  } else {
699  delta = templ.dytwo();
700  sigma = templ.sytwo();
701  }
702 
703  yrec = 0.5f*(fypix+lypix-2*shifty+2.f*originy)*ysize-delta;
704  if(sigma <= 0.) {
705  sigmay = 43.3f;
706  } else {
707  sigmay = sigma;
708  }
709 
710  // Do probability calculation for one-pixel clusters
711 
712  chi21max = fmax(chi21min, (double)templ.chi2yminone());
713  chi2ymin -=chi21max;
714  if(chi2ymin < 0.) {chi2ymin = 0.;}
715  // proby = gsl_cdf_chisq_Q(chi2ymin, mean1pix);
716  meany = fmax(mean1pix, (double)templ.chi2yavgone());
717  hchi2 = chi2ymin/2.; hndof = meany/2.;
718  proby = 1. - TMath::Gamma(hndof, hchi2);
719 
720  } else {
721 
722  // For cluster > 1 pix, make the second, interpolating pass with the templates
723 
724  binl = minbin - 1;
725  binh = binl + 2;
726  if(binl < fybin) { binl = fybin;}
727  if(binh > lybin) { binh = lybin;}
728  ssa = 0.;
729  sa2 = 0.;
730  ssba = 0.;
731  saba = 0.;
732  sba2 = 0.;
733  for(i=fypix-2; i<=lypix+2; ++i) {
734  ssa += ysw[i]*ytemp[binl][i];
735  sa2 += ytemp[binl][i]*ytemp[binl][i]*yw2[i];
736  ssba += ysw[i]*(ytemp[binh][i] - ytemp[binl][i]);
737  saba += ytemp[binl][i]*(ytemp[binh][i] - ytemp[binl][i])*yw2[i];
738  sba2 += (ytemp[binh][i] - ytemp[binl][i])*(ytemp[binh][i] - ytemp[binl][i])*yw2[i];
739  }
740 
741  // rat is the fraction of the "distance" from template a to template b
742 
743  rat=(ssba*ssa-ss2*saba)/(ss2*sba2-ssba*ssba);
744  if(rat < 0.) {rat=0.;}
745  if(rat > 1.) {rat=1.0;}
746  rnorm = (ssa+rat*ssba)/ss2;
747 
748  // Calculate the charges in the first and last pixels
749 
750  qfy = ysum[fypix];
751  if(yd[fypix]) {qfy+=ysum[fypix+1];}
752  if(logypx > 1) {
753  qly=ysum[lypix];
754  if(yd[lypix-1]) {qly+=ysum[lypix-1];}
755  } else {
756  qly = qfy;
757  }
758 
759  // Now calculate the mean bias correction and uncertainties
760 
761  float qyfrac = (qfy-qly)/(qfy+qly);
762  bias = templ.yflcorr(binq,qyfrac)+templ.yavg(binq);
763 
764  // uncertainty and final correction depend upon charge bin
765 
766  yrec = (0.125f*binl+BHY-2.5f+rat*(binh-binl)*0.125f-(float)shifty+originy)*ysize - bias;
767  sigmay = templ.yrms(binq);
768 
769  // Do goodness of fit test in y
770 
771  if(rnorm <= 0.) {LOGERROR("SiPixelTemplateReco") << "illegal chi2y normalization (2) = " << rnorm << ENDL; rnorm = 1.;}
772  chi2y=ss2-2./rnorm*ssa-2./rnorm*rat*ssba+(sa2+2.*rat*saba+rat*rat*sba2)/(rnorm*rnorm)-templ.chi2ymin(binq);
773  if(chi2y < 0.0) {chi2y = 0.0;}
774  meany = templ.chi2yavg(binq);
775  if(meany < 0.01) {meany = 0.01;}
776  // gsl function that calculates the chi^2 tail prob for non-integral dof
777  // proby = gsl_cdf_chisq_Q(chi2y, meany);
778  // proby = ROOT::Math::chisquared_cdf_c(chi2y, meany);
779  hchi2 = chi2y/2.; hndof = meany/2.;
780  proby = 1. - TMath::Gamma(hndof, hchi2);
781  }
782 
783  // Do the x-reconstruction next
784 
785  // Apply the first-pass template algorithm to all clusters
786 
787  // Modify the template if double pixels are present
788 
789  if(nxpix > logxpx) {
790  i=fxpix;
791  while(i < lxpix) {
792  if(xd[i] && !xd[i+1]) {
793  for(j=fxbin; j<=lxbin; ++j) {
794 
795  // Sum the adjacent cells and put the average signal in both
796 
797  sigavg = (xtemp[j][i] + xtemp[j][i+1])/2.;
798  xtemp[j][i] = sigavg;
799  xtemp[j][i+1] = sigavg;
800  }
801  i += 2;
802  } else {
803  ++i;
804  }
805  }
806  }
807 
808  // Define the maximum signal to allow before de-weighting a pixel
809 
810  sxthr = 1.1f*templ.sxmax();
811 
812  // Make sure that there will be at least two pixels that are not de-weighted
813  std::sort(&xsort[0], &xsort[logxpx]);
814  if(logxpx == 1) {sxthr = 1.01f*xsort[0];} else {
815  if (xsort[1] > sxthr) { sxthr = 1.01f*xsort[1]; }
816  }
817 
818  // Evaluate pixel-by-pixel uncertainties (weights) for the templ analysis
819 
820  // for(i=0; i<BXSIZE; ++i) { xsig2[i] = 0.; }
821  templ.xsigma2(fxpix, lxpix, sxthr, xsum, xsig2);
822 
823  // Find the template bin that minimizes the Chi^2
824 
825  chi2xmin = 1.e15;
826  for(i=fxbin; i<=lxbin; ++i) { chi2xbin[i] = -1.e15f;}
827  ss2 = 0.f;
828  for(i=fxpix-2; i<=lxpix+2; ++i) {
829  xw2[i] = 1.f/xsig2[i];
830  xsw[i] = xsum[i]*xw2[i];
831  ss2 += xsum[i]*xsw[i];
832  }
833 
834  minbin = -1;
835  deltaj = djx;
836  jmin = fxbin;
837  jmax = lxbin;
838  while(deltaj > 0) {
839  for(j=jmin; j<=jmax; j+=deltaj) {
840  if(chi2xbin[j] < -100.f) {
841  ssa = 0.f;
842  sa2 = 0.f;
843  for(i=fxpix-2; i<=lxpix+2; ++i) {
844  ssa += xsw[i]*xtemp[j][i];
845  sa2 += xtemp[j][i]*xtemp[j][i]*xw2[i];
846  }
847  rat=ssa/ss2;
848  if(rat <= 0.f) {LOGERROR("SiPixelTemplateReco") << "illegal chi2xmin normalization (1) = " << rat << ENDL; rat = 1.;}
849  chi2xbin[j]=ss2-2.f*ssa/rat+sa2/(rat*rat);
850  }
851  if(chi2xbin[j] < chi2xmin) {
852  chi2xmin = chi2xbin[j];
853  minbin = j;
854  }
855  }
856  deltaj /= 2;
857  if(minbin > fxbin) {jmin = minbin - deltaj;} else {jmin = fxbin;}
858  if(minbin < lxbin) {jmax = minbin + deltaj;} else {jmax = lxbin;}
859  }
860 
861  if (theVerboseLevel > 9) {
862  LOGDEBUG("SiPixelTemplateReco") <<
863  "minbin " << minbin << " chi2xmin = " << chi2xmin << ENDL;
864  }
865 
866  // Do not apply final template pass to 1-pixel clusters (use calibrated offset)
867 
868  if(logxpx == 1) {
869 
870  if(nxpix ==1) {
871  delta = templ.dxone();
872  sigma = templ.sxone();
873  } else {
874  delta = templ.dxtwo();
875  sigma = templ.sxtwo();
876  }
877  xrec = 0.5*(fxpix+lxpix-2*shiftx+2.*originx)*xsize-delta;
878  if(sigma <= 0.) {
879  sigmax = 28.9;
880  } else {
881  sigmax = sigma;
882  }
883 
884  // Do probability calculation for one-pixel clusters
885 
886  chi21max = fmax(chi21min, (double)templ.chi2xminone());
887  chi2xmin -=chi21max;
888  if(chi2xmin < 0.) {chi2xmin = 0.;}
889  meanx = fmax(mean1pix, (double)templ.chi2xavgone());
890  hchi2 = chi2xmin/2.; hndof = meanx/2.;
891  probx = 1. - TMath::Gamma(hndof, hchi2);
892 
893  } else {
894 
895  // Now make the second, interpolating pass with the templates
896 
897  binl = minbin - 1;
898  binh = binl + 2;
899  if(binl < fxbin) { binl = fxbin;}
900  if(binh > lxbin) { binh = lxbin;}
901  ssa = 0.;
902  sa2 = 0.;
903  ssba = 0.;
904  saba = 0.;
905  sba2 = 0.;
906  for(i=fxpix-2; i<=lxpix+2; ++i) {
907  ssa += xsw[i]*xtemp[binl][i];
908  sa2 += xtemp[binl][i]*xtemp[binl][i]*xw2[i];
909  ssba += xsw[i]*(xtemp[binh][i] - xtemp[binl][i]);
910  saba += xtemp[binl][i]*(xtemp[binh][i] - xtemp[binl][i])*xw2[i];
911  sba2 += (xtemp[binh][i] - xtemp[binl][i])*(xtemp[binh][i] - xtemp[binl][i])*xw2[i];
912  }
913 
914  // rat is the fraction of the "distance" from template a to template b
915 
916  rat=(ssba*ssa-ss2*saba)/(ss2*sba2-ssba*ssba);
917  if(rat < 0.f) {rat=0.f;}
918  if(rat > 1.f) {rat=1.0f;}
919  rnorm = (ssa+rat*ssba)/ss2;
920 
921  // Calculate the charges in the first and last pixels
922 
923  qfx = xsum[fxpix];
924  if(xd[fxpix]) {qfx+=xsum[fxpix+1];}
925  if(logxpx > 1) {
926  qlx=xsum[lxpix];
927  if(xd[lxpix-1]) {qlx+=xsum[lxpix-1];}
928  } else {
929  qlx = qfx;
930  }
931 
932  // Now calculate the mean bias correction and uncertainties
933 
934  float qxfrac = (qfx-qlx)/(qfx+qlx);
935  bias = templ.xflcorr(binq,qxfrac)+templ.xavg(binq);
936 
937  // uncertainty and final correction depend upon charge bin
938 
939  xrec = (0.125f*binl+BHX-2.5f+rat*(binh-binl)*0.125f-(float)shiftx+originx)*xsize - bias;
940  sigmax = templ.xrms(binq);
941 
942  // Do goodness of fit test in x
943 
944  if(rnorm <= 0.) {LOGERROR("SiPixelTemplateReco") << "illegal chi2x normalization (2) = " << rnorm << ENDL; rnorm = 1.;}
945  chi2x=ss2-2./rnorm*ssa-2./rnorm*rat*ssba+(sa2+2.*rat*saba+rat*rat*sba2)/(rnorm*rnorm)-templ.chi2xmin(binq);
946  if(chi2x < 0.0) {chi2x = 0.0;}
947  meanx = templ.chi2xavg(binq);
948  if(meanx < 0.01) {meanx = 0.01;}
949  // gsl function that calculates the chi^2 tail prob for non-integral dof
950  // probx = gsl_cdf_chisq_Q(chi2x, meanx);
951  // probx = ROOT::Math::chisquared_cdf_c(chi2x, meanx, trx0);
952  hchi2 = chi2x/2.; hndof = meanx/2.;
953  probx = 1. - TMath::Gamma(hndof, hchi2);
954  }
955 
956  // Don't return exact zeros for the probability
957 
958  if(proby < probmin) {proby = probmin;}
959  if(probx < probmin) {probx = probmin;}
960 
961  // Decide whether to generate a cluster charge probability
962 
963  if(calc_probQ) {
964 
965  // Calculate the Vavilov probability that the cluster charge is OK
966 
967  templ.vavilov_pars(mpv, sigmaQ, kappa);
968 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
969  if((sigmaQ <=0.) || (mpv <= 0.) || (kappa < 0.01) || (kappa > 9.9)) {
970  throw cms::Exception("DataCorrupt") << "SiPixelTemplateReco::Vavilov parameters mpv/sigmaQ/kappa = " << mpv << "/" << sigmaQ << "/" << kappa << std::endl;
971  }
972 #else
973  assert((sigmaQ > 0.) && (mpv > 0.) && (kappa > 0.01) && (kappa < 10.));
974 #endif
975  xvav = ((double)qtotal-mpv)/sigmaQ;
976  beta2 = 1.;
977  if(use_VVIObj) {
978  // VVIObj is a private port of CERNLIB VVIDIS
979  VVIObj vvidist(kappa, beta2, 1);
980  prvav = vvidist.fcn(xvav);
981 
982  // std::cout << "vav: " << kappa << " " << xvav << " " << prvav
983  // << " " << TMath::VavilovI(xvav, kappa, beta2) << std::endl;
984 
985  } else {
986  // Use faster but less accurate TMath Vavilov distribution function
987  prvav = TMath::VavilovI(xvav, kappa, beta2);
988  }
989  // Change to upper tail probability
990  // if(prvav > 0.5) prvav = 1. - prvav;
991  // probQ = (float)(2.*prvav);
992  probQ = 1. - prvav;
993  if(probQ < probQmin) {probQ = probQmin;}
994  } else {
995  probQ = -1;
996  }
997 
998  return 0;
999 } // PixelTempReco2D
ExpressReco_HICollisions_FallBack.sigma
tuple sigma
Definition: ExpressReco_HICollisions_FallBack.py:4426
delta
dbl * delta
Definition: mlp_gen.cc:36
ExpressReco_HICollisions_FallBack.speed
tuple speed
Definition: ExpressReco_HICollisions_FallBack.py:23342
i
int i
Definition: DBlmapReader.cc:9
BXSIZE
#define BXSIZE
Definition: SiPixelTemplateDefs.h:26
SiPixelTemplate::chi2xminone
float chi2xminone()
//!&lt; minimum of x chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:441
SiPixelTemplate::symax
float symax()
average pixel signal for y-projection of cluster
Definition: SiPixelTemplate.h:275
SiPixelTemplate::yavg
float yavg(int i)
average y-bias of reconstruction binned in 4 charge bins
Definition: SiPixelTemplate.h:298
BYSIZE
#define BYSIZE
Definition: SiPixelTemplateDefs.h:18
TXSIZE
#define TXSIZE
Definition: SiPixelTemplateDefs.h:23
SiPixelTemplate::chi2xmin
float chi2xmin(int i)
minimum y chi^2 in 4 charge bins
Definition: SiPixelTemplate.h:375
BXM1
#define BXM1
Definition: SiPixelTemplateDefs.h:28
SiPixelTemplate::chi2ymin
float chi2ymin(int i)
minimum y chi^2 in 4 charge bins
Definition: SiPixelTemplate.h:361
SiPixelTemplate::qmin
float qmin()
minimum cluster charge for valid hit (keeps 99.9% of simulated hits)
Definition: SiPixelTemplate.h:285
SiPixelTemplate::xflcorr
float xflcorr(int binq, float qflx)
Definition: SiPixelTemplate.cc:1573
SiPixelTemplate::sytwo
float sytwo()
rms for one double-pixel y-clusters
Definition: SiPixelTemplate.h:279
SiPixelTemplate::chi2yminone
float chi2yminone()
//!&lt; minimum of y chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:439
SiPixelTemplate::interpolate
bool interpolate(int id, float cotalpha, float cotbeta, float locBz)
Definition: SiPixelTemplate.cc:844
SiPixelTemplate::sxone
float sxone()
rms for one pixel x-clusters
Definition: SiPixelTemplate.h:282
BYM1
#define BYM1
Definition: SiPixelTemplateDefs.h:20
BXM2
#define BXM2
Definition: SiPixelTemplateDefs.h:29
SiPixelTemplateReco::array_2d
boost::multi_array< float, 2 > array_2d
Definition: SiPixelTemplateReco.h:69
SiPixelTemplate::qscale
float qscale()
charge scaling factor
Definition: SiPixelTemplate.h:273
SiPixelTemplate::dxone
float dxone()
mean offset/correction for one pixel x-clusters
Definition: SiPixelTemplate.h:281
theVerboseLevel
static int theVerboseLevel
Definition: SiPixelTemplateReco.cc:65
SiPixelTemplate::chi2yavg
float chi2yavg(int i)
average y chi^2 in 4 charge bins
Definition: SiPixelTemplate.h:354
SiPixelTemplate::vavilov_pars
void vavilov_pars(double &mpv, double &sigma, double &kappa)
Definition: SiPixelTemplate.cc:2905
SiPixelTemplate::xtemp
void xtemp(int fxbin, int lxbin, float xtemplate[41][BXSIZE])
Definition: SiPixelTemplate.cc:1703
BHX
#define BHX
Definition: SiPixelTemplateDefs.h:27
SiPixelTemplate::yflcorr
float yflcorr(int binq, float qfly)
Definition: SiPixelTemplate.cc:1521
j
int j
Definition: DBlmapReader.cc:9
SiPixelTemplate::xsize
float xsize()
pixel x-size (microns)
Definition: SiPixelTemplate.h:447
BYM2
#define BYM2
Definition: SiPixelTemplateDefs.h:21
f
double f[11][100]
Definition: MuScleFitUtils.cc:79
SiPixelTemplate::ytemp
void ytemp(int fybin, int lybin, float ytemplate[41][BYSIZE])
Definition: SiPixelTemplate.cc:1625
SiPixelTemplate::sxtwo
float sxtwo()
rms for one double-pixel x-clusters
Definition: SiPixelTemplate.h:284
BYM3
#define BYM3
Definition: SiPixelTemplateDefs.h:22
SiPixelTemplate::dytwo
float dytwo()
mean offset/correction for one double-pixel y-clusters
Definition: SiPixelTemplate.h:278
LOGDEBUG
#define LOGDEBUG(x)
Definition: SiPixelTemplateReco.cc:64
SiPixelTemplate::s50
float s50()
1/2 of the pixel threshold signal in electrons
Definition: SiPixelTemplate.h:274
gen::k
int k[5][pyjets_maxn]
Definition: Pythia6Hadronizer.cc:86
TYSIZE
#define TYSIZE
Definition: SiPixelTemplateDefs.h:14
maxpix
static const int maxpix
Definition: SiPixelLorentzAngle.h:49
LOGERROR
#define LOGERROR(x)
Definition: SiPixelTemplateReco.cc:63
SiPixelTemplate::syone
float syone()
rms for one pixel y-clusters
Definition: SiPixelTemplate.h:277
SiPixelTemplate::chi2yavgone
float chi2yavgone()
//!&lt; average y chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:438
SiPixelTemplate::qavg
float qavg()
average cluster charge for this set of track angles
Definition: SiPixelTemplate.h:271
BHY
#define BHY
Definition: SiPixelTemplateDefs.h:19
SiPixelTemplate::sxmax
float sxmax()
average pixel signal for x-projection of cluster
Definition: SiPixelTemplate.h:280
SiPixelTemplate::xsigma2
void xsigma2(int fxpix, int lxpix, float sxthr, float xsum[BXSIZE], float xsig2[BXSIZE])
Definition: SiPixelTemplate.cc:1414
VVIObj
Definition: VVIObj.h:23
SiPixelTemplate::chi2xavgone
float chi2xavgone()
//!&lt; average x chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:440
SiPixelTemplate::pixmax
float pixmax()
maximum pixel charge
Definition: SiPixelTemplate.h:272
SiPixelTemplate::ysigma2
void ysigma2(int fypix, int lypix, float sythr, float ysum[BYSIZE], float ysig2[BYSIZE])
Gamma
dbl * Gamma
Definition: mlp_gen.cc:38
if
perl if(1 lt scalar(@::datatypes))
Definition: edlooper.cc:31
ENDL
#define ENDL
Definition: SiPixelTemplateReco.cc:66
SiPixelTemplate::chi2xavg
float chi2xavg(int i)
averaage x chi^2 in 4 charge bins
Definition: SiPixelTemplate.h:368
SiPixelTemplate::dyone
float dyone()
mean offset/correction for one pixel y-clusters
Definition: SiPixelTemplate.h:276
SiPixelTemplate::xavg
float xavg(int i)
average x-bias of reconstruction binned in 4 charge bins
Definition: SiPixelTemplate.h:326
SiPixelTemplate::dxtwo
float dxtwo()
mean offset/correction for one double-pixel x-clusters
Definition: SiPixelTemplate.h:283
SiPixelTemplate::ysize
float ysize()
pixel y-size (microns)
Definition: SiPixelTemplate.h:448
int SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
float  locBz,
array_2d &  cluster,
std::vector< bool > &  ydouble,
std::vector< bool > &  xdouble,
SiPixelTemplate &  templ,
float &  yrec,
float &  sigmay,
float &  proby,
float &  xrec,
float &  sigmax,
float &  probx,
int &  qbin,
int  speed,
float &  probQ 
)

Reconstruct the best estimate of the hit position for pixel clusters.

Parameters
id- (input) identifier of the template to use
cotalpha- (input) the cotangent of the alpha track angle (see CMS IN 2004/014)
cotbeta- (input) the cotangent of the beta track angle (see CMS IN 2004/014)
locBz- (input) the sign of the local B_z field for FPix (usually B_z<0 when cot(beta)>0 and B_z>0 when cot(beta)<0
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].
ydouble- (input) STL vector of 21 element array to flag a double-pixel
xdouble- (input) STL vector of 7 element array to flag a double-pixel
templ- (input) the template used in the reconstruction
yrec- (output) best estimate of y-coordinate of hit in microns
sigmay- (output) best estimate of uncertainty on yrec in microns
proby- (output) probability describing goodness-of-fit for y-reco
xrec- (output) best estimate of x-coordinate of hit in microns
sigmax- (output) best estimate of uncertainty on xrec in microns
probx- (output) probability describing goodness-of-fit for x-reco
qbin- (output) index (0-4) describing the charge of the cluster [0: 1.5<Q/Qavg, 1: 1<Q/Qavg<1.5, 2: 0.85<Q/Qavg<1, 3: 0.95Qmin<Q<0.85Qavg, 4: Q<0.95Qmin]
speed- (input) switch (-1-4) trading speed vs robustness -1 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4), calculates Q probability 0 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4) 1 faster, searches reduced 25 bin range (no big pix) + 33 bins (big pix at ends) at full density 2 faster yet, searches same range as 1 but at 1/2 density 3 fastest, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster) 4 fastest w/ Q prob, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster), calculates Q probability
probQ- (output) the Vavilov-distribution-based cluster charge probability

Definition at line 1030 of file SiPixelTemplateReco.cc.

References PixelTempReco2D().

1036 {
1037  // Local variables
1038  const bool deadpix = false;
1039  std::vector<std::pair<int, int> > zeropix;
1040 
1041  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, ydouble, xdouble, templ,
1042  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1043 
1044 } // PixelTempReco2D
ExpressReco_HICollisions_FallBack.speed
tuple speed
Definition: ExpressReco_HICollisions_FallBack.py:23342
SiPixelTemplateReco::PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ)
Definition: SiPixelTemplateReco.cc:120
int SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
array_2d &  cluster,
std::vector< bool > &  ydouble,
std::vector< bool > &  xdouble,
SiPixelTemplate &  templ,
float &  yrec,
float &  sigmay,
float &  proby,
float &  xrec,
float &  sigmax,
float &  probx,
int &  qbin,
int  speed,
float &  probQ 
)

Reconstruct the best estimate of the hit position for pixel clusters.

Parameters
id- (input) identifier of the template to use
cotalpha- (input) the cotangent of the alpha track angle (see CMS IN 2004/014)
cotbeta- (input) the cotangent of the beta track angle (see CMS IN 2004/014)
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].
ydouble- (input) STL vector of 21 element array to flag a double-pixel
xdouble- (input) STL vector of 7 element array to flag a double-pixel
templ- (input) the template used in the reconstruction
yrec- (output) best estimate of y-coordinate of hit in microns
sigmay- (output) best estimate of uncertainty on yrec in microns
proby- (output) probability describing goodness-of-fit for y-reco
xrec- (output) best estimate of x-coordinate of hit in microns
sigmax- (output) best estimate of uncertainty on xrec in microns
probx- (output) probability describing goodness-of-fit for x-reco
qbin- (output) index (0-4) describing the charge of the cluster [0: 1.5<Q/Qavg, 1: 1<Q/Qavg<1.5, 2: 0.85<Q/Qavg<1, 3: 0.95Qmin<Q<0.85Qavg, 4: Q<0.95Qmin]
speed- (input) switch (-1-4) trading speed vs robustness -1 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4), calculates Q probability 0 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4) 1 faster, searches reduced 25 bin range (no big pix) + 33 bins (big pix at ends) at full density 2 faster yet, searches same range as 1 but at 1/2 density 3 fastest, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster) 4 fastest w/ Q prob, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster), calculates Q probability
probQ- (output) the Vavilov-distribution-based cluster charge probability

Definition at line 1074 of file SiPixelTemplateReco.cc.

References PixelTempReco2D().

1080 {
1081  // Local variables
1082  const bool deadpix = false;
1083  std::vector<std::pair<int, int> > zeropix;
1084  float locBz = -1.;
1085  if(cotbeta < 0.) {locBz = -locBz;}
1086 
1087  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, ydouble, xdouble, templ,
1088  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1089 
1090 } // PixelTempReco2D
ExpressReco_HICollisions_FallBack.speed
tuple speed
Definition: ExpressReco_HICollisions_FallBack.py:23342
SiPixelTemplateReco::PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ)
Definition: SiPixelTemplateReco.cc:120
int SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
array_2d &  cluster,
std::vector< bool > &  ydouble,
std::vector< bool > &  xdouble,
SiPixelTemplate &  templ,
float &  yrec,
float &  sigmay,
float &  proby,
float &  xrec,
float &  sigmax,
float &  probx,
int &  qbin,
int  speed 
)

Reconstruct the best estimate of the hit position for pixel clusters.

Parameters
id- (input) identifier of the template to use
cotalpha- (input) the cotangent of the alpha track angle (see CMS IN 2004/014)
cotbeta- (input) the cotangent of the beta track angle (see CMS IN 2004/014)
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].
ydouble- (input) STL vector of 21 element array to flag a double-pixel
xdouble- (input) STL vector of 7 element array to flag a double-pixel
templ- (input) the template used in the reconstruction
yrec- (output) best estimate of y-coordinate of hit in microns
sigmay- (output) best estimate of uncertainty on yrec in microns
proby- (output) probability describing goodness-of-fit for y-reco
xrec- (output) best estimate of x-coordinate of hit in microns
sigmax- (output) best estimate of uncertainty on xrec in microns
probx- (output) probability describing goodness-of-fit for x-reco
qbin- (output) index (0-4) describing the charge of the cluster [0: 1.5<Q/Qavg, 1: 1<Q/Qavg<1.5, 2: 0.85<Q/Qavg<1, 3: 0.95Qmin<Q<0.85Qavg, 4: Q<0.95Qmin]
speed- (input) switch (0-3) trading speed vs robustness 0 totally bombproof, searches the entire 41 bin range at full density (equiv to V2_4) 1 faster, searches reduced 25 bin range (no big pix) + 33 bins (big pix at ends) at full density 2 faster yet, searches same range as 1 but at 1/2 density 3 fastest, searches same range as 1 but at 1/4 density (no big pix) and 1/2 density (big pix in cluster)

Definition at line 1118 of file SiPixelTemplateReco.cc.

References PixelTempReco2D().

1123 {
1124  // Local variables
1125  const bool deadpix = false;
1126  std::vector<std::pair<int, int> > zeropix;
1127  float locBz = -1.;
1128  if(cotbeta < 0.) {locBz = -locBz;}
1129  float probQ;
1130  if(speed < 0) speed = 0;
1131  if(speed > 3) speed = 3;
1132 
1133  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, ydouble, xdouble, templ,
1134  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1135 
1136 } // PixelTempReco2D
ExpressReco_HICollisions_FallBack.speed
tuple speed
Definition: ExpressReco_HICollisions_FallBack.py:23342
SiPixelTemplateReco::PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, array_2d &cluster, std::vector< bool > &ydouble, std::vector< bool > &xdouble, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ)
Definition: SiPixelTemplateReco.cc:120
int SiPixelTemplateReco::PixelTempSplit ( int  id,
bool  fpix,
float  cotalpha,
float  cotbeta,
array_2d  cluster,
std::vector< bool >  ydouble,
std::vector< bool >  xdouble,
SiPixelTemplate &  templ,
float &  yrec1,
float &  yrec2,
float &  sigmay,
float &  proby,
float &  xrec1,
float &  xrec2,
float &  sigmax,
float &  probx,
int &  qbin,
bool  deadpix,
std::vector< std::pair< int, int > >  zeropix 
)

Reconstruct the best estimate of the hit positions for pixel clusters.

Parameters
id- (input) identifier of the template to use
fpix- (input) logical input indicating whether to use FPix templates (true) or Barrel templates (false)
cotalpha- (input) the cotangent of the alpha track angle (see CMS IN 2004/014)
cotbeta- (input) the cotangent of the beta track angle (see CMS IN 2004/014)
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].
ydouble- (input) STL vector of 21 element array to flag a double-pixel
xdouble- (input) STL vector of 7 element array to flag a double-pixel
templ- (input) the template used in the reconstruction
yrec1- (output) best estimate of first y-coordinate of hit in microns
yrec2- (output) best estimate of second y-coordinate of hit in microns
sigmay- (output) best estimate of uncertainty on yrec in microns
proby- (output) probability describing goodness-of-fit for y-reco
xrec1- (output) best estimate of first x-coordinate of hit in microns
xrec2- (output) best estimate of second x-coordinate of hit in microns
sigmax- (output) best estimate of uncertainty on xrec in microns
probx- (output) probability describing goodness-of-fit for x-reco
qbin- (output) index (0-4) describing the charge of the cluster [0: 1.5<Q/Qavg, 1: 1<Q/Qavg<1.5, 2: 0.85<Q/Qavg<1, 3: 0.95Qmin<Q<0.85Qavg, 4: Q<0.95Qmin]
deadpix- (input) bool to indicate that there are dead pixels to be included in the analysis
zeropix- (input) vector of index pairs pointing to the dead pixels

Definition at line 70 of file SiPixelTemplateSplit.cc.

References BHX, BHY, BXM1, BXM2, BXSIZE, BYM1, BYM2, BYM3, BYSIZE, SiPixelTemplate::chi2xavg(), SiPixelTemplate::chi2xavgone(), SiPixelTemplate::chi2xminone(), SiPixelTemplate::chi2yavg(), SiPixelTemplate::chi2yavgone(), SiPixelTemplate::chi2yminone(), delta, SiPixelTemplate::dxone(), SiPixelTemplate::dxtwo(), SiPixelTemplate::dyone(), SiPixelTemplate::dytwo(), ENDL, Gamma, i, SiPixelTemplate::interpolate(), j, gen::k, LOGDEBUG, LOGERROR, max(), maxpix, min, SiPixelTemplate::pixmax(), SiPixelTemplate::qavg(), SiPixelTemplate::qmin(), SiPixelTemplate::qscale(), SiPixelTemplate::s50(), ExpressReco_HICollisions_FallBack::sigma, python.multivaluedict::sort(), SiPixelTemplate::sxmax(), SiPixelTemplate::sxone(), SiPixelTemplate::sxtwo(), SiPixelTemplate::symax(), SiPixelTemplate::syone(), SiPixelTemplate::sytwo(), theVerboseLevel, TXSIZE, TYSIZE, SiPixelTemplate::xavgc2m(), SiPixelTemplate::xrmsc2m(), SiPixelTemplate::xsigma2(), SiPixelTemplate::xtemp3d(), SiPixelTemplate::yavgc2m(), SiPixelTemplate::yrmsc2m(), SiPixelTemplate::ysigma2(), and SiPixelTemplate::ytemp3d().

Referenced by PixelCPETemplateReco::localPosition(), and PixelTempSplit().

76 {
77  // Local variables
78  static int i, j, k, minbin, binq, midpix, fypix, nypix, lypix, logypx;
79  static int fxpix, nxpix, lxpix, logxpx, shifty, shiftx, nyzero[TYSIZE];
80  static int nclusx, nclusy;
81  static int nybin, ycbin, nxbin, xcbin, minbinj, minbink;
82  static float sythr, sxthr, delta, sigma, sigavg, pseudopix, qscale;
83  static float ss2, ssa, sa2, rat, fq, qtotal, qpixel;
84  static float originx, originy, bias, maxpix, minmax;
85  static double chi2x, meanx, chi2y, meany, chi2ymin, chi2xmin, chi21max;
86  static double hchi2, hndof;
87  static float ysum[BYSIZE], xsum[BXSIZE], ysort[BYSIZE], xsort[BXSIZE];
88  static float ysig2[BYSIZE], xsig2[BXSIZE];
89  static bool yd[BYSIZE], xd[BXSIZE], anyyd, anyxd;
90  const float ysize={150.}, xsize={100.}, sqrt2={2.};
91  const float probmin={1.110223e-16};
92 // const float sqrt2={1.41421356};
93 
94 // The minimum chi2 for a valid one pixel cluster = pseudopixel contribution only
95 
96  const double mean1pix={0.100}, chi21min={0.160};
97 
98 // First, interpolate the template needed to analyze this cluster
99 // check to see of the track direction is in the physical range of the loaded template
100 
101  if(!templ.interpolate(id, fpix, cotalpha, cotbeta)) {
102  LOGDEBUG("SiPixelTemplateReco") << "input cluster direction cot(alpha) = " << cotalpha << ", cot(beta) = " << cotbeta << " is not within the acceptance of fpix = "
103  << fpix << ", template ID = " << id << ", no reconstruction performed" << ENDL;
104  return 20;
105  }
106 
107 // Define size of pseudopixel
108 
109  pseudopix = templ.s50();
110 
111 // Get charge scaling factor
112 
113  qscale = templ.qscale();
114 
115 // Check that the cluster container is (up to) a 7x21 matrix and matches the dimensions of the double pixel flags
116 
117  if(cluster.num_dimensions() != 2) {
118  LOGERROR("SiPixelTemplateReco") << "input cluster container (BOOST Multiarray) has wrong number of dimensions" << ENDL;
119  return 3;
120  }
121  nclusx = (int)cluster.shape()[0];
122  nclusy = (int)cluster.shape()[1];
123  if(nclusx != (int)xdouble.size()) {
124  LOGERROR("SiPixelTemplateReco") << "input cluster container x-size is not equal to double pixel flag container size" << ENDL;
125  return 4;
126  }
127  if(nclusy != (int)ydouble.size()) {
128  LOGERROR("SiPixelTemplateReco") << "input cluster container y-size is not equal to double pixel flag container size" << ENDL;
129  return 5;
130  }
131 
132 // enforce maximum size
133 
134  if(nclusx > TXSIZE) {nclusx = TXSIZE;}
135  if(nclusy > TYSIZE) {nclusy = TYSIZE;}
136 
137 // First, rescale all pixel charges
138 
139  for(i=0; i<nclusy; ++i) {
140  for(j=0; j<nclusx; ++j) {
141  if(cluster[j][i] > 0) {cluster[j][i] *= qscale;}
142  }
143  }
144 
145 // Next, sum the total charge and "decapitate" big pixels
146 
147  qtotal = 0.;
148  minmax = 2.*templ.pixmax();
149  for(i=0; i<nclusy; ++i) {
150  maxpix = minmax;
151  if(ydouble[i]) {maxpix *=2.;}
152  for(j=0; j<nclusx; ++j) {
153  qtotal += cluster[j][i];
154  if(cluster[j][i] > maxpix) {cluster[j][i] = maxpix;}
155  }
156  }
157 
158 // Do the cluster repair here
159 
160  if(deadpix) {
161  fypix = BYM3; lypix = -1;
162  for(i=0; i<nclusy; ++i) {
163  ysum[i] = 0.; nyzero[i] = 0;
164 // Do preliminary cluster projection in y
165  for(j=0; j<nclusx; ++j) {
166  ysum[i] += cluster[j][i];
167  }
168  if(ysum[i] > 0.) {
169 // identify ends of cluster to determine what the missing charge should be
170  if(i < fypix) {fypix = i;}
171  if(i > lypix) {lypix = i;}
172  }
173  }
174 
175 // Now loop over dead pixel list and "fix" everything
176 
177 //First see if the cluster ends are redefined and that we have only one dead pixel per column
178 
179  std::vector<std::pair<int, int> >::const_iterator zeroIter = zeropix.begin(), zeroEnd = zeropix.end();
180  for ( ; zeroIter != zeroEnd; ++zeroIter ) {
181  i = zeroIter->second;
182  if(i<0 || i>TYSIZE-1) {LOGERROR("SiPixelTemplateReco") << "dead pixel column y-index " << i << ", no reconstruction performed" << ENDL;
183  return 11;}
184 
185 // count the number of dead pixels in each column
186  ++nyzero[i];
187 // allow them to redefine the cluster ends
188  if(i < fypix) {fypix = i;}
189  if(i > lypix) {lypix = i;}
190  }
191 
192  nypix = lypix-fypix+1;
193 
194 // Now adjust the charge in the dead pixels to sum to 0.5*truncation value in the end columns and the truncation value in the interior columns
195 
196  for (zeroIter = zeropix.begin(); zeroIter != zeroEnd; ++zeroIter ) {
197  i = zeroIter->second; j = zeroIter->first;
198  if(j<0 || j>TXSIZE-1) {LOGERROR("SiPixelTemplateReco") << "dead pixel column x-index " << j << ", no reconstruction performed" << ENDL;
199  return 12;}
200  if((i == fypix || i == lypix) && nypix > 1) {maxpix = templ.symax()/2.;} else {maxpix = templ.symax();}
201  if(ydouble[i]) {maxpix *=2.;}
202  if(nyzero[i] > 0 && nyzero[i] < 3) {qpixel = (maxpix - ysum[i])/(float)nyzero[i];} else {qpixel = 1.;}
203  if(qpixel < 1.) {qpixel = 1.;}
204  cluster[j][i] = qpixel;
205 // Adjust the total cluster charge to reflect the charge of the "repaired" cluster
206  qtotal += qpixel;
207  }
208 // End of cluster repair section
209  }
210 
211 // Next, make y-projection of the cluster and copy the double pixel flags into a 25 element container
212 
213  for(i=0; i<BYSIZE; ++i) { ysum[i] = 0.; yd[i] = false;}
214  k=0;
215  anyyd = false;
216  for(i=0; i<nclusy; ++i) {
217  for(j=0; j<nclusx; ++j) {
218  ysum[k] += cluster[j][i];
219  }
220 
221 // If this is a double pixel, put 1/2 of the charge in 2 consective single pixels
222 
223  if(ydouble[i]) {
224  ysum[k] /= 2.;
225  ysum[k+1] = ysum[k];
226  yd[k] = true;
227  yd[k+1] = false;
228  k=k+2;
229  anyyd = true;
230  } else {
231  yd[k] = false;
232  ++k;
233  }
234  if(k > BYM1) {break;}
235  }
236 
237 // Next, make x-projection of the cluster and copy the double pixel flags into an 11 element container
238 
239  for(i=0; i<BXSIZE; ++i) { xsum[i] = 0.; xd[i] = false;}
240  k=0;
241  anyxd = false;
242  for(j=0; j<nclusx; ++j) {
243  for(i=0; i<nclusy; ++i) {
244  xsum[k] += cluster[j][i];
245  }
246 
247 // If this is a double pixel, put 1/2 of the charge in 2 consective single pixels
248 
249  if(xdouble[j]) {
250  xsum[k] /= 2.;
251  xsum[k+1] = xsum[k];
252  xd[k]=true;
253  xd[k+1]=false;
254  k=k+2;
255  anyxd = true;
256  } else {
257  xd[k]=false;
258  ++k;
259  }
260  if(k > BXM1) {break;}
261  }
262 
263 // next, identify the y-cluster ends, count total pixels, nypix, and logical pixels, logypx
264 
265  fypix=-1;
266  nypix=0;
267  lypix=0;
268  logypx=0;
269  for(i=0; i<BYSIZE; ++i) {
270  if(ysum[i] > 0.) {
271  if(fypix == -1) {fypix = i;}
272  if(!yd[i]) {
273  ysort[logypx] = ysum[i];
274  ++logypx;
275  }
276  ++nypix;
277  lypix = i;
278  }
279  }
280 
281 // Make sure cluster is continuous
282 
283  if((lypix-fypix+1) != nypix || nypix == 0) {
284  LOGDEBUG("SiPixelTemplateReco") << "y-length of pixel cluster doesn't agree with number of pixels above threshold" << ENDL;
285  if (theVerboseLevel > 2) {
286  LOGDEBUG("SiPixelTemplateReco") << "ysum[] = ";
287  for(i=0; i<BYSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << ysum[i] << ", ";}
288  LOGDEBUG("SiPixelTemplateReco") << ysum[BYSIZE-1] << ENDL;
289  }
290 
291  return 1;
292  }
293 
294 // If cluster is longer than max template size, technique fails
295 
296  if(nypix > TYSIZE) {
297  LOGDEBUG("SiPixelTemplateReco") << "y-length of pixel cluster is larger than maximum template size" << ENDL;
298  if (theVerboseLevel > 2) {
299  LOGDEBUG("SiPixelTemplateReco") << "ysum[] = ";
300  for(i=0; i<BYSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << ysum[i] << ", ";}
301  LOGDEBUG("SiPixelTemplateReco") << ysum[BYSIZE-1] << ENDL;
302  }
303 
304  return 6;
305  }
306 
307 // next, center the cluster on pixel 12 if necessary
308 
309  midpix = (fypix+lypix)/2;
310  shifty = BHY - midpix;
311  if(shifty > 0) {
312  for(i=lypix; i>=fypix; --i) {
313  ysum[i+shifty] = ysum[i];
314  ysum[i] = 0.;
315  yd[i+shifty] = yd[i];
316  yd[i] = false;
317  }
318  } else if (shifty < 0) {
319  for(i=fypix; i<=lypix; ++i) {
320  ysum[i+shifty] = ysum[i];
321  ysum[i] = 0.;
322  yd[i+shifty] = yd[i];
323  yd[i] = false;
324  }
325  }
326  lypix +=shifty;
327  fypix +=shifty;
328 
329 // If the cluster boundaries are OK, add pesudopixels, otherwise quit
330 
331  if(fypix > 1 && fypix < BYM2) {
332  ysum[fypix-1] = pseudopix;
333  ysum[fypix-2] = 0.2*pseudopix;
334  } else {return 8;}
335  if(lypix > 1 && lypix < BYM2) {
336  ysum[lypix+1] = pseudopix;
337  ysum[lypix+2] = 0.2*pseudopix;
338  } else {return 8;}
339 
340 // finally, determine if pixel[0] is a double pixel and make an origin correction if it is
341 
342  if(ydouble[0]) {
343  originy = -0.5;
344  } else {
345  originy = 0.;
346  }
347 
348 // next, identify the x-cluster ends, count total pixels, nxpix, and logical pixels, logxpx
349 
350  fxpix=-1;
351  nxpix=0;
352  lxpix=0;
353  logxpx=0;
354  for(i=0; i<BXSIZE; ++i) {
355  if(xsum[i] > 0.) {
356  if(fxpix == -1) {fxpix = i;}
357  if(!xd[i]) {
358  xsort[logxpx] = xsum[i];
359  ++logxpx;
360  }
361  ++nxpix;
362  lxpix = i;
363  }
364  }
365 
366 // Make sure cluster is continuous
367 
368  if((lxpix-fxpix+1) != nxpix) {
369 
370  LOGDEBUG("SiPixelTemplateReco") << "x-length of pixel cluster doesn't agree with number of pixels above threshold" << ENDL;
371  if (theVerboseLevel > 2) {
372  LOGDEBUG("SiPixelTemplateReco") << "xsum[] = ";
373  for(i=0; i<BXSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << xsum[i] << ", ";}
374  LOGDEBUG("SiPixelTemplateReco") << ysum[BXSIZE-1] << ENDL;
375  }
376 
377  return 2;
378  }
379 
380 // If cluster is longer than max template size, technique fails
381 
382  if(nxpix > TXSIZE) {
383 
384  LOGDEBUG("SiPixelTemplateReco") << "x-length of pixel cluster is larger than maximum template size" << ENDL;
385  if (theVerboseLevel > 2) {
386  LOGDEBUG("SiPixelTemplateReco") << "xsum[] = ";
387  for(i=0; i<BXSIZE-1; ++i) {LOGDEBUG("SiPixelTemplateReco") << xsum[i] << ", ";}
388  LOGDEBUG("SiPixelTemplateReco") << ysum[BXSIZE-1] << ENDL;
389  }
390 
391  return 7;
392  }
393 
394 // next, center the cluster on pixel 5 if necessary
395 
396  midpix = (fxpix+lxpix)/2;
397  shiftx = BHX - midpix;
398  if(shiftx > 0) {
399  for(i=lxpix; i>=fxpix; --i) {
400  xsum[i+shiftx] = xsum[i];
401  xsum[i] = 0.;
402  xd[i+shiftx] = xd[i];
403  xd[i] = false;
404  }
405  } else if (shiftx < 0) {
406  for(i=fxpix; i<=lxpix; ++i) {
407  xsum[i+shiftx] = xsum[i];
408  xsum[i] = 0.;
409  xd[i+shiftx] = xd[i];
410  xd[i] = false;
411  }
412  }
413  lxpix +=shiftx;
414  fxpix +=shiftx;
415 
416 // If the cluster boundaries are OK, add pesudopixels, otherwise quit
417 
418  if(fxpix > 1 && fxpix <BXM2) {
419  xsum[fxpix-1] = pseudopix;
420  xsum[fxpix-2] = 0.2*pseudopix;
421  } else {return 9;}
422  if(lxpix > 1 && lxpix < BXM2) {
423  xsum[lxpix+1] = pseudopix;
424  xsum[lxpix+2] = 0.2*pseudopix;
425  } else {return 9;}
426 
427 // finally, determine if pixel[0] is a double pixel and make an origin correction if it is
428 
429  if(xdouble[0]) {
430  originx = -0.5;
431  } else {
432  originx = 0.;
433  }
434 
435 // uncertainty and final corrections depend upon total charge bin
436 
437  fq = qtotal/templ.qavg();
438  if(fq > 3.0) {
439  binq=0;
440  } else {
441  if(fq > 2.0) {
442  binq=1;
443  } else {
444  if(fq > 1.70) {
445  binq=2;
446  } else {
447  binq=3;
448  }
449  }
450  }
451 
452 // Return the charge bin via the parameter list unless the charge is too small (then flag it)
453 
454  qbin = binq;
455  if(!deadpix && qtotal < 1.9*templ.qmin()) {qbin = 5;} else {
456  if(!deadpix && qtotal < 1.9*templ.qmin(1)) {qbin = 4;}
457  }
458 
459  if (theVerboseLevel > 9) {
460  LOGDEBUG("SiPixelTemplateReco") <<
461  "ID = " << id << " FPix = " << fpix <<
462  " cot(alpha) = " << cotalpha << " cot(beta) = " << cotbeta <<
463  " nclusx = " << nclusx << " nclusy = " << nclusy << ENDL;
464  }
465 
466 
467 // Next, generate the 3d y- and x-templates
468 
469  array_3d ytemp;
470  templ.ytemp3d(nypix, ytemp);
471 
472  array_3d xtemp;
473  templ.xtemp3d(nxpix, xtemp);
474 
475 // Do the y-reconstruction first
476 
477 // retrieve the number of y-bins
478 
479  nybin = ytemp.shape()[0]; ycbin = nybin/2;
480 
481 // Modify the template if double pixels are present
482 
483  if(nypix > logypx) {
484  i=fypix;
485  while(i < lypix) {
486  if(yd[i] && !yd[i+1]) {
487  for(j=0; j<nybin; ++j) {
488  for(k=0; k<=j; ++k) {
489 
490 // Sum the adjacent cells and put the average signal in both
491 
492  sigavg = (ytemp[j][k][i] + ytemp[j][k][i+1])/2.;
493  ytemp[j][k][i] = sigavg;
494  ytemp[j][k][i+1] = sigavg;
495  }
496  }
497  i += 2;
498  } else {
499  ++i;
500  }
501  }
502  }
503 
504 // Define the maximum signal to allow before de-weighting a pixel
505 
506  sythr = 2.1*(templ.symax());
507 
508 // Make sure that there will be at least two pixels that are not de-weighted
509 
510  std::sort(&ysort[0], &ysort[logypx]);
511  if(logypx == 1) {sythr = 1.01*ysort[0];} else {
512  if (ysort[1] > sythr) { sythr = 1.01*ysort[1]; }
513  }
514 
515 // Evaluate pixel-by-pixel uncertainties (weights) for the templ analysis
516 
517  for(i=0; i<BYSIZE; ++i) { ysig2[i] = 0.;}
518  templ.ysigma2(fypix, lypix, sythr, ysum, ysig2);
519 
520 // Find the template bin that minimizes the Chi^2
521 
522  chi2ymin = 1.e15;
523  minbinj = -1;
524  minbink = -1;
525  for(j=0; j<nybin; ++j) {
526  for(k=0; k<=j; ++k) {
527  ss2 = 0.;
528  ssa = 0.;
529  sa2 = 0.;
530  for(i=fypix-2; i<=lypix+2; ++i) {
531  ss2 += ysum[i]*ysum[i]/ysig2[i];
532  ssa += ysum[i]*ytemp[j][k][i]/ysig2[i];
533  sa2 += ytemp[j][k][i]*ytemp[j][k][i]/ysig2[i];
534  }
535  rat=ssa/ss2;
536  if(rat <= 0.) {LOGERROR("SiPixelTemplateReco") << "illegal chi2ymin normalization = " << rat << ENDL; rat = 1.;}
537  chi2y=ss2-2.*ssa/rat+sa2/(rat*rat);
538  if(chi2y < chi2ymin) {
539  chi2ymin = chi2y;
540  minbinj = j;
541  minbink = k;
542  }
543  }
544  }
545 
546  if (theVerboseLevel > 9) {
547  LOGDEBUG("SiPixelTemplateReco") <<
548  "minbins " << minbinj << "," << minbink << " chi2ymin = " << chi2ymin << ENDL;
549  }
550 
551 // Do not apply final template pass to 1-pixel clusters (use calibrated offset)
552 
553  if(logypx == 1) {
554 
555  if(nypix ==1) {
556  delta = templ.dyone();
557  sigma = templ.syone();
558  } else {
559  delta = templ.dytwo();
560  sigma = templ.sytwo();
561  }
562 
563  yrec1 = 0.5*(fypix+lypix-2*shifty+2.*originy)*ysize-delta;
564  yrec2 = yrec1;
565 
566  if(sigma <= 0.) {
567  sigmay = 43.3;
568  } else {
569  sigmay = sigma;
570  }
571 
572 // Do probability calculation for one-pixel clusters
573 
574  chi21max = fmax(chi21min, (double)templ.chi2yminone());
575  chi2ymin -=chi21max;
576  if(chi2ymin < 0.) {chi2ymin = 0.;}
577  // proby = gsl_cdf_chisq_Q(chi2ymin, mean1pix);
578  meany = fmax(mean1pix, (double)templ.chi2yavgone());
579  hchi2 = chi2ymin/2.; hndof = meany/2.;
580  proby = 1. - TMath::Gamma(hndof, hchi2);
581 
582  } else {
583 
584 // For cluster > 1 pix, use chi^2 minimm to recontruct the two y-positions
585 
586 // at small eta, the templates won't actually work on two pixel y-clusters so just return the pixel centers
587 
588  if(logypx == 2 && fabsf(cotbeta) < 0.25) {
589  switch(nypix) {
590  case 2:
591 // Both pixels are small
592  yrec1 = (fypix-shifty+originy)*ysize;
593  yrec2 = (lypix-shifty+originy)*ysize;
594  sigmay = 43.3;
595  break;
596  case 3:
597 // One big pixel and one small pixel
598  if(yd[fypix]) {
599  yrec1 = (fypix+0.5-shifty+originy)*ysize;
600  yrec2 = (lypix-shifty+originy)*ysize;
601  sigmay = 43.3;
602  } else {
603  yrec1 = (fypix-shifty+originy)*ysize;
604  yrec2 = (lypix-0.5-shifty+originy)*ysize;
605  sigmay = 65.;
606  }
607  break;
608  case 4:
609 // Two big pixels
610  yrec1 = (fypix+0.5-shifty+originy)*ysize;
611  yrec2 = (lypix-0.5-shifty+originy)*ysize;
612  sigmay = 86.6;
613  break;
614  default:
615 // Something is screwy ...
616  LOGERROR("SiPixelTemplateReco") << "weird problem: logical y-pixels = " << logypx << ", total ysize in normal pixels = " << nypix << ENDL;
617  return 10;
618  }
619  } else {
620 
621 // uncertainty and final correction depend upon charge bin
622 
623  bias = templ.yavgc2m(binq);
624  yrec1 = (0.125*(minbink-ycbin)+BHY-(float)shifty+originy)*ysize - bias;
625  yrec2 = (0.125*(minbinj-ycbin)+BHY-(float)shifty+originy)*ysize - bias;
626  sigmay = sqrt2*templ.yrmsc2m(binq);
627 
628  }
629 
630 // Do goodness of fit test in y
631 
632  if(chi2ymin < 0.0) {chi2ymin = 0.0;}
633  meany = 2.*templ.chi2yavg(binq);
634  if(meany < 0.01) {meany = 0.01;}
635 // gsl function that calculates the chi^2 tail prob for non-integral dof
636 // proby = gsl_cdf_chisq_Q(chi2y, meany);
637 // proby = ROOT::Math::chisquared_cdf_c(chi2y, meany);
638  hchi2 = chi2ymin/2.; hndof = meany/2.;
639  proby = 1. - TMath::Gamma(hndof, hchi2);
640  }
641 
642 // Do the x-reconstruction next
643 
644 // retrieve the number of x-bins
645 
646  nxbin = xtemp.shape()[0]; xcbin = nxbin/2;
647 
648 // Modify the template if double pixels are present
649 
650  if(nxpix > logxpx) {
651  i=fxpix;
652  while(i < lxpix) {
653  if(xd[i] && !xd[i+1]) {
654  for(j=0; j<nxbin; ++j) {
655  for(k=0; k<=j; ++k) {
656 
657 // Sum the adjacent cells and put the average signal in both
658 
659  sigavg = (xtemp[j][k][i] + xtemp[j][k][i+1])/2.;
660  xtemp[j][k][i] = sigavg;
661  xtemp[j][k][i+1] = sigavg;
662  }
663  }
664  i += 2;
665  } else {
666  ++i;
667  }
668  }
669  }
670 
671 // Define the maximum signal to allow before de-weighting a pixel
672 
673  sxthr = 2.1*templ.sxmax();
674 
675 // Make sure that there will be at least two pixels that are not de-weighted
676  std::sort(&xsort[0], &xsort[logxpx]);
677  if(logxpx == 1) {sxthr = 1.01*xsort[0];} else {
678  if (xsort[1] > sxthr) { sxthr = 1.01*xsort[1]; }
679  }
680 
681 // Evaluate pixel-by-pixel uncertainties (weights) for the templ analysis
682 
683  for(i=0; i<BYSIZE; ++i) { xsig2[i] = 0.; }
684  templ.xsigma2(fxpix, lxpix, sxthr, xsum, xsig2);
685 
686 // Find the template bin that minimizes the Chi^2
687 
688  chi2xmin = 1.e15;
689  minbinj = -1;
690  minbink = -1;
691  for(j=0; j<nxbin; ++j) {
692  for(k=0; k<=j; ++k) {
693  ss2 = 0.;
694  ssa = 0.;
695  sa2 = 0.;
696  for(i=fxpix-2; i<=lxpix+2; ++i) {
697  ss2 += xsum[i]*xsum[i]/xsig2[i];
698  ssa += xsum[i]*xtemp[j][k][i]/xsig2[i];
699  sa2 += xtemp[j][k][i]*xtemp[j][k][i]/xsig2[i];
700  }
701  rat=ssa/ss2;
702  if(rat <= 0.) {LOGERROR("SiPixelTemplateReco") << "illegal chi2ymin normalization = " << rat << ENDL; rat = 1.;}
703  chi2x=ss2-2.*ssa/rat+sa2/(rat*rat);
704  if(chi2x < chi2xmin) {
705  chi2xmin = chi2x;
706  minbinj = j;
707  minbink = k;
708  }
709  }
710  }
711 
712  if (theVerboseLevel > 9) {
713  LOGDEBUG("SiPixelTemplateReco") <<
714  "minbin " << minbin << " chi2xmin = " << chi2xmin << ENDL;
715  }
716 
717 // Do not apply final template pass to 1-pixel clusters (use calibrated offset)
718 
719  if(logxpx == 1) {
720 
721  if(nxpix ==1) {
722  delta = templ.dxone();
723  sigma = templ.sxone();
724  } else {
725  delta = templ.dxtwo();
726  sigma = templ.sxtwo();
727  }
728  xrec1 = 0.5*(fxpix+lxpix-2*shiftx+2.*originx)*xsize-delta;
729  xrec2 = xrec1;
730  if(sigma <= 0.) {
731  sigmax = 28.9;
732  } else {
733  sigmax = sigma;
734  }
735 
736 // Do probability calculation for one-pixel clusters
737 
738  chi21max = fmax(chi21min, (double)templ.chi2xminone());
739  chi2xmin -=chi21max;
740  if(chi2xmin < 0.) {chi2xmin = 0.;}
741  meanx = fmax(mean1pix, (double)templ.chi2xavgone());
742  hchi2 = chi2xmin/2.; hndof = meanx/2.;
743  probx = 1. - TMath::Gamma(hndof, hchi2);
744 
745  } else {
746 
747 // For cluster > 1 pix, use chi^2 minimm to recontruct the two x-positions
748 
749 // uncertainty and final correction depend upon charge bin
750 
751  bias = templ.xavgc2m(binq);
752  k = std::min(minbink, minbinj);
753  j = std::max(minbink, minbinj);
754  xrec1 = (0.125*(minbink-xcbin)+BHX-(float)shiftx+originx)*xsize - bias;
755  xrec2 = (0.125*(minbinj-xcbin)+BHX-(float)shiftx+originx)*xsize - bias;
756  sigmax = sqrt2*templ.xrmsc2m(binq);
757 
758 // Do goodness of fit test in y
759 
760  if(chi2xmin < 0.0) {chi2xmin = 0.0;}
761  meanx = 2.*templ.chi2xavg(binq);
762  if(meanx < 0.01) {meanx = 0.01;}
763  hchi2 = chi2xmin/2.; hndof = meanx/2.;
764  probx = 1. - TMath::Gamma(hndof, hchi2);
765  }
766 
767  // Don't return exact zeros for the probability
768 
769  if(proby < probmin) {proby = probmin;}
770  if(probx < probmin) {probx = probmin;}
771 
772  return 0;
773 } // PixelTempSplit
ExpressReco_HICollisions_FallBack.sigma
tuple sigma
Definition: ExpressReco_HICollisions_FallBack.py:4426
SiPixelTemplate::yrmsc2m
float yrmsc2m(int i)
1st pass chi2 min search: average y-rms of reconstruction binned in 4 charge bins ...
Definition: SiPixelTemplate.h:389
delta
dbl * delta
Definition: mlp_gen.cc:36
i
int i
Definition: DBlmapReader.cc:9
BXSIZE
#define BXSIZE
Definition: SiPixelTemplateDefs.h:26
SiPixelTemplate::chi2xminone
float chi2xminone()
//!&lt; minimum of x chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:441
SiPixelTemplate::ytemp3d
void ytemp3d(int nypix, array_3d &ytemplate)
Definition: SiPixelTemplate.cc:1779
theVerboseLevel
static int theVerboseLevel
Definition: SiPixelTemplateSplit.cc:33
SiPixelTemplate::symax
float symax()
average pixel signal for y-projection of cluster
Definition: SiPixelTemplate.h:275
BYSIZE
#define BYSIZE
Definition: SiPixelTemplateDefs.h:18
TXSIZE
#define TXSIZE
Definition: SiPixelTemplateDefs.h:23
LOGERROR
#define LOGERROR(x)
Definition: SiPixelTemplateSplit.cc:31
BXM1
#define BXM1
Definition: SiPixelTemplateDefs.h:28
SiPixelTemplate::qmin
float qmin()
minimum cluster charge for valid hit (keeps 99.9% of simulated hits)
Definition: SiPixelTemplate.h:285
min
#define min(a, b)
Definition: mlp_lapack.h:161
SiPixelTemplate::xavgc2m
float xavgc2m(int i)
1st pass chi2 min search: average x-bias of reconstruction binned in 4 charge bins ...
Definition: SiPixelTemplate.h:410
SiPixelTemplate::sytwo
float sytwo()
rms for one double-pixel y-clusters
Definition: SiPixelTemplate.h:279
ENDL
#define ENDL
Definition: SiPixelTemplateSplit.cc:34
SiPixelTemplate::chi2yminone
float chi2yminone()
//!&lt; minimum of y chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:439
SiPixelTemplate::interpolate
bool interpolate(int id, float cotalpha, float cotbeta, float locBz)
Definition: SiPixelTemplate.cc:844
SiPixelTemplate::sxone
float sxone()
rms for one pixel x-clusters
Definition: SiPixelTemplate.h:282
BYM1
#define BYM1
Definition: SiPixelTemplateDefs.h:20
BXM2
#define BXM2
Definition: SiPixelTemplateDefs.h:29
SiPixelTemplate::qscale
float qscale()
charge scaling factor
Definition: SiPixelTemplate.h:273
SiPixelTemplate::dxone
float dxone()
mean offset/correction for one pixel x-clusters
Definition: SiPixelTemplate.h:281
SiPixelTemplate::xtemp3d
void xtemp3d(int nxpix, array_3d &xtemplate)
Definition: SiPixelTemplate.cc:1871
max
const T & max(const T &a, const T &b)
Definition: MaterialBudgetTrackerHistos.cc:4
SiPixelTemplate::chi2yavg
float chi2yavg(int i)
average y chi^2 in 4 charge bins
Definition: SiPixelTemplate.h:354
BHX
#define BHX
Definition: SiPixelTemplateDefs.h:27
j
int j
Definition: DBlmapReader.cc:9
LOGDEBUG
#define LOGDEBUG(x)
Definition: SiPixelTemplateSplit.cc:32
BYM2
#define BYM2
Definition: SiPixelTemplateDefs.h:21
SiPixelTemplate::sxtwo
float sxtwo()
rms for one double-pixel x-clusters
Definition: SiPixelTemplate.h:284
BYM3
#define BYM3
Definition: SiPixelTemplateDefs.h:22
SiPixelTemplate::dytwo
float dytwo()
mean offset/correction for one double-pixel y-clusters
Definition: SiPixelTemplate.h:278
SiPixelTemplate::s50
float s50()
1/2 of the pixel threshold signal in electrons
Definition: SiPixelTemplate.h:274
gen::k
int k[5][pyjets_maxn]
Definition: Pythia6Hadronizer.cc:86
SiPixelTemplate::yavgc2m
float yavgc2m(int i)
1st pass chi2 min search: average y-bias of reconstruction binned in 4 charge bins ...
Definition: SiPixelTemplate.h:382
TYSIZE
#define TYSIZE
Definition: SiPixelTemplateDefs.h:14
maxpix
static const int maxpix
Definition: SiPixelLorentzAngle.h:49
SiPixelTemplate::syone
float syone()
rms for one pixel y-clusters
Definition: SiPixelTemplate.h:277
SiPixelTemplate::chi2yavgone
float chi2yavgone()
//!&lt; average y chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:438
SiPixelTemplate::qavg
float qavg()
average cluster charge for this set of track angles
Definition: SiPixelTemplate.h:271
BHY
#define BHY
Definition: SiPixelTemplateDefs.h:19
SiPixelTemplate::sxmax
float sxmax()
average pixel signal for x-projection of cluster
Definition: SiPixelTemplate.h:280
SiPixelTemplate::xsigma2
void xsigma2(int fxpix, int lxpix, float sxthr, float xsum[BXSIZE], float xsig2[BXSIZE])
Definition: SiPixelTemplate.cc:1414
SiPixelTemplate::chi2xavgone
float chi2xavgone()
//!&lt; average x chi^2 for 1 pixel clusters
Definition: SiPixelTemplate.h:440
SiPixelTemplate::pixmax
float pixmax()
maximum pixel charge
Definition: SiPixelTemplate.h:272
SiPixelTemplate::ysigma2
void ysigma2(int fypix, int lypix, float sythr, float ysum[BYSIZE], float ysig2[BYSIZE])
Gamma
dbl * Gamma
Definition: mlp_gen.cc:38
SiPixelTemplate::chi2xavg
float chi2xavg(int i)
averaage x chi^2 in 4 charge bins
Definition: SiPixelTemplate.h:368
SiPixelTemplate::dyone
float dyone()
mean offset/correction for one pixel y-clusters
Definition: SiPixelTemplate.h:276
SiPixelTemplate::dxtwo
float dxtwo()
mean offset/correction for one double-pixel x-clusters
Definition: SiPixelTemplate.h:283
array_3d
boost::multi_array< float, 3 > array_3d
Definition: SiPixelTemplate.h:79
SiPixelTemplate::xrmsc2m
float xrmsc2m(int i)
1st pass chi2 min search: average x-rms of reconstruction binned in 4 charge bins ...
Definition: SiPixelTemplate.h:417
int SiPixelTemplateReco::PixelTempSplit ( int  id,
bool  fpix,
float  cotalpha,
float  cotbeta,
array_2d  cluster,
std::vector< bool >  ydouble,
std::vector< bool >  xdouble,
SiPixelTemplate &  templ,
float &  yrec1,
float &  yrec2,
float &  sigmay,
float &  proby,
float &  xrec1,
float &  xrec2,
float &  sigmax,
float &  probx,
int &  qbin 
)

Reconstruct the best estimate of the hit positions for pixel clusters.

Parameters
id- (input) identifier of the template to use
fpix- (input) logical input indicating whether to use FPix templates (true) or Barrel templates (false)
cotalpha- (input) the cotangent of the alpha track angle (see CMS IN 2004/014)
cotbeta- (input) the cotangent of the beta track angle (see CMS IN 2004/014)
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].
ydouble- (input) STL vector of 21 element array to flag a double-pixel
xdouble- (input) STL vector of 7 element array to flag a double-pixel
templ- (input) the template used in the reconstruction
yrec1- (output) best estimate of first y-coordinate of hit in microns
yrec2- (output) best estimate of second y-coordinate of hit in microns
sigmay- (output) best estimate of uncertainty on yrec in microns
proby- (output) probability describing goodness-of-fit for y-reco
xrec1- (output) best estimate of first x-coordinate of hit in microns
xrec2- (output) best estimate of second x-coordinate of hit in microns
sigmax- (output) best estimate of uncertainty on xrec in microns
probx- (output) probability describing goodness-of-fit for x-reco
qbin- (output) index (0-4) describing the charge of the cluster [0: 1.5<Q/Qavg, 1: 1<Q/Qavg<1.5, 2: 0.85<Q/Qavg<1, 3: 0.95Qmin<Q<0.85Qavg, 4: Q<0.95Qmin]

Definition at line 800 of file SiPixelTemplateSplit.cc.

References PixelTempSplit().

805 {
806  // Local variables
807  const bool deadpix = false;
808  std::vector<std::pair<int, int> > zeropix;
809 
810  return SiPixelTemplateReco::PixelTempSplit(id, fpix, cotalpha, cotbeta, cluster, ydouble, xdouble, templ,
811  yrec1, yrec2, sigmay, proby, xrec1, xrec2, sigmax, probx, qbin, deadpix, zeropix);
812 
813 } // PixelTempSplit
SiPixelTemplateReco::PixelTempSplit
int PixelTempSplit(int id, bool fpix, float cotalpha, float cotbeta, array_2d cluster, std::vector< bool > ydouble, std::vector< bool > xdouble, SiPixelTemplate &templ, float &yrec1, float &yrec2, float &sigmay, float &proby, float &xrec1, float &xrec2, float &sigmax, float &probx, int &qbin, bool deadpix, std::vector< std::pair< int, int > > zeropix)
Definition: SiPixelTemplateSplit.cc:70
Generated for CMSSW Reference Manual by doxygen 1.8.5