401 {
402
404
405
406
407
408 double meanCalib[4];
409 meanCalib[0] = 1.295;
410 meanCalib[1] = 81.454;
411
412 meanCalib[2] = 1.26282;
413 meanCalib[3] = 40.9303;
414
415
416 HWIdentifier hwid;
417
418
423
424
426
427
429 TGraphErrors* grrms;
432 int nevt, ndigerr = 0;
433 int badPts;
434 double maxRMS;
435 double maxPointInFitRange;
436
437
438 TF1 *fslope = new TF1("fslope", "[0]*x", 0, 1000);
439
442
443
446
447
448 for (; adcIter != adcIterE; ++adcIter) {
449 hwid = (adcIter)->first;
450 MeanDacMap = (adcIter)->second;
454
459
460 int fragId = (
ros << 8) | drawer;
462 int gain_offset = (demonstrator) ? 2 : 0;
463 int gain_ind =
gain + gain_offset;
466 double prevMean = 0.0, deltaMean = 0.0, averSlope = 0.0;
467 double prevCharge = 0.0, deltaCharge = 0.0, deltaChargeMin = 999999.;
468 double prevSlope = 0.0, slope = 0.0, maxGoodAmp=0.0, minGoodAmp = 999999.;
469 double maxGoodCharge =
m_linfitMax[gain_ind], minGoodCharge = 999999.;
470 double maxCharge = 0.0, minCharge = 999999., maxAmp = 0.0, minAmp = 999999.;
471 bool signalInRange = true;
472 int nptGood = 0;
473
474
475 npt = MeanDacMap->size();
477 gr =
new TGraphErrors(npt);
478 grrms = new TGraphErrors(npt);
479
480 if (npt == 0) {
484 << ros << "/" << drawer << "/" << chan << "/" << gain );
485 } else {
486
487
489
490
493
494
496 badPts = 0;
497 maxPointInFitRange = 0.0;
498 maxRMS = 0.0;
499 for (; dacIter != dacIterE; ++dacIter) {
500 dac = (dacIter)->first;
501 mean = (dacIter)->second;
502 meansq = (*MeanSqDacMap)[
dac];
503 nevt = (*NEvtDacMap)[
dac];
504 ndigerr = (*NDigitalErrorsDacMap)[
dac];
505
508 meansq = meansq / nevt;
509
510 rms = (meansq <= mean2) ? 0. : sqrt(meansq - mean2);
511 err = sqrt(rms * rms / nevt + 0.5 * 0.5);
512
513
515
521 }
522
523 deltaMean =
mean - prevMean;
524 deltaCharge =
charge - prevCharge;
525 if (deltaCharge != 0) {
526 slope = deltaMean / deltaCharge;
527 if (deltaCharge<deltaChargeMin) deltaChargeMin = deltaCharge;
528 }
529
530
533 signalInRange = false;
535 << ros << "/" << drawer
536 << "/" << chan << "/" << gain
539 << " => Ignoring this point ");
540 } else {
541 if (prevSlope != 0) {
542 double R = std::abs(slope/prevSlope);
543 double R1 = (nptGood>1) ? (slope / averSlope) : 1.;
544 if (R<0.025 || (nptGood>2 && R1>2.5 && R>2.5) || (nptGood==2 && R1>4.1) ) {
545
547 << ros << "/" << drawer
548 << "/" << chan << "/" << gain
551 << " dC " << deltaCharge
552 << " dA " << deltaMean
553 << " expected dA " << prevSlope*deltaCharge
554 << " aver dA " << averSlope*deltaCharge
555 << " rms " << rms
556 << " => Removing this point ");
557 continue;
558 }
559 }
560 if (signalInRange) {
563 if (
mean>maxGoodAmp) maxGoodAmp =
mean;
564 if (
mean<minGoodAmp) minGoodAmp =
mean;
565 }
568 if (slope > 0) {
569 prevSlope = slope;
570 }
571 ++nptGood;
572 averSlope += (slope - averSlope) / nptGood;
573 if (slope < 0) {
574 --nptGood;
575 }
576 if (averSlope < 0) {
577 nptGood = 0;
578 averSlope = 0.;
579 }
580 }
581
582 if (rms < 0.01) badPts++;
583 if (
mean > maxPointInFitRange) maxPointInFitRange =
mean;
584 if (rms > maxRMS) maxRMS =
rms;
585
587
589 if (prevSlope !=0 ) {
590 double R = slope / prevSlope;
591 if (R>0.4 && R<2.5) {
594 prevSlope = (prevSlope + slope)/2.;
595 }
596 } else {
598 if (R>0.4 && R<2.5) {
599 prevSlope = slope;
600 }
601 }
602 }
else if (nptGood < 3 && signalInRange && charge >=
m_linfitMax[gain_ind] &&
mean <
m_maxAmp[gain_ind] && prevSlope != 0) {
603 double R = slope/prevSlope;
604 if (R>0.9 && R<1.1) {
607 ++nptGood;
608 }
609 }
610
611
613 gr->SetPointError(pt, 0.0, err);
614
616 grrms->SetPointError(pt, 0.0, rms);
617
619 }
620
621
622 for (
int i=npt-1;
i>=
pt; --
i)
624
625 slope = (prevMean > 0 && prevCharge > 0) ? prevMean/prevCharge :
m_defaultCalib[gain_ind];
626 fslope->SetParameter(0, slope);
629 << ros << "/" << drawer
630 << "/" << chan << "/" << gain
631 << " up to " << maxGoodCharge
632 << " pC; N good points " << nptGood);
633 }
634
635 if (deltaChargeMin>999.) deltaChargeMin = 0.1;
636 else deltaChargeMin *= 0.5;
637 gr->Fit(
"fslope",
"q",
"", minGoodCharge-deltaChargeMin, maxGoodCharge+deltaChargeMin);
639 << ros << "/" << drawer
640 << "/" << chan << "/" << gain
641 << " minC " << minGoodCharge
642 << " maxC " << maxGoodCharge
643 << " deltaC " << deltaChargeMin*2
644 << " minA " << minGoodAmp
645 << " maxA " << maxGoodAmp
646 << " nptG " << nptGood
647 << " slope " << slope
648 << " fit " << fslope->GetParameter(0));
649 if (nptGood < 3) {
650 const char * bms[3] = {"No good points to fit ",
651 "Only one point to fit ",
652 "Only two points to fit "};
653 const char * ems[3] = {" => Put zero slope",
654 "",
655 ""};
657 << ros << "/" << drawer
658 << "/" << chan << "/" << gain
661 << " slope " << slope
662 << " fit " << fslope->GetParameter(0)
663 << ems[nptGood]);
664 }
665 slope = (nptGood > 0) ? fslope->GetParameter(0) : 0.;
666
667 if (slope > 0.) {
668 if (maxGoodCharge > minGoodCharge) {
669 averSlope = (maxGoodAmp-minGoodAmp)/(maxGoodCharge-minGoodCharge);
670 if (slope > 2.5 * averSlope) {
672 << ros << "/" << drawer
673 << "/" << chan << "/" << gain
676 << " avslope " << averSlope
677 << " fit " << slope
678 << " => Put zero slope");
679 slope = 0.;
680 }
681 } else if (maxCharge > minCharge) {
682 averSlope = (maxAmp-minAmp)/(maxCharge-minCharge);
683 if (slope > 10. * averSlope) {
685 << ros << "/" << drawer
686 << "/" << chan << "/" << gain
689 << " AVslope " << averSlope
690 << " fit " << slope
691 << " => Put zero slope");
692 slope = 0.;
693 }
694 }
695 }
696
698
699
700 if (ndigerr == 0) {
702 }
703
705 if (fslope->GetNDF() == 0)
707 else
709
710
711 if (TMath::Prob(fslope->GetChisquare(), fslope->GetNDF()) > 2.e-6) {
713 }
714
715
717
718
719
722
723
724
725
726
727
728
729 ratio = (slope / meanCalib[gain_ind]);
731
732
733
734 if (maxPointInFitRange > 600) {
736 }
737
738
739
740 if (maxRMS < 5.0) {
742 }
743
744
745
746
747
748
751 }
752
753
756 }
757
758
759
760
761 int NoStuckBit = 1;
762 for(
int i = 0;
i <
NBITS;
i++) {
763
765
766
767
769 NoStuckBit = 0;
771 << ros << " " << drawer << " " << chan << " " << gain << " " << i << "\n");
772
773 }
774
777 NoStuckBit = 0;
779 << ros << " " << drawer << " " << chan << " " << gain << " " << i << "\n");
780 }
781 }
782
783
784 if(NoStuckBit) {
786 }
787
789 grrms->SetName(
"scan_" +
arrayString(ros, drawer, chan, gain));
790
793 }
794 }
795 return StatusCode::SUCCESS;
796}
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
void mean(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")
double R(const INavigable4Momentum *p1, const double v_eta, const double v_phi)
1.14.0