27 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
39 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
43 #define LOGERROR(x) LogError(x)
44 #define LOGINFO(x) LogInfo(x)
50 #define LOGERROR(x) std::cout << x << ": "
51 #define LOGINFO(x) std::cout << x << ": "
52 #define ENDL std::endl
65 const int code_version = {21};
72 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
74 int nzeros = 4 - tempfile.length();
79 tempfile =
dir +
"template_summary2D_zp" + tempfile +
".out";
81 tempfile =
file.fullPath();
85 std::ostringstream tout;
86 tout <<
"template_summary2D_zp" << std::setw(4) << std::setfill(
'0') << std::right << filenum <<
".out" << std::ends;
87 tempfile = tout.str();
93 std::ifstream
in_file(tempfile);
100 for (
int i = 0; (
c =
in_file.get()) !=
'\n'; ++
i) {
107 LOGINFO(
"SiPixelTemplate2D") <<
"Loading Pixel Template File - " << theCurrentTemp.
head.
title <<
ENDL;
118 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 0A, no template load" <<
ENDL;
127 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 0B, no template load" <<
ENDL;
136 theCurrentTemp.
head.
fbin[1] = 1.00f;
137 theCurrentTemp.
head.
fbin[2] = 0.85f;
140 LOGINFO(
"SiPixelTemplate2D") <<
"Template ID = " << theCurrentTemp.
head.
ID <<
", Template Version "
142 <<
", NTy = " << theCurrentTemp.
head.
NTy <<
", NTyx = " << theCurrentTemp.
head.
NTyx
143 <<
", NTxx = " << theCurrentTemp.
head.
NTxx <<
", Dtype = " << theCurrentTemp.
head.
Dtype
144 <<
", Bias voltage " << theCurrentTemp.
head.
Vbias <<
", temperature "
146 <<
", Q-scaling factor " << theCurrentTemp.
head.
qscale <<
", 1/2 multi dcol threshold "
147 << theCurrentTemp.
head.
s50 <<
", 1/2 single dcol threshold "
149 <<
", y Lorentz Bias " << theCurrentTemp.
head.
lorybias <<
", x Lorentz width "
151 <<
", Q/Q_avg fractions for Qbin defs " << theCurrentTemp.
head.
fbin[0] <<
", "
153 <<
", pixel x-size " << theCurrentTemp.
head.
xsize <<
", y-size "
157 LOGERROR(
"SiPixelTemplate2D") <<
"code expects version " << code_version <<
", no template load" <<
ENDL;
161 if (theCurrentTemp.
head.
NTy != 0) {
163 <<
"Trying to load 1-d template info into the 2-d template object, check your DB/global tag!" <<
ENDL;
175 for (
int iy = 0; iy < theCurrentTemp.
head.
NTyx; ++iy) {
176 for (
int jx = 0; jx < theCurrentTemp.
head.
NTxx; ++jx) {
177 in_file >> theCurrentTemp.
entry[iy][jx].runnum >> theCurrentTemp.
entry[iy][jx].costrk[0] >>
178 theCurrentTemp.
entry[iy][jx].costrk[1] >> theCurrentTemp.
entry[iy][jx].costrk[2];
181 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 1, no template load, run # "
182 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
188 theCurrentTemp.
entry[iy][jx].cotalpha =
189 theCurrentTemp.
entry[iy][jx].costrk[0] / theCurrentTemp.
entry[iy][jx].costrk[2];
191 theCurrentTemp.
entry[iy][jx].cotbeta =
192 theCurrentTemp.
entry[iy][jx].costrk[1] / theCurrentTemp.
entry[iy][jx].costrk[2];
194 in_file >> theCurrentTemp.
entry[iy][jx].qavg >> theCurrentTemp.
entry[iy][jx].pixmax >>
195 theCurrentTemp.
entry[iy][jx].sxymax >> theCurrentTemp.
entry[iy][jx].iymin >>
196 theCurrentTemp.
entry[iy][jx].iymax >> theCurrentTemp.
entry[iy][jx].jxmin >>
197 theCurrentTemp.
entry[iy][jx].jxmax;
200 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 2, no template load, run # "
201 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
205 for (
int k = 0;
k < 2; ++
k) {
206 in_file >> theCurrentTemp.
entry[iy][jx].xypar[
k][0] >> theCurrentTemp.
entry[iy][jx].xypar[
k][1] >>
207 theCurrentTemp.
entry[iy][jx].xypar[
k][2] >> theCurrentTemp.
entry[iy][jx].xypar[
k][3] >>
208 theCurrentTemp.
entry[iy][jx].xypar[
k][4];
212 <<
"Error reading file 3, no template load, run # " << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
217 for (
int k = 0;
k < 2; ++
k) {
218 in_file >> theCurrentTemp.
entry[iy][jx].lanpar[
k][0] >> theCurrentTemp.
entry[iy][jx].lanpar[
k][1] >>
219 theCurrentTemp.
entry[iy][jx].lanpar[
k][2] >> theCurrentTemp.
entry[iy][jx].lanpar[
k][3] >>
220 theCurrentTemp.
entry[iy][jx].lanpar[
k][4];
224 <<
"Error reading file 4, no template load, run # " << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
232 for (
int l = 0;
l < 7; ++
l) {
233 for (
int k = 0;
k < 7; ++
k) {
240 <<
"Error reading file 5, no template load, run # " << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
250 in_file >> theCurrentTemp.
entry[iy][jx].chi2ppix >> theCurrentTemp.
entry[iy][jx].chi2scale >>
251 theCurrentTemp.
entry[iy][jx].offsetx[0] >> theCurrentTemp.
entry[iy][jx].offsetx[1] >>
252 theCurrentTemp.
entry[iy][jx].offsetx[2] >> theCurrentTemp.
entry[iy][jx].offsetx[3] >>
253 theCurrentTemp.
entry[iy][jx].offsety[0] >> theCurrentTemp.
entry[iy][jx].offsety[1] >>
254 theCurrentTemp.
entry[iy][jx].offsety[2] >> theCurrentTemp.
entry[iy][jx].offsety[3];
257 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 6, no template load, run # "
258 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
262 in_file >> theCurrentTemp.
entry[iy][jx].clsleny >> theCurrentTemp.
entry[iy][jx].clslenx >>
263 theCurrentTemp.
entry[iy][jx].mpvvav >> theCurrentTemp.
entry[iy][jx].sigmavav >>
264 theCurrentTemp.
entry[iy][jx].kappavav >> theCurrentTemp.
entry[iy][jx].scalexavg >>
265 theCurrentTemp.
entry[iy][jx].scaleyavg >> theCurrentTemp.
entry[iy][jx].delyavg >>
266 theCurrentTemp.
entry[iy][jx].delysig >> theCurrentTemp.
entry[iy][jx].spare[0];
269 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 7, no template load, run # "
270 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
274 in_file >> theCurrentTemp.
entry[iy][jx].scalex[0] >> theCurrentTemp.
entry[iy][jx].scalex[1] >>
275 theCurrentTemp.
entry[iy][jx].scalex[2] >> theCurrentTemp.
entry[iy][jx].scalex[3] >>
276 theCurrentTemp.
entry[iy][jx].scaley[0] >> theCurrentTemp.
entry[iy][jx].scaley[1] >>
277 theCurrentTemp.
entry[iy][jx].scaley[2] >> theCurrentTemp.
entry[iy][jx].scaley[3] >>
278 theCurrentTemp.
entry[iy][jx].spare[1] >> theCurrentTemp.
entry[iy][jx].spare[2];
281 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 8, no template load, run # "
282 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
292 pixelTemp.push_back(theCurrentTemp);
299 LOGERROR(
"SiPixelTemplate2D") <<
"Error opening File" << tempfile <<
ENDL;
305 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
313 std::vector<SiPixelTemplateStore2D>& pixelTemp) {
316 const int code_version = {21};
325 for (
int m = 0;
m <
db.numOfTempl(); ++
m) {
329 for (
int i = 0;
i < 20; ++
i) {
335 db.incrementIndex(1);
338 LOGINFO(
"SiPixelTemplate2D") <<
"Loading Pixel Template File - " << theCurrentTemp.
head.
title <<
ENDL;
350 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 0A, no template load" <<
ENDL;
354 LOGINFO(
"SiPixelTemplate2D") <<
"Loading Pixel Template File - " << theCurrentTemp.
head.
title
355 <<
" code version = " << code_version <<
" object version "
363 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 0B, no template load" <<
ENDL;
371 theCurrentTemp.
head.
fbin[0] = 1.50f;
372 theCurrentTemp.
head.
fbin[1] = 1.00f;
373 theCurrentTemp.
head.
fbin[2] = 0.85f;
376 LOGINFO(
"SiPixelTemplate2D") <<
"Template ID = " << theCurrentTemp.
head.
ID <<
", Template Version "
378 <<
", NTy = " << theCurrentTemp.
head.
NTy <<
", NTyx = " << theCurrentTemp.
head.
NTyx
379 <<
", NTxx = " << theCurrentTemp.
head.
NTxx <<
", Dtype = " << theCurrentTemp.
head.
Dtype
380 <<
", Bias voltage " << theCurrentTemp.
head.
Vbias <<
", temperature "
382 <<
", Q-scaling factor " << theCurrentTemp.
head.
qscale <<
", 1/2 multi dcol threshold "
383 << theCurrentTemp.
head.
s50 <<
", 1/2 single dcol threshold "
385 <<
", y Lorentz Bias " << theCurrentTemp.
head.
lorybias <<
", x Lorentz width "
387 <<
", Q/Q_avg fractions for Qbin defs " << theCurrentTemp.
head.
fbin[0] <<
", "
389 <<
", pixel x-size " << theCurrentTemp.
head.
xsize <<
", y-size "
393 LOGINFO(
"SiPixelTemplate2D") <<
"code expects version " << code_version <<
" finds "
397 if (theCurrentTemp.
head.
NTy != 0) {
399 <<
"Trying to load 1-d template info into the 2-d template object, check your DB/global tag!" <<
ENDL;
411 for (
int iy = 0; iy < theCurrentTemp.
head.
NTyx; ++iy) {
412 for (
int jx = 0; jx < theCurrentTemp.
head.
NTxx; ++jx) {
413 db >> theCurrentTemp.
entry[iy][jx].runnum >> theCurrentTemp.
entry[iy][jx].costrk[0] >>
414 theCurrentTemp.
entry[iy][jx].costrk[1] >> theCurrentTemp.
entry[iy][jx].costrk[2];
417 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 1, no template load, run # "
418 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
424 theCurrentTemp.
entry[iy][jx].cotalpha =
425 theCurrentTemp.
entry[iy][jx].costrk[0] / theCurrentTemp.
entry[iy][jx].costrk[2];
427 theCurrentTemp.
entry[iy][jx].cotbeta =
428 theCurrentTemp.
entry[iy][jx].costrk[1] / theCurrentTemp.
entry[iy][jx].costrk[2];
430 db >> theCurrentTemp.
entry[iy][jx].qavg >> theCurrentTemp.
entry[iy][jx].pixmax >>
431 theCurrentTemp.
entry[iy][jx].sxymax >> theCurrentTemp.
entry[iy][jx].iymin >>
432 theCurrentTemp.
entry[iy][jx].iymax >> theCurrentTemp.
entry[iy][jx].jxmin >>
433 theCurrentTemp.
entry[iy][jx].jxmax;
436 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 2, no template load, run # "
437 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
441 for (
int k = 0;
k < 2; ++
k) {
442 db >> theCurrentTemp.
entry[iy][jx].xypar[
k][0] >> theCurrentTemp.
entry[iy][jx].xypar[
k][1] >>
443 theCurrentTemp.
entry[iy][jx].xypar[
k][2] >> theCurrentTemp.
entry[iy][jx].xypar[
k][3] >>
444 theCurrentTemp.
entry[iy][jx].xypar[
k][4];
448 <<
"Error reading file 3, no template load, run # " << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
453 for (
int k = 0;
k < 2; ++
k) {
454 db >> theCurrentTemp.
entry[iy][jx].lanpar[
k][0] >> theCurrentTemp.
entry[iy][jx].lanpar[
k][1] >>
455 theCurrentTemp.
entry[iy][jx].lanpar[
k][2] >> theCurrentTemp.
entry[iy][jx].lanpar[
k][3] >>
456 theCurrentTemp.
entry[iy][jx].lanpar[
k][4];
460 <<
"Error reading file 4, no template load, run # " << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
468 for (
int l = 0;
l < 7; ++
l) {
469 for (
int k = 0;
k < 7; ++
k) {
476 <<
"Error reading file 5, no template load, run # " << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
486 db >> theCurrentTemp.
entry[iy][jx].chi2ppix >> theCurrentTemp.
entry[iy][jx].chi2scale >>
487 theCurrentTemp.
entry[iy][jx].offsetx[0] >> theCurrentTemp.
entry[iy][jx].offsetx[1] >>
488 theCurrentTemp.
entry[iy][jx].offsetx[2] >> theCurrentTemp.
entry[iy][jx].offsetx[3] >>
489 theCurrentTemp.
entry[iy][jx].offsety[0] >> theCurrentTemp.
entry[iy][jx].offsety[1] >>
490 theCurrentTemp.
entry[iy][jx].offsety[2] >> theCurrentTemp.
entry[iy][jx].offsety[3];
493 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 6, no template load, run # "
494 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
498 db >> theCurrentTemp.
entry[iy][jx].clsleny >> theCurrentTemp.
entry[iy][jx].clslenx >>
499 theCurrentTemp.
entry[iy][jx].mpvvav >> theCurrentTemp.
entry[iy][jx].sigmavav >>
500 theCurrentTemp.
entry[iy][jx].kappavav >> theCurrentTemp.
entry[iy][jx].scalexavg >>
501 theCurrentTemp.
entry[iy][jx].scaleyavg >> theCurrentTemp.
entry[iy][jx].delyavg >>
502 theCurrentTemp.
entry[iy][jx].delysig >> theCurrentTemp.
entry[iy][jx].spare[0];
505 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 7, no template load, run # "
506 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
510 db >> theCurrentTemp.
entry[iy][jx].scalex[0] >> theCurrentTemp.
entry[iy][jx].scalex[1] >>
511 theCurrentTemp.
entry[iy][jx].scalex[2] >> theCurrentTemp.
entry[iy][jx].scalex[3] >>
512 theCurrentTemp.
entry[iy][jx].scaley[0] >> theCurrentTemp.
entry[iy][jx].scaley[1] >>
513 theCurrentTemp.
entry[iy][jx].scaley[2] >> theCurrentTemp.
entry[iy][jx].scaley[3] >>
514 theCurrentTemp.
entry[iy][jx].spare[1] >> theCurrentTemp.
entry[iy][jx].spare[2];
517 LOGERROR(
"SiPixelTemplate2D") <<
"Error reading file 8, no template load, run # "
518 << theCurrentTemp.
entry[iy][jx].runnum <<
ENDL;
526 pixelTemp.push_back(theCurrentTemp);
536 if (
id != id_current_) {
540 for (
int i = 0;
i < (
int)thePixelTemp_.size(); ++
i) {
541 if (
id == thePixelTemp_[
i].head.ID) {
547 Dtype_ = thePixelTemp_[index_id_].head.Dtype;
551 qscale_ = thePixelTemp_[index_id_].head.qscale;
555 s50_ = thePixelTemp_[index_id_].head.s50;
559 for (
int j = 0;
j < 3; ++
j) {
560 fbin_[
j] = thePixelTemp_[index_id_].head.fbin[
j];
565 lorywidth_ = thePixelTemp_[index_id_].head.lorywidth;
566 lorxwidth_ = thePixelTemp_[index_id_].head.lorxwidth;
570 xsize_ = thePixelTemp_[index_id_].head.xsize;
571 ysize_ = thePixelTemp_[index_id_].head.ysize;
572 zsize_ = thePixelTemp_[index_id_].head.zsize;
576 Nyx_ = thePixelTemp_[index_id_].head.NTyx;
577 Nxx_ = thePixelTemp_[index_id_].head.NTxx;
578 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
579 if (Nyx_ < 2 || Nxx_ < 2) {
580 throw cms::Exception(
"DataCorrupt") <<
"template ID = " << id_current_
581 <<
"has too few entries: Nyx/Nxx = " << Nyx_ <<
"/" << Nxx_ << std::endl;
584 assert(Nyx_ > 1 && Nxx_ > 1);
586 int imidx = Nxx_ / 2;
588 cotalpha0_ = thePixelTemp_[index_id_].entry[0][0].cotalpha;
589 cotalpha1_ = thePixelTemp_[index_id_].entry[0][Nxx_ - 1].cotalpha;
590 deltacota_ = (cotalpha1_ - cotalpha0_) / (
float)(Nxx_ - 1);
592 cotbeta0_ = thePixelTemp_[index_id_].entry[0][imidx].cotbeta;
593 cotbeta1_ = thePixelTemp_[index_id_].entry[Nyx_ - 1][imidx].cotbeta;
594 deltacotb_ = (cotbeta1_ - cotbeta0_) / (
float)(Nyx_ - 1);
601 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
602 if (index_id_ < 0 || index_id_ >= (
int)thePixelTemp_.size()) {
603 throw cms::Exception(
"DataCorrupt") <<
"SiPixelTemplate2D::interpolate can't find needed template ID = " <<
id
604 <<
", Are you using the correct global tag?" << std::endl;
607 assert(index_id_ >= 0 && index_id_ < (
int)thePixelTemp_.size());
630 float acotb, dcota, dcotb;
634 if (
id != id_current_ || cotalpha != cota_current_ || cotbeta != cotb_current_) {
635 cota_current_ = cotalpha;
636 cotb_current_ = cotbeta;
638 success_ = getid(
id);
641 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
642 if (index_id_ < 0 || index_id_ >= (
int)thePixelTemp_.size()) {
643 throw cms::Exception(
"DataCorrupt") <<
"SiPixelTemplate2D::interpolate can't find needed template ID = " <<
id
644 <<
", Are you using the correct global tag?" << std::endl;
647 assert(index_id_ >= 0 && index_id_ < (
int)thePixelTemp_.size());
652 float cota = cotalpha;
670 if (locBx * locBz < 0.
f) {
679 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
681 <<
"SiPixelTemplate2D::illegal subdetector ID = " << thePixelTemp_[index_id_].head.Dtype << std::endl;
683 std::cout <<
"SiPixelTemplate:2D:illegal subdetector ID = " << thePixelTemp_[index_id_].head.Dtype << std::endl;
687 if (cota < cotalpha0_) {
692 }
else if (cota > cotalpha1_) {
698 jx0_ = (
int)((cota - cotalpha0_) / deltacota_ + 0.5f);
699 dcota = (cota - (cotalpha0_ + jx0_ * deltacota_)) / deltacota_;
700 adcota_ = fabs(dcota);
716 if (acotb < cotbeta0_) {
721 }
else if (acotb > cotbeta1_) {
727 iy0_ = (
int)((acotb - cotbeta0_) / deltacotb_ + 0.5f);
728 dcotb = (acotb - (cotbeta0_ + iy0_ * deltacotb_)) / deltacotb_;
729 adcotb_ = fabs(dcotb);
743 lorydrift_ = lorywidth_ / 2.;
745 lorydrift_ = -lorydrift_;
746 lorxdrift_ = lorxwidth_ / 2.;
748 lorxdrift_ = -lorxdrift_;
752 entry00_ = &thePixelTemp_[index_id_].entry[iy0_][jx0_];
753 entry10_ = &thePixelTemp_[index_id_].entry[iy1_][jx0_];
754 entry01_ = &thePixelTemp_[index_id_].entry[iy0_][jx1_];
758 qavg_ = entry00_->qavg + adcota_ * (entry01_->qavg - entry00_->qavg) + adcotb_ * (entry10_->qavg - entry00_->qavg);
760 pixmax_ = entry00_->pixmax + adcota_ * (entry01_->pixmax - entry00_->pixmax) +
761 adcotb_ * (entry10_->pixmax - entry00_->pixmax);
763 sxymax_ = entry00_->sxymax + adcota_ * (entry01_->sxymax - entry00_->sxymax) +
764 adcotb_ * (entry10_->sxymax - entry00_->sxymax);
766 chi2avgone_ = entry00_->chi2avgone + adcota_ * (entry01_->chi2avgone - entry00_->chi2avgone) +
767 adcotb_ * (entry10_->chi2avgone - entry00_->chi2avgone);
769 chi2minone_ = entry00_->chi2minone + adcota_ * (entry01_->chi2minone - entry00_->chi2minone) +
770 adcotb_ * (entry10_->chi2minone - entry00_->chi2minone);
772 clsleny_ = entry00_->clsleny + adcota_ * (entry01_->clsleny - entry00_->clsleny) +
773 adcotb_ * (entry10_->clsleny - entry00_->clsleny);
775 clslenx_ = entry00_->clslenx + adcota_ * (entry01_->clslenx - entry00_->clslenx) +
776 adcotb_ * (entry10_->clslenx - entry00_->clslenx);
778 chi2ppix_ = entry00_->chi2ppix + adcota_ * (entry01_->chi2ppix - entry00_->chi2ppix) +
779 adcotb_ * (entry10_->chi2ppix - entry00_->chi2ppix);
781 chi2scale_ = entry00_->chi2scale + adcota_ * (entry01_->chi2scale - entry00_->chi2scale) +
782 adcotb_ * (entry10_->chi2scale - entry00_->chi2scale);
784 scaleyavg_ = entry00_->scaleyavg + adcota_ * (entry01_->scaleyavg - entry00_->scaleyavg) +
785 adcotb_ * (entry10_->scaleyavg - entry00_->scaleyavg);
787 scalexavg_ = entry00_->scalexavg + adcota_ * (entry01_->scalexavg - entry00_->scalexavg) +
788 adcotb_ * (entry10_->scalexavg - entry00_->scalexavg);
790 delyavg_ = entry00_->delyavg + adcota_ * (entry01_->delyavg - entry00_->delyavg) +
791 adcotb_ * (entry10_->delyavg - entry00_->delyavg);
793 delysig_ = entry00_->delysig + adcota_ * (entry01_->delysig - entry00_->delysig) +
794 adcotb_ * (entry10_->delysig - entry00_->delysig);
796 mpvvav_ = entry00_->mpvvav + adcota_ * (entry01_->mpvvav - entry00_->mpvvav) +
797 adcotb_ * (entry10_->mpvvav - entry00_->mpvvav);
799 sigmavav_ = entry00_->sigmavav + adcota_ * (entry01_->sigmavav - entry00_->sigmavav) +
800 adcotb_ * (entry10_->sigmavav - entry00_->sigmavav);
802 kappavav_ = entry00_->kappavav + adcota_ * (entry01_->kappavav - entry00_->kappavav) +
803 adcotb_ * (entry10_->kappavav - entry00_->kappavav);
805 for (
int i = 0;
i < 4; ++
i) {
806 scalex_[
i] = entry00_->scalex[
i] + adcota_ * (entry01_->scalex[
i] - entry00_->scalex[
i]) +
807 adcotb_ * (entry10_->scalex[
i] - entry00_->scalex[
i]);
809 scaley_[
i] = entry00_->scaley[
i] + adcota_ * (entry01_->scaley[
i] - entry00_->scaley[
i]) +
810 adcotb_ * (entry10_->scaley[
i] - entry00_->scaley[
i]);
812 offsetx_[
i] = entry00_->offsetx[
i] + adcota_ * (entry01_->offsetx[
i] - entry00_->offsetx[
i]) +
813 adcotb_ * (entry10_->offsetx[
i] - entry00_->offsetx[
i]);
815 offsetx_[
i] = -offsetx_[
i];
817 offsety_[
i] = entry00_->offsety[
i] + adcota_ * (entry01_->offsety[
i] - entry00_->offsety[
i]) +
818 adcotb_ * (entry10_->offsety[
i] - entry00_->offsety[
i]);
820 offsety_[
i] = -offsety_[
i];
823 for (
int i = 0;
i < 2; ++
i) {
824 for (
int j = 0;
j < 5; ++
j) {
827 xypary0x0_[1 -
i][
j] = (
float)entry00_->xypar[
i][
j];
828 xypary1x0_[1 -
i][
j] = (
float)entry10_->xypar[
i][
j];
829 xypary0x1_[1 -
i][
j] = (
float)entry01_->xypar[
i][
j];
830 lanpar_[1 -
i][
j] = entry00_->lanpar[
i][
j] + adcota_ * (entry01_->lanpar[
i][
j] - entry00_->lanpar[
i][
j]) +
831 adcotb_ * (entry10_->lanpar[
i][
j] - entry00_->lanpar[
i][
j]);
833 xypary0x0_[
i][
j] = (
float)entry00_->xypar[
i][
j];
834 xypary1x0_[
i][
j] = (
float)entry10_->xypar[
i][
j];
835 xypary0x1_[
i][
j] = (
float)entry01_->xypar[
i][
j];
836 lanpar_[
i][
j] = entry00_->lanpar[
i][
j] + adcota_ * (entry01_->lanpar[
i][
j] - entry00_->lanpar[
i][
j]) +
837 adcotb_ * (entry10_->lanpar[
i][
j] - entry00_->lanpar[
i][
j]);
853 #ifdef SI_PIXEL_TEMPLATE_STANDALONE
878 for (
int i = 0;
i < 3; ++
i) {
889 qavg_ = entry00_->qavg;
891 pixmax_ = entry00_->pixmax;
893 sxymax_ = entry00_->sxymax;
895 clsleny_ = entry00_->clsleny;
897 clslenx_ = entry00_->clslenx;
907 for (
int i = 0;
i < 4; ++
i) {
918 float cotbeta = entry00_->cotbeta;
934 if (locBx * locBz < 0.
f) {
942 std::cout <<
"SiPixelTemplate:2D:illegal subdetector ID = " << iDtype << std::endl;
947 lorydrift_ = lorywidth_ / 2.;
949 lorydrift_ = -lorydrift_;
950 lorxdrift_ = lorxwidth_ / 2.;
952 lorxdrift_ = -lorxdrift_;
954 for (
int i = 0;
i < 2; ++
i) {
955 for (
int j = 0;
j < 5; ++
j) {
958 xypary0x0_[1 -
i][
j] = (
float)entry00_->xypar[
i][
j];
959 xypary1x0_[1 -
i][
j] = (
float)entry00_->xypar[
i][
j];
960 xypary0x1_[1 -
i][
j] = (
float)entry00_->xypar[
i][
j];
961 lanpar_[1 -
i][
j] = entry00_->lanpar[
i][
j];
963 xypary0x0_[
i][
j] = (
float)entry00_->xypar[
i][
j];
964 xypary1x0_[
i][
j] = (
float)entry00_->xypar[
i][
j];
965 xypary0x1_[
i][
j] = (
float)entry00_->xypar[
i][
j];
966 lanpar_[
i][
j] = entry00_->lanpar[
i][
j];
993 int pixx, pixy,
k0, k1, l0, l1, deltax, deltay, iflipy, jflipx, imin, imax, jmin, jmax;
995 float dx,
dy, ddx, ddy, adx, ady;
997 const float deltaxy[2] = {16.67f, 25.0f};
1006 pixy = (
int)floorf(yhit / ysize_);
1007 dy = yhit - (pixy + 0.5f) * ysize_;
1016 ddy = 6.f *
dy / ysize_ - (
k0 - 3);
1027 pixx = (
int)floorf(xhit / xsize_);
1028 dx = xhit - (pixx + 0.5f) * xsize_;
1032 l0 = (
int)(
dx / xsize_ * 6.
f + 3.5
f);
1037 ddx = 6.f *
dx / xsize_ - (l0 - 3);
1053 imin =
std::min(entry00_->iymin, entry10_->iymin);
1054 imin_ =
std::min(imin, entry01_->iymin);
1056 jmin =
std::min(entry00_->jxmin, entry10_->jxmin);
1057 jmin_ =
std::min(jmin, entry01_->jxmin);
1059 imax =
std::max(entry00_->iymax, entry10_->iymax);
1060 imax_ =
std::max(imax, entry01_->iymax);
1062 jmax =
std::max(entry00_->jxmax, entry10_->jxmax);
1063 jmax_ =
std::max(jmax, entry01_->jxmax);
1074 deltax = pixx -
T2HX;
1075 deltay = pixy -
T2HY;
1079 for (
int j = 0;
j <
BXM2; ++
j) {
1080 for (
int i = 0;
i <
BYM2; ++
i) {
1081 xytemp_[
j][
i] = 0.f;
1087 for (
int j = jmin_;
j <= jmax_; ++
j) {
1091 m = deltax + jflipx;
1095 for (
int i = imin_;
i <= imax_; ++
i) {
1098 n = deltay + iflipy;
1102 if (
m >= 0 && m <= BXM3 && n >= 0 &&
n <=
BYM3) {
1103 xytemp_[
m][
n] = (
float)entry00_->xytemp[
k0][l0][
i][
j] +
1104 adx * (
float)(entry00_->xytemp[
k0][l1][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1105 ady * (
float)(entry00_->xytemp[k1][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1106 adcota_ * (
float)(entry01_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1107 adcotb_ * (
float)(entry10_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]);
1114 for (
int n = 1;
n <
BYM3; ++
n) {
1117 for (
int m = 1;
m <
BXM3; ++
m) {
1118 xytemp_[
m][
n] += xytemp_[
m][
n + 1];
1121 for (
int i =
n + 1;
i <
BYM3; ++
i) {
1122 for (
int m = 1;
m <
BXM3; ++
m) {
1123 xytemp_[
m][
i] = xytemp_[
m][
i + 1];
1131 for (
int m = 1;
m <
BXM3; ++
m) {
1134 for (
int n = 1;
n <
BYM3; ++
n) {
1135 xytemp_[
m][
n] += xytemp_[
m + 1][
n];
1138 for (
int j =
m + 1;
j <
BXM3; ++
j) {
1140 xytemp_[
j][
n] = xytemp_[
j + 1][
n];
1148 float qtemptot = 0.f;
1150 for (
int n = 1;
n <
BYM3; ++
n) {
1151 for (
int m = 1;
m <
BXM3; ++
m) {
1152 if (xytemp_[
m][
n] != 0.
f) {
1153 template2d[
m][
n] += xytemp_[
m][
n];
1154 qtemptot += xytemp_[
m][
n];
1159 QTemplate = qtemptot;
1164 for (
int k = 0;
k < 2; ++
k) {
1165 for (
int i = 0;
i <
BXM2; ++
i) {
1166 for (
int j = 0;
j <
BYM2; ++
j) {
1167 dxytempdx[
k][
i][
j] = 0.f;
1168 dxytempdy[
k][
i][
j] = 0.f;
1169 dpdx2d[
k][
i][
j] = 0.f;
1176 pixx = (
int)floorf((xhit + deltaxy[0]) / xsize_);
1177 dx = (xhit + deltaxy[0]) - (pixx + 0.5
f) * xsize_;
1181 l0 = (
int)(
dx / xsize_ * 6.
f + 3.5
f);
1186 ddx = 6.f *
dx / xsize_ - (l0 - 3);
1208 deltax = pixx -
T2HX;
1212 for (
int j = jmin_;
j <= jmax_; ++
j) {
1216 m = deltax + jflipx;
1220 for (
int i = imin_;
i <= imax_; ++
i) {
1223 n = deltay + iflipy;
1227 if (
m >= 0 && m <= BXM3 && n >= 0 &&
n <=
BYM3) {
1228 dxytempdx[1][
m][
n] = (
float)entry00_->xytemp[
k0][l0][
i][
j] +
1229 adx * (
float)(entry00_->xytemp[
k0][l1][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1230 ady * (
float)(entry00_->xytemp[k1][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1231 adcota_ * (
float)(entry01_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1232 adcotb_ * (
float)(entry10_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]);
1239 for (
int n = 1;
n <
BYM3; ++
n) {
1242 for (
int m = 1;
m <
BXM3; ++
m) {
1243 dxytempdx[1][
m][
n] += dxytempdx[1][
m][
n + 1];
1246 for (
int i =
n + 1;
i <
BYM3; ++
i) {
1247 for (
int m = 1;
m <
BXM3; ++
m) {
1248 dxytempdx[1][
m][
i] = dxytempdx[1][
m][
i + 1];
1256 for (
int m = 1;
m <
BXM3; ++
m) {
1259 for (
int n = 1;
n <
BYM3; ++
n) {
1260 dxytempdx[1][
m][
n] += dxytempdx[1][
m + 1][
n];
1263 for (
int j =
m + 1;
j <
BXM3; ++
j) {
1264 for (
int n = 1;
n <
BYM3; ++
n) {
1265 dxytempdx[1][
j][
n] = dxytempdx[1][
j + 1][
n];
1273 pixx = (
int)floorf((xhit - deltaxy[0]) / xsize_);
1274 dx = (xhit - deltaxy[0]) - (pixx + 0.5
f) * xsize_;
1278 l0 = (
int)(
dx / xsize_ * 6.
f + 3.5
f);
1283 ddx = 6.f *
dx / xsize_ - (l0 - 3);
1305 deltax = pixx -
T2HX;
1309 for (
int j = jmin_;
j <= jmax_; ++
j) {
1313 m = deltax + jflipx;
1317 for (
int i = imin_;
i <= imax_; ++
i) {
1320 n = deltay + iflipy;
1324 if (
m >= 0 && m <= BXM3 && n >= 0 &&
n <=
BYM3) {
1325 dxytempdx[0][
m][
n] = (
float)entry00_->xytemp[
k0][l0][
i][
j] +
1326 adx * (
float)(entry00_->xytemp[
k0][l1][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1327 ady * (
float)(entry00_->xytemp[k1][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1328 adcota_ * (
float)(entry01_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1329 adcotb_ * (
float)(entry10_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]);
1336 for (
int n = 1;
n <
BYM3; ++
n) {
1339 for (
int m = 1;
m <
BXM3; ++
m) {
1340 dxytempdx[0][
m][
n] += dxytempdx[0][
m][
n + 1];
1343 for (
int i =
n + 1;
i <
BYM3; ++
i) {
1344 for (
int m = 1;
m <
BXM3; ++
m) {
1345 dxytempdx[0][
m][
i] = dxytempdx[0][
m][
i + 1];
1353 for (
int m = 1;
m <
BXM3; ++
m) {
1356 for (
int n = 1;
n <
BYM3; ++
n) {
1357 dxytempdx[0][
m][
n] += dxytempdx[0][
m + 1][
n];
1360 for (
int j =
m + 1;
j <
BXM3; ++
j) {
1361 for (
int n = 1;
n <
BYM3; ++
n) {
1362 dxytempdx[0][
j][
n] = dxytempdx[0][
j + 1][
n];
1370 for (
int n = 1;
n <
BYM3; ++
n) {
1371 for (
int m = 1;
m <
BXM3; ++
m) {
1372 dpdx2d[0][
m][
n] = (dxytempdx[1][
m][
n] - dxytempdx[0][
m][
n]) / (2. * deltaxy[0]);
1378 pixy = (
int)floorf((yhit + deltaxy[1]) / ysize_);
1379 dy = (yhit + deltaxy[1]) - (pixy + 0.5
f) * ysize_;
1388 ddy = 6.f *
dy / ysize_ - (
k0 - 3);
1399 pixx = (
int)floorf(xhit / xsize_);
1400 dx = xhit - (pixx + 0.5f) * xsize_;
1404 l0 = (
int)(
dx / xsize_ * 6.
f + 3.5
f);
1409 ddx = 6.f *
dx / xsize_ - (l0 - 3);
1432 deltax = pixx -
T2HX;
1433 deltay = pixy -
T2HY;
1437 for (
int j = jmin_;
j <= jmax_; ++
j) {
1441 m = deltax + jflipx;
1445 for (
int i = imin_;
i <= imax_; ++
i) {
1448 n = deltay + iflipy;
1452 if (
m >= 0 && m <= BXM3 && n >= 0 &&
n <=
BYM3) {
1453 dxytempdy[1][
m][
n] = (
float)entry00_->xytemp[
k0][l0][
i][
j] +
1454 adx * (
float)(entry00_->xytemp[
k0][l1][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1455 ady * (
float)(entry00_->xytemp[k1][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1456 adcota_ * (
float)(entry01_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1457 adcotb_ * (
float)(entry10_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]);
1464 for (
int n = 1;
n <
BYM3; ++
n) {
1467 for (
int m = 1;
m <
BXM3; ++
m) {
1468 dxytempdy[1][
m][
n] += dxytempdy[1][
m][
n + 1];
1471 for (
int i =
n + 1;
i <
BYM3; ++
i) {
1472 for (
int m = 1;
m <
BXM3; ++
m) {
1473 dxytempdy[1][
m][
i] = dxytempdy[1][
m][
i + 1];
1481 for (
int m = 1;
m <
BXM3; ++
m) {
1484 for (
int n = 1;
n <
BYM3; ++
n) {
1485 dxytempdy[1][
m][
n] += dxytempdy[1][
m + 1][
n];
1488 for (
int j =
m + 1;
j <
BXM3; ++
j) {
1489 for (
int n = 1;
n <
BYM3; ++
n) {
1490 dxytempdy[1][
j][
n] = dxytempdy[1][
j + 1][
n];
1498 pixy = (
int)floorf((yhit - deltaxy[1]) / ysize_);
1499 dy = (yhit - deltaxy[1]) - (pixy + 0.5
f) * ysize_;
1508 ddy = 6.f *
dy / ysize_ - (
k0 - 3);
1530 deltay = pixy -
T2HY;
1534 for (
int j = jmin_;
j <= jmax_; ++
j) {
1538 m = deltax + jflipx;
1542 for (
int i = imin_;
i <= imax_; ++
i) {
1545 n = deltay + iflipy;
1549 if (
m >= 0 && m <= BXM3 && n >= 0 &&
n <=
BYM3) {
1550 dxytempdy[0][
m][
n] = (
float)entry00_->xytemp[
k0][l0][
i][
j] +
1551 adx * (
float)(entry00_->xytemp[
k0][l1][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1552 ady * (
float)(entry00_->xytemp[k1][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1553 adcota_ * (
float)(entry01_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]) +
1554 adcotb_ * (
float)(entry10_->xytemp[
k0][l0][
i][
j] - entry00_->xytemp[
k0][l0][
i][
j]);
1561 for (
int n = 1;
n <
BYM3; ++
n) {
1564 for (
int m = 1;
m <
BXM3; ++
m) {
1565 dxytempdy[0][
m][
n] += dxytempdy[0][
m][
n + 1];
1568 for (
int i =
n + 1;
i <
BYM3; ++
i) {
1569 for (
int m = 1;
m <
BXM3; ++
m) {
1570 dxytempdy[0][
m][
i] = dxytempdy[0][
m][
i + 1];
1578 for (
int m = 1;
m <
BXM3; ++
m) {
1581 for (
int n = 1;
n <
BYM3; ++
n) {
1582 dxytempdy[0][
m][
n] += dxytempdy[0][
m + 1][
n];
1585 for (
int j =
m + 1;
j <
BXM3; ++
j) {
1586 for (
int n = 1;
n <
BYM3; ++
n) {
1587 dxytempdy[0][
j][
n] = dxytempdy[0][
j + 1][
n];
1595 for (
int n = 1;
n <
BYM3; ++
n) {
1596 for (
int m = 1;
m <
BXM3; ++
m) {
1597 dpdx2d[1][
m][
n] = (dxytempdy[1][
m][
n] - dxytempdy[0][
m][
n]) / (2. * deltaxy[1]);
1615 float xhit,
float yhit,
bool ydouble[
BYM2],
bool xdouble[
BXM2],
float template2d[
BXM2][
BYM2]) {
1618 bool derivatives =
false;
1643 std::vector<bool>& ydouble,
1644 std::vector<bool>& xdouble,
1648 bool derivatives =
false;
1652 if (cotbeta < 0.
f) {
1655 float locBz = locBx;
1656 if (cotalpha < 0.
f) {
1663 yd[
BYM2 - 1] =
false;
1665 yd[
i + 1] = ydouble[
i];
1668 xd[
BXM2 - 1] =
false;
1670 xd[
j + 1] = xdouble[
j];
1696 float sigi, sigi2, sigi3, sigi4, qscale, err2, err00;
1700 #ifndef SI_PIXEL_TEMPLATE_STANDALONE
1701 if (index < 1 || index >=
BYM2) {
1702 throw cms::Exception(
"DataCorrupt") <<
"SiPixelTemplate2D::ysigma2 called with index = " <<
index << std::endl;
1712 if (qpixel < sxymax_) {
1717 qscale = qpixel / sxymax_;
1719 sigi2 = sigi * sigi;
1720 sigi3 = sigi2 * sigi;
1721 sigi4 = sigi3 * sigi;
1723 err00 = xypary0x0_[0][0] + xypary0x0_[0][1] * sigi + xypary0x0_[0][2] * sigi2 + xypary0x0_[0][3] * sigi3 +
1724 xypary0x0_[0][4] * sigi4;
1726 adcota_ * (xypary0x1_[0][0] + xypary0x1_[0][1] * sigi + xypary0x1_[0][2] * sigi2 + xypary0x1_[0][3] * sigi3 +
1727 xypary0x1_[0][4] * sigi4 - err00) +
1728 adcotb_ * (xypary1x0_[0][0] + xypary1x0_[0][1] * sigi + xypary1x0_[0][2] * sigi2 + xypary1x0_[0][3] * sigi3 +
1729 xypary1x0_[0][4] * sigi4 - err00);
1731 err00 = xypary0x0_[1][0] + xypary0x0_[1][1] * sigi + xypary0x0_[1][2] * sigi2 + xypary0x0_[1][3] * sigi3 +
1732 xypary0x0_[1][4] * sigi4;
1734 adcota_ * (xypary0x1_[1][0] + xypary0x1_[1][1] * sigi + xypary0x1_[1][2] * sigi2 + xypary0x1_[1][3] * sigi3 +
1735 xypary0x1_[1][4] * sigi4 - err00) +
1736 adcotb_ * (xypary1x0_[1][0] + xypary1x0_[1][1] * sigi + xypary1x0_[1][2] * sigi2 + xypary1x0_[1][3] * sigi3 +
1737 xypary1x0_[1][4] * sigi4 - err00);
1739 xysig2 = qscale * err2;
1740 if (xysig2 <= 0.
f) {
1741 xysig2 = s50_ * s50_;
1756 for (
int i = 0;
i < 2; ++
i) {
1757 for (
int j = 0;
j < 5; ++
j) {
1758 lanpar[
i][
j] = lanpar_[
i][
j];