CMS 3D CMS Logo

Classes | Functions
SiPixelTemplateReco Namespace Reference

Classes

struct  ClusMatrix
 

Functions

int PixelTempReco2D (int id, float cotalpha, float cotbeta, float locBz, float locBx, ClusMatrix &cluster, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ, int &nypix, int &nxpix)
 
int PixelTempReco2D (int id, float cotalpha, float cotbeta, float locBz, float locBx, ClusMatrix &cluster, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, float &probQ)
 
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,
float  locBx,
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 &  nypix,
int &  nxpix 
)

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)
locBz- (input) the sign of this quantity is used to determine whether to flip cot(beta)<0 quantities from cot(beta)>0 (FPix only) for Phase 0 FPix IP-related tracks, locBz < 0 for cot(beta) > 0 and locBz > 0 for cot(beta) < 0 for Phase 1 FPix IP-related tracks, see next comment
locBx- (input) the sign of this quantity is used to determine whether to flip cot(alpha/beta)<0 quantities from cot(alpha/beta)>0 (FPix only) for Phase 1 FPix IP-related tracks, locBx/locBz > 0 for cot(alpha) > 0 and locBx/locBz < 0 for cot(alpha) < 0 for Phase 1 FPix IP-related tracks, locBx > 0 for cot(beta) > 0 and locBx < 0 for cot(beta) < 0//!
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]. 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
nypix- (output) the projected y-size of the cluster
nxpix- (output) the projected x-size of the cluster

Definition at line 134 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(), SiPixelTemplate::cxtemp(), SiPixelTemplate::cytemp(), delta, SiPixelTemplate::dxone(), SiPixelTemplate::dxtwo(), SiPixelTemplate::dyone(), SiPixelTemplate::dytwo(), ENDL, Exception, f, SiPixelTemplate::fbin(), VVIObjF::fcn(), Gamma, mps_fire::i, createfilelist::int, SiPixelTemplate::interpolate(), gen::k, kappa, LOGDEBUG, LOGERROR, SiPixelTemplateReco::ClusMatrix::matrix, maxpix, SiPixelTemplateReco::ClusMatrix::mcol, SiPixelTemplateReco::ClusMatrix::mrow, SiPixelTemplate::pixmax(), SiPixelTemplate::qavg(), SiPixelTemplate::qmin(), SiPixelTemplate::qscale(), SiPixelTemplate::s50(), ctppsDiamondLocalTracks_cfi::sigma, 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().

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

References PixelTempReco2D().

1033 {
1034  // Local variables
1035  const bool deadpix = false;
1036  std::vector<std::pair<int, int> > zeropix;
1037  int nypix, nxpix;
1038 
1039  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, locBx, cluster, templ,
1040  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ, nypix, nxpix);
1041 
1042 } // PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, float locBx, ClusMatrix &cluster, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ, int &nypix, int &nxpix)
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 1073 of file SiPixelTemplateReco.cc.

References f, and PixelTempReco2D().

1078 {
1079  // Local variables
1080  const bool deadpix = false;
1081  std::vector<std::pair<int, int> > zeropix;
1082  int nypix, nxpix;
1083  float locBx, locBz;
1084  locBx = 1.f;
1085  if(cotbeta < 0.f) {locBx = -1.f;}
1086  locBz = locBx;
1087  if(cotalpha < 0.f) {locBz = -locBx;}
1088 
1089  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, locBx, cluster, templ,
1090  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ, nypix, nxpix);
1091 
1092 } // PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, float locBx, ClusMatrix &cluster, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ, int &nypix, int &nxpix)
double f[11][100]
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 1120 of file SiPixelTemplateReco.cc.

References f, and PixelTempReco2D().

1124 {
1125  // Local variables
1126  const bool deadpix = false;
1127  std::vector<std::pair<int, int> > zeropix;
1128  int nypix, nxpix;
1129  float locBx, locBz;
1130  locBx = 1.f;
1131  if(cotbeta < 0.f) {locBx = -1.f;}
1132  locBz = locBx;
1133  if(cotalpha < 0.f) {locBz = -locBx;}
1134  float probQ;
1135  if(speed < 0) speed = 0;
1136  if(speed > 3) speed = 3;
1137 
1138  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, locBx, cluster, templ,
1139  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ, nypix, nxpix);
1140 
1141 } // PixelTempReco2D
int PixelTempReco2D(int id, float cotalpha, float cotbeta, float locBz, float locBx, ClusMatrix &cluster, SiPixelTemplate &templ, float &yrec, float &sigmay, float &proby, float &xrec, float &sigmax, float &probx, int &qbin, int speed, bool deadpix, std::vector< std::pair< int, int > > &zeropix, float &probQ, int &nypix, int &nxpix)
double f[11][100]