CMS 3D CMS Logo

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

Classes

struct  ClusMatrix
 

Functions

int PixelTempReco2D (int id, float cotalpha, float cotbeta, float locBz, ClusMatrix &cluster, 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, ClusMatrix &cluster, 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, ClusMatrix &cluster, 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, ClusMatrix &cluster, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed)
 

Function Documentation

int SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
float  locBz,
ClusMatrix cluster,
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 125 of file SiPixelTemplateReco.cc.

References assert(), 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(), SiPixelTemplate::cxtemp(), SiPixelTemplate::cytemp(), delta, SiPixelTemplate::dxone(), SiPixelTemplate::dxtwo(), SiPixelTemplate::dyone(), SiPixelTemplate::dytwo(), ENDL, Exception, f, VVIObjF::fcn(), Gamma, i, SiPixelTemplate::interpolate(), j, relval_2017::k, kappa, LOGDEBUG, LOGERROR, SiPixelTemplateReco::ClusMatrix::matrix, maxpix, SiPixelTemplateReco::ClusMatrix::mcol, SiPixelTemplateReco::ClusMatrix::mrow, SiPixelTemplate::pixmax(), SiPixelTemplate::qavg(), SiPixelTemplate::qmin(), SiPixelTemplate::qscale(), SiPixelTemplate::s50(), SiPixelTemplate::sxmax(), SiPixelTemplate::sxone(), SiPixelTemplate::sxtwo(), SiPixelTemplate::symax(), SiPixelTemplate::syone(), SiPixelTemplate::sytwo(), theVerboseLevel, TXSIZE, TYSIZE, SiPixelTemplate::vavilov_pars(), SiPixelTemplate::xavg(), SiPixelTemplateReco::ClusMatrix::xdouble, SiPixelTemplate::xflcorr(), SiPixelTemplate::xrms(), SiPixelTemplate::xsigma2(), xsize, SiPixelTemplate::xsize(), SiPixelTemplate::xtemp(), SiPixelTemplate::yavg(), SiPixelTemplateReco::ClusMatrix::ydouble, SiPixelTemplate::yflcorr(), SiPixelTemplate::yrms(), SiPixelTemplate::ysigma2(), ysize, SiPixelTemplate::ysize(), and SiPixelTemplate::ytemp().

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

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

References PixelTempReco2D().

1010 {
1011  // Local variables
1012  const bool deadpix = false;
1013  std::vector<std::pair<int, int> > zeropix;
1014 
1015  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, templ,
1016  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1017 
1018 } // PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, ClusMatrix &cluster, 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 SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
ClusMatrix cluster,
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 1048 of file SiPixelTemplateReco.cc.

References PixelTempReco2D().

1053 {
1054  // Local variables
1055  const bool deadpix = false;
1056  std::vector<std::pair<int, int> > zeropix;
1057  float locBz = -1.f;
1058  if(cotbeta < 0.) {locBz = -locBz;}
1059 
1060  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, templ,
1061  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1062 
1063 } // PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, ClusMatrix &cluster, 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 SiPixelTemplateReco::PixelTempReco2D ( int  id,
float  cotalpha,
float  cotbeta,
ClusMatrix cluster,
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 1091 of file SiPixelTemplateReco.cc.

References PixelTempReco2D().

1095 {
1096  // Local variables
1097  const bool deadpix = false;
1098  std::vector<std::pair<int, int> > zeropix;
1099  float locBz = -1.f;
1100  if(cotbeta < 0.) {locBz = -locBz;}
1101  float probQ;
1102  if(speed < 0) speed = 0;
1103  if(speed > 3) speed = 3;
1104 
1105  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, templ,
1106  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1107 
1108 } // PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, ClusMatrix &cluster, 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)