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
 

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)
 

Typedef Documentation

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

Definition at line 70 of file SiPixelTemplateReco.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 125 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, edm::hlt::Exception, f, VVIObjF::fcn(), Gamma, i, SiPixelTemplate::interpolate(), j, gen::k, LOGDEBUG, LOGERROR, maxpix, SiPixelTemplate::pixmax(), SiPixelTemplate::qavg(), SiPixelTemplate::qmin(), SiPixelTemplate::qscale(), SiPixelTemplate::s50(), 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().

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

References PixelTempReco2D().

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

References PixelTempReco2D().

1078 {
1079  // Local variables
1080  const bool deadpix = false;
1081  std::vector<std::pair<int, int> > zeropix;
1082  float locBz = -1.f;
1083  if(cotbeta < 0.) {locBz = -locBz;}
1084 
1085  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, ydouble, xdouble, templ,
1086  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1087 
1088 } // PixelTempReco2D
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:125
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 1116 of file SiPixelTemplateReco.cc.

References PixelTempReco2D().

1121 {
1122  // Local variables
1123  const bool deadpix = false;
1124  std::vector<std::pair<int, int> > zeropix;
1125  float locBz = -1.f;
1126  if(cotbeta < 0.) {locBz = -locBz;}
1127  float probQ;
1128  if(speed < 0) speed = 0;
1129  if(speed > 3) speed = 3;
1130 
1131  return SiPixelTemplateReco::PixelTempReco2D(id, cotalpha, cotbeta, locBz, cluster, ydouble, xdouble, templ,
1132  yrec, sigmay, proby, xrec, sigmax, probx, qbin, speed, deadpix, zeropix, probQ);
1133 
1134 } // PixelTempReco2D
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:125
Generated for CMSSW Reference Manual by doxygen 1.8.5