1 /**************************************************************************
2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ////////////////////////////////////////////////////////////
19 // Author: Henrik Tydesjo //
20 // This class is used in the detector algorithm framework //
21 // to process the data stored in special container files //
22 // (see AliITSOnlineSPDscan). For instance, minimum //
23 // threshold values can be extracted. //
24 ////////////////////////////////////////////////////////////
26 #include "AliITSOnlineSPDscanAnalyzer.h"
27 #include "AliITSOnlineSPDscan.h"
28 #include "AliITSOnlineSPDscanSingle.h"
29 #include "AliITSOnlineSPDscanMultiple.h"
30 #include "AliITSOnlineSPDscanMeanTh.h"
31 #include "AliITSOnlineCalibrationSPDhandler.h"
32 #include "AliITSRawStreamSPD.h"
43 Double_t itsSpdErrorf(Double_t *x, Double_t *par){
44 if (par[2]<0) par[2]=0;
45 Double_t val = par[2]+(0.12*256*32-par[2])*(0.5+0.5*TMath::Erf((x[0]-par[0])/par[1]/sqrt(2.)));
48 //Double_t itsSpdErrorfOrig(Double_t *x, Double_t *par){
49 // return 0.5+0.5*TMath::Erf((x[0]-par[0])/par[1]/sqrt(2.));
51 //_________________________________________________________________________
52 Double_t itsSpdScurveForMeanTh(Double_t *x, Double_t *par){
53 if (par[2]<0) par[2]=0;
54 Double_t val = 1.- par[2]*(1.-TMath::Erf((x[0]-par[0])/par[1]/sqrt(2.)));
55 // Double_t val = par[2]+(0.12*256*32-par[2])*(0.5+0.5*TMath::Erf((x[0]-par[0])/par[1]/sqrt(2.)));
59 //_________________________________________________________________________
60 AliITSOnlineSPDscanAnalyzer::AliITSOnlineSPDscanAnalyzer(const Char_t *fileName, AliITSOnlineCalibrationSPDhandler *handler, Bool_t readFromGridFile) :
61 fType(99),fDacId(99),fFileName(fileName),fScanObj(NULL),fHandler(handler),fTriggers(NULL),fTPeff(0),fTPeffHS(NULL),fDeadPixel(0),fDeadPixelHS(NULL),fNoisyPixel(0),fNoisyPixelHS(NULL),
62 fOverWrite(kFALSE),fNoiseThreshold(0.01),fNoiseMinimumEvents(100),
63 fMinNrStepsBeforeIncrease(5),fMinIncreaseFromBaseLine(2),fStepDownDacSafe(5),fMaxBaseLineLevel(10)
66 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
67 for (UInt_t hs=0; hs<6; hs++) {
68 fMeanMultiplicity[hs][chipNr]=NULL;
69 fHitEventEfficiency[hs][chipNr]=NULL;
72 for (UInt_t hs=0; hs<6; hs++) {
74 fDeadPixelChip[hs]=NULL;
75 fNoisyPixelChip[hs]=NULL;
78 for (UInt_t mod=0; mod<240; mod++) {
79 fbModuleScanned[mod]=kFALSE;
82 Init(readFromGridFile);
84 //_________________________________________________________________________
85 AliITSOnlineSPDscanAnalyzer::AliITSOnlineSPDscanAnalyzer(const AliITSOnlineSPDscanAnalyzer& handle) :
86 fType(99),fDacId(99),fFileName("."),fScanObj(NULL),fHandler(NULL),fTriggers(NULL),fTPeff(0),fTPeffHS(NULL),fDeadPixel(0),fDeadPixelHS(NULL),fNoisyPixel(0),fNoisyPixelHS(NULL),
87 fOverWrite(kFALSE),fNoiseThreshold(0.01),fNoiseMinimumEvents(100),
88 fMinNrStepsBeforeIncrease(5),fMinIncreaseFromBaseLine(2),fStepDownDacSafe(5),fMaxBaseLineLevel(10)
90 // copy constructor, only copies the filename and params (not the processed data)
91 fFileName=handle.fFileName;
92 fOverWrite=handle.fOverWrite;
93 fNoiseThreshold=handle.fNoiseThreshold;
94 fNoiseMinimumEvents=handle.fNoiseMinimumEvents;
95 fMinNrStepsBeforeIncrease=handle.fMinNrStepsBeforeIncrease;
96 fMinIncreaseFromBaseLine=handle.fMinIncreaseFromBaseLine;
97 fStepDownDacSafe=handle.fStepDownDacSafe;
98 fMaxBaseLineLevel=handle.fMaxBaseLineLevel;
100 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
101 for (UInt_t hs=0; hs<6; hs++) {
102 fMeanMultiplicity[hs][chipNr]=NULL;
103 fHitEventEfficiency[hs][chipNr]=NULL;
106 for (UInt_t hs=0; hs<6; hs++) {
108 fDeadPixelChip[hs]=NULL;
109 fNoisyPixelChip[hs]=NULL;
112 for (UInt_t mod=0; mod<240; mod++) {
113 fbModuleScanned[mod]=kFALSE;
118 //_________________________________________________________________________
119 AliITSOnlineSPDscanAnalyzer::~AliITSOnlineSPDscanAnalyzer() {
121 for (UInt_t hs=0; hs<6; hs++) {
122 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
123 if (fMeanMultiplicity[hs][chipNr]!=NULL) {
124 delete fMeanMultiplicity[hs][chipNr];
125 fMeanMultiplicity[hs][chipNr]=NULL;
127 if (fHitEventEfficiency[hs][chipNr]!=NULL) {
128 delete fHitEventEfficiency[hs][chipNr];
129 fHitEventEfficiency[hs][chipNr]=NULL;
134 if (fTriggers!=NULL) {
139 DeleteUniformityHistograms();
141 if (fScanObj!=NULL) {
146 //_________________________________________________________________________
147 AliITSOnlineSPDscanAnalyzer& AliITSOnlineSPDscanAnalyzer::operator=(const AliITSOnlineSPDscanAnalyzer& handle) {
148 // assignment operator, only copies the filename and params (not the processed data)
150 for (UInt_t hs=0; hs<6; hs++) {
151 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
152 if (fMeanMultiplicity[hs][chipNr]!=NULL) {
153 delete fMeanMultiplicity[hs][chipNr];
155 if (fHitEventEfficiency[hs][chipNr]!=NULL) {
156 delete fHitEventEfficiency[hs][chipNr];
160 if (fTriggers!=NULL) {
165 DeleteUniformityHistograms();
167 if (fScanObj!=NULL) {
172 fFileName=handle.fFileName;
173 fOverWrite=handle.fOverWrite;
174 fNoiseThreshold=handle.fNoiseThreshold;
175 fNoiseMinimumEvents=handle.fNoiseMinimumEvents;
176 fMinNrStepsBeforeIncrease=handle.fMinNrStepsBeforeIncrease;
177 fMinIncreaseFromBaseLine=handle.fMinIncreaseFromBaseLine;
178 fStepDownDacSafe=handle.fStepDownDacSafe;
179 fMaxBaseLineLevel=handle.fMaxBaseLineLevel;
181 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
182 for (UInt_t hs=0; hs<6; hs++) {
183 fMeanMultiplicity[hs][chipNr]=NULL;
184 fHitEventEfficiency[hs][chipNr]=NULL;
187 for (UInt_t mod=0; mod<240; mod++) {
188 fbModuleScanned[mod]=kFALSE;
200 //_________________________________________________________________________
201 void AliITSOnlineSPDscanAnalyzer::Init(Bool_t readFromGridFile) {
202 // first checks type of container and then initializes container obj
203 if (!readFromGridFile) {
204 FILE* fp0 = fopen(fFileName.Data(), "r");
213 fScanObj = new AliITSOnlineSPDscan(fFileName.Data(),readFromGridFile);
214 fType = fScanObj->GetType();
221 fScanObj = new AliITSOnlineSPDscanSingle(fFileName.Data(),readFromGridFile);
226 fScanObj = new AliITSOnlineSPDscanMultiple(fFileName.Data(),readFromGridFile);
227 fDacId = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacId();
230 fScanObj = new AliITSOnlineSPDscanMeanTh(fFileName.Data(),readFromGridFile);
231 fDacId = ((AliITSOnlineSPDscanMeanTh*)fScanObj)->GetDacId();
234 Error("AliITSOnlineSPDscanAnalyzer::Init","Type %d not defined!",fType);
241 //_________________________________________________________________________
242 void AliITSOnlineSPDscanAnalyzer::SetParam(const Char_t *pname, const Char_t *pval) {
244 TString name = pname;
246 if (name.CompareTo("fOverWrite")==0) {
247 if (val.CompareTo("YES")==0 || val.CompareTo("1")==0) {
250 else fOverWrite = kFALSE;
252 else if (name.CompareTo("fNoiseThreshold")==0) {
253 fNoiseThreshold = val.Atof();
255 else if (name.CompareTo("fNoiseMinimumEvents")==0) {
256 fNoiseMinimumEvents = val.Atoi();
258 else if (name.CompareTo("fMinNrStepsBeforeIncrease")==0) {
259 fMinNrStepsBeforeIncrease = val.Atoi();
261 else if (name.CompareTo("fMinIncreaseFromBaseLine")==0) {
262 fMinIncreaseFromBaseLine = val.Atof();
264 else if (name.CompareTo("fStepDownDacSafe")==0) {
265 fStepDownDacSafe = val.Atoi();
267 else if (name.CompareTo("fMaxBaseLineLevel")==0) {
268 fMaxBaseLineLevel = val.Atof();
271 Error("AliITSOnlineSPDscanAnalyzer::SetParam","Parameter %s in configuration file unknown.",name.Data());
274 //_________________________________________________________________________
275 void AliITSOnlineSPDscanAnalyzer::ReadParamsFromLocation(const Char_t *dirName) {
276 // opens file (default name) in dir dirName and reads parameters from it
277 TString paramsFileName = Form("%s/standal_params.txt",dirName);
279 paramsFile.open(paramsFileName, ifstream::in);
280 if (paramsFile.fail()) {
281 printf("No config file (%s) present. Using default tuning parameters.\n",paramsFileName.Data());
287 paramsFile >> paramN;
288 if (paramsFile.eof()) break;
289 paramsFile >> paramV;
290 SetParam(paramN,paramV);
291 if (paramsFile.eof()) break;
296 //_________________________________________________________________________
297 Bool_t AliITSOnlineSPDscanAnalyzer::IsChipPresent(UInt_t hs, UInt_t chipNr) {
298 // is the chip present?
299 if (fScanObj==NULL) {
300 Warning("AliITSOnlineSPDscanAnalyzer::IsChipPresent","No data!");
303 return fScanObj->GetChipPresent(hs,chipNr);
305 //_________________________________________________________________________
306 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels(/*Char_t *oldcalibDir*/) {
307 // process dead pixel data, for uniformity scan,
308 // NB: This will not be the general way of finding dead pixels.
309 if (fScanObj==NULL) {
310 Warning("AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels","No data!");
313 // should be type kUNIMA
315 Warning("AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels","Dead pixels only for scan type %d.",kUNIMA);
318 // handler should be initialized
319 if (fHandler==NULL) {
320 Error("AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels","Calibration handler is not initialized!");
324 UInt_t routerNr = fScanObj->GetRouterNr();
325 UInt_t rowStart = fScanObj->GetRowStart();
326 UInt_t rowEnd = fScanObj->GetRowEnd();
327 for (UInt_t hs=0; hs<6; hs++) {
328 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
329 if (fScanObj->GetChipPresent(hs,chipNr) && fScanObj->GetAverageMultiplicity(0,hs,chipNr)>0) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
330 if (fOverWrite) {fHandler->ResetDeadForChip(routerNr,hs,chipNr);}
331 for (UInt_t col=0; col<32; col++) {
332 for (UInt_t row=rowStart; row<=rowEnd; row++) {
333 if (col!=1 && col!=9 && col!=17 && col!=25) { //exclude test columns!!!
334 if (fScanObj->GetHits(0,hs,chipNr,col,row)==0) {
335 fHandler->SetDeadPixel(routerNr,hs,chipNr,col,row);
345 //_________________________________________________________________________
346 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessUniformity() {
347 // process uniformity scan data (thanks to Roberta Ferretti for providing this method)
348 if (fScanObj==NULL) {
349 Warning("AliITSOnlineSPDscanAnalyzer::ProcessUniformity","No data!");
352 // should be type kUNIMA
354 Warning("AliITSOnlineSPDscanAnalyzer::ProcessUniformity","Only for scan type %d.",kUNIMA);
358 CreateUniformityHistograms(); // create all histograms that will be filled here
360 // UInt_t routerNr = fScanObj->GetRouterNr();
361 UInt_t rowStart = fScanObj->GetRowStart();
362 UInt_t rowEnd = fScanObj->GetRowEnd();
363 UInt_t nrTriggers = fScanObj->GetTriggers(0)/(rowEnd-rowStart+1);
368 UInt_t numChipsActive=0;
370 for (UInt_t hs=0; hs<6; hs++) {
371 Float_t pixel100hs=0;
374 UInt_t numChipsActiveHS=0;
376 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
377 Float_t pixel100chip=0;
379 Float_t pixelNchip=0;
381 if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
385 for (UInt_t col=0; col<32; col++) {
386 for (UInt_t row=rowStart; row<=rowEnd; row++) {
387 if (col!=1 && col!=9 && col!=17 && col!=25) { //exclude test columns!!!
389 if (fScanObj->GetHits(0,hs,chipNr,col,row)==nrTriggers) {
394 if (fScanObj->GetHits(0,hs,chipNr,col,row)==0) {
399 if (fScanObj->GetHits(0,hs,chipNr,col,row)>nrTriggers) {
408 Float_t tPeffChip=(pixel100chip/(28*(rowEnd-rowStart+1)))*100;
409 fTPeffChip[hs]->Fill(chipNr,tPeffChip);
411 Float_t deadPixelChip=(zerichip/(28*(rowEnd-rowStart+1)))*100;
412 fDeadPixelChip[hs]->Fill(chipNr,deadPixelChip);
414 Float_t noisyPixelChip=(pixelNchip/(28*(rowEnd-rowStart+1)))*100;
415 fNoisyPixelChip[hs]->Fill(chipNr,noisyPixelChip);
419 Float_t tPeffHS=(pixel100hs/(28*numChipsActiveHS*(rowEnd-rowStart+1)))*100;
420 fTPeffHS->Fill(hs,tPeffHS);
422 Float_t deadPixelHS=(zerihs/(28*numChipsActiveHS*(rowEnd-rowStart+1)))*100;
423 fDeadPixelHS->Fill(hs,deadPixelHS);
425 Float_t noisyPixelHS=(pixelNhs/(28*numChipsActiveHS*(rowEnd-rowStart+1)))*100;
426 fNoisyPixelHS->Fill(hs,noisyPixelHS);
429 fTPeff=(pixel100/(28*numChipsActive*(rowEnd-rowStart+1)))*100;
430 fDeadPixel=(zeri/(28*numChipsActive*(rowEnd-rowStart+1)))*100;
431 fNoisyPixel=(pixelN/(28*numChipsActive*(rowEnd-rowStart+1)))*100;
435 //_________________________________________________________________________
436 void AliITSOnlineSPDscanAnalyzer::CreateUniformityHistograms() {
437 // create uniformity histograms to be filled by "ProcessUniformity" method
438 DeleteUniformityHistograms(); // make sure no old histograms are lying around...
439 UInt_t eq = GetRouterNr();
442 label = Form("Ratio of 'Good' Pixels Per HS (eq %d)",eq);
443 fTPeffHS = new TH1F(label.Data(),label.Data(),6,-0.5,5.5);
444 fTPeffHS->SetXTitle("hs");
445 fTPeffHS->SetYTitle("ratio [%]");
446 fTPeffHS->SetFillColor(kBlue);
447 fTPeffHS->SetStats(0);
449 label = Form("Ratio of 'Dead' Pixels Per HS (eq %d)",eq);
450 fDeadPixelHS = new TH1F(label.Data(),label.Data(),6,-0.5,5.5);
451 fDeadPixelHS->SetXTitle("hs");
452 fDeadPixelHS->SetYTitle("ratio [%]");
453 fDeadPixelHS->SetFillColor(kBlue);
454 fDeadPixelHS->SetStats(0);
456 label = Form("Ratio of 'Noisy' Pixels Per HS (eq %d)",eq);
457 fNoisyPixelHS = new TH1F(label.Data(),label.Data(),6,-0.5,5.5);
458 fNoisyPixelHS->SetXTitle("hs");
459 fNoisyPixelHS->SetYTitle("ratio [%]");
460 fNoisyPixelHS->SetFillColor(kBlue);
461 fNoisyPixelHS->SetStats(0);
463 for (UInt_t hs=0; hs<6; hs++) {
464 label = Form("Ratio of 'Good' Pixels Per Chip (eq %d, hs %d)",eq,hs);
465 fTPeffChip[hs] = new TH1F(label.Data(),label.Data(),10,-0.5,9.5);
466 fTPeffChip[hs]->SetXTitle("chip");
467 fTPeffChip[hs]->SetYTitle("ratio [%]");
468 fTPeffChip[hs]->SetFillColor(kBlue);
469 fTPeffChip[hs]->SetStats(0);
471 label = Form("Ratio of 'Dead' Pixels Per Chip (eq %d, hs %d)",eq,hs);
472 fDeadPixelChip[hs] = new TH1F(label.Data(),label.Data(),10,-0.5,9.5);
473 fDeadPixelChip[hs]->SetXTitle("chip");
474 fDeadPixelChip[hs]->SetYTitle("ratio [%]");
475 fDeadPixelChip[hs]->SetFillColor(kBlue);
476 fDeadPixelChip[hs]->SetStats(0);
478 label = Form("Ratio of 'Noisy' Pixels Per Chip (eq %d, hs %d)",eq,hs);
479 fNoisyPixelChip[hs] = new TH1F(label.Data(),label.Data(),10,-0.5,9.5);
480 fNoisyPixelChip[hs]->SetXTitle("chip");
481 fNoisyPixelChip[hs]->SetYTitle("ratio [%]");
482 fNoisyPixelChip[hs]->SetFillColor(kBlue);
483 fNoisyPixelChip[hs]->SetStats(0);
487 //_________________________________________________________________________
488 void AliITSOnlineSPDscanAnalyzer::DeleteUniformityHistograms() {
489 // remove uniformity histograms if they are created
490 if (fTPeffHS!=NULL) {
494 if (fDeadPixelHS!=NULL) {
498 if (fNoisyPixelHS!=NULL) {
499 delete fNoisyPixelHS;
502 for (UInt_t hs=0; hs<6; hs++) {
503 if (fTPeffChip[hs]!=NULL) {
504 delete fTPeffChip[hs];
507 if (fDeadPixelChip[hs]!=NULL) {
508 delete fDeadPixelChip[hs];
509 fDeadPixelChip[hs]=NULL;
511 if (fNoisyPixelChip[hs]!=NULL) {
512 delete fNoisyPixelChip[hs];
513 fNoisyPixelChip[hs]=NULL;
517 //_________________________________________________________________________
518 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels(/*Char_t *oldcalibDir*/) {
519 // process noisy pixel data
520 if (fScanObj==NULL) {
521 Warning("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","No data!");
524 // should be type kNOISE
525 if (fType != kNOISE) {
526 Warning("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","Noisy pixels only for scan type %d.",kNOISE);
529 // handler should be initialized
530 if (fHandler==NULL) {
531 Error("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","Calibration handler is not initialized!");
534 // check if enough statistics
535 if (fScanObj->GetTriggers(0)<fNoiseMinimumEvents) {
536 Warning("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","Process noisy: Too few events.");
540 UInt_t routerNr = fScanObj->GetRouterNr();
541 for (UInt_t hs=0; hs<6; hs++) {
542 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
543 if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
544 if (fOverWrite) {fHandler->ResetNoisyForChip(routerNr,hs,chipNr);}
545 for (UInt_t col=0; col<32; col++) {
546 for (UInt_t row=0; row<256; row++) {
547 if (fScanObj->GetHitsEfficiency(0,hs,chipNr,col,row)>fNoiseThreshold) {
548 fHandler->SetNoisyPixel(routerNr,hs,chipNr,col,row);
557 //_________________________________________________________________________
558 Int_t AliITSOnlineSPDscanAnalyzer::GetDelay(UInt_t hs, UInt_t chipNr) {
560 if (hs>=6 || chipNr>10) return -1;
561 if (fScanObj==NULL) {
562 Warning("AliITSOnlineSPDscanAnalyzer::GetDelay","No data!");
565 // should be type kDELAY or kDAC with id 42 (delay_ctrl)
566 if (fType!=kDELAY && (fType!=kDAC || fDacId!=42)) {
567 Warning("AliITSOnlineSPDscanAnalyzer::GetDelay","Delay only for scan type %d or %d and dac_id 42.",kDELAY,kDAC);
570 if (fMeanMultiplicity[hs][chipNr]==NULL) {
571 if (!ProcessMeanMultiplicity()) {
578 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
581 fMeanMultiplicity[hs][chipNr]->GetPoint(step,thisDac,thisMult);
582 if (thisMult > maxVal) {
589 return ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(maxStep);
596 //_________________________________________________________________________
597 Int_t AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima(UInt_t hs, UInt_t chipNr) {
598 // in case of a uniformity scan, returns the nr of noisy pixels, (here > 200 hits)
599 if (hs>=6 || chipNr>10) return -1;
600 if (fScanObj==NULL) {
601 Error("AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima","No data!");
604 // should be type kUNIMA
605 if (fType != kUNIMA) {
606 Error("AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima","Noisy pixels Unima only for scan type %d.",kUNIMA);
609 if (fScanObj->GetTriggers(0)!=25600) {
610 Error("AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima","Process noisy unima: Incorrect number of events (!=25600.");
615 if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
616 for (UInt_t col=0; col<32; col++) {
617 for (UInt_t row=0; row<256; row++) {
618 if (fScanObj->GetHits(0,hs,chipNr,col,row)>200) {
629 //_________________________________________________________________________
630 Int_t AliITSOnlineSPDscanAnalyzer::FindLastMinThDac(UInt_t hs, UInt_t chipNr) {
631 // returns dac value where fMinIncreaseFromBaseLine reached
632 if (hs>=6 || chipNr>10) return -1;
633 if (fMeanMultiplicity[hs][chipNr]==NULL) {
634 if (!ProcessMeanMultiplicity()) {
638 Double_t firstVal, dummy1;
639 fMeanMultiplicity[hs][chipNr]->GetPoint(0,dummy1,firstVal);
641 while (step<fScanObj->GetNSteps()-1) {
642 Double_t graphVal, dummy2;
643 fMeanMultiplicity[hs][chipNr]->GetPoint(step+1,dummy2,graphVal);
644 if (graphVal>firstVal+fMinIncreaseFromBaseLine) break;
647 if (step==fScanObj->GetNSteps()-1) return -1;
648 return ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step);
651 Int_t AliITSOnlineSPDscanAnalyzer::FindClosestLowerStep(Float_t dacValueInput) {
652 // returns step closest (lower) to a dacvalue
654 while (step<fScanObj->GetNSteps()-1) {
655 Int_t dacVal = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step+1);
656 if (dacVal>=dacValueInput) break;
661 //_________________________________________________________________________
662 Float_t AliITSOnlineSPDscanAnalyzer::GetCompareLine(UInt_t step, UInt_t hs, UInt_t chipNr, Float_t basePar2) {
663 // returns value to compare mean mult with (when finding min th)
664 if (hs>=6 || chipNr>10) return -1;
665 if (step<fMinNrStepsBeforeIncrease) return -1;
666 Float_t baseLine = basePar2;
667 if (baseLine<0) baseLine=0;
672 for (UInt_t st=1;st<2*step/3;st++) { // skip first point...
673 fMeanMultiplicity[hs][chipNr]->GetPoint(st,d,m);
675 baseS+=(m-baseLine)*(m-baseLine);
677 baseAdd=2*sqrt( baseS/(2*step/3-1) - (baseM/(2*step/3-1))*(baseM/(2*step/3-1)) );
678 baseAdd+=0.03; // magic number
679 if (baseAdd>fMinIncreaseFromBaseLine) baseAdd=fMinIncreaseFromBaseLine;
680 return baseLine + baseAdd;
683 Int_t AliITSOnlineSPDscanAnalyzer::GetMinTh(UInt_t hs, UInt_t chipNr) {
684 // calculates and returns the minimum threshold
685 if (hs>=6 || chipNr>10) return -1;
686 if (fScanObj==NULL) {
687 Error("AliITSOnlineSPDscanAnalyzer::GetMinTh","No data!");
690 // should be type kMINTH or kDAC with id 39 (pre_vth)
691 if (fType!=kMINTH && (fType!=kDAC || fDacId!=39)) {
692 Error("AliITSOnlineSPDscanAnalyzer::GetMinTh","MinTh only for scan type %d OR %d with dac_id 39.",kMINTH,kDAC);
695 if (fMeanMultiplicity[hs][chipNr]==NULL) {
696 if (!ProcessMeanMultiplicity()) {
701 Int_t lastDac = FindLastMinThDac(hs,chipNr);
703 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: Increase of Mean Multiplicity by %1.2f never reached.",hs,chipNr,fMinIncreaseFromBaseLine);
707 Int_t minDac = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(0);
708 TString funcName = Form("Fit minth func HS%d CHIP%d",hs,chipNr);
709 TF1 *minThFunc = new TF1(funcName.Data(),itsSpdErrorf,100,500,3);
710 minThFunc->SetParameter(0,lastDac+10);
711 minThFunc->SetParameter(1,2);
712 minThFunc->SetParameter(2,0);
713 minThFunc->SetParName(0,"Mean");
714 minThFunc->SetParName(1,"Sigma");
715 minThFunc->SetParName(2,"BaseLine");
716 minThFunc->SetLineWidth(1);
717 if (fMeanMultiplicity[hs][chipNr]==NULL) {
718 if (!ProcessMeanMultiplicity()) {
722 fMeanMultiplicity[hs][chipNr]->Fit(funcName,"Q0","",minDac,lastDac);
724 // Double_t mean = fMinThFunc[hs][chipNr]->GetParameter(0);
725 // Double_t sigma = fMinThFunc[hs][chipNr]->GetParameter(1);
726 Double_t baseLine = minThFunc->GetParameter(2);
729 if (baseLine>fMaxBaseLineLevel) {
730 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: BaseLine too large (%1.2f>%1.2f).",hs,chipNr,baseLine,fMaxBaseLineLevel);
733 UInt_t step=FindClosestLowerStep(lastDac);
734 Float_t compareLine = GetCompareLine(step,hs,chipNr,baseLine);
735 if (compareLine==-1) {
736 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: Not enough steps (%d<%d) before increase to get a compare line.",hs,chipNr,step,fMinNrStepsBeforeIncrease);
740 Double_t mult, dummy;
742 while (mult > compareLine && step>0) {
743 fMeanMultiplicity[hs][chipNr]->GetPoint(step,dummy,mult);
746 Int_t minth = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step+1)-fStepDownDacSafe;
752 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: Did not find a point below the compare line (%f).",hs,chipNr,compareLine);
756 //_________________________________________________________________________
757 TArrayI AliITSOnlineSPDscanAnalyzer::GetMeanTh(UInt_t hs, UInt_t chipNr) {
758 // calculates and returns the mean threshold
759 TArrayI fitresults(4);
760 if (hs>=6 || chipNr>10) return fitresults;
761 if (fScanObj==NULL) {
762 Error("AliITSOnlineSPDscanAnalyzer::GetMeanTh","No data!");
765 // should be type kMEANTH or kDAC with id 39
766 if (fType!=kMEANTH && (fType!=kDAC || fDacId!=105)) {
767 Error("AliITSOnlineSPDscanAnalyzer::GetMeanTh","MeanTh only for scan type %d OR %d with dac_id 105.",kMEANTH,kDAC);
770 if (fHitEventEfficiency[hs][chipNr]==NULL) {
771 printf("processing hit efficiency \n");
772 if (!ProcessHitEventEfficiency()) {
773 printf("...not processed!!\n");
778 fHitEventEfficiency[hs][chipNr]->GetPoint(fHitEventEfficiency[hs][chipNr]->GetN()-1,x,y);
780 fHitEventEfficiency[hs][chipNr]->GetPoint(0,x,y);
783 Double_t mean = 0.5*(min+max);
784 TString funcName = Form("Fit meanth func HS%d CHIP%d",hs,chipNr);
785 TF1 *minThFunc = new TF1(funcName.Data(),itsSpdScurveForMeanTh,min,max,3);
786 minThFunc->SetParameter(0,mean);
787 minThFunc->SetParameter(1,264); // 4 (mV) * 66 (el)
788 minThFunc->SetParLimits(1,100,1000); // minimum value is 1 mV (-> 100 in TPAmplitude)
789 minThFunc->SetParameter(2,0.4);
790 minThFunc->SetParName(0,"Mean");
791 minThFunc->SetParName(1,"Sigma");
792 minThFunc->SetParName(2,"Half");
793 minThFunc->SetLineWidth(1);
795 fHitEventEfficiency[hs][chipNr]->Fit(funcName,"Q","",min,max);
797 fitresults.AddAt((Int_t)minThFunc->GetParameter(0),0);
798 fitresults.AddAt((Int_t)minThFunc->GetParError(0),1);
799 fitresults.AddAt((Int_t)minThFunc->GetParameter(1),2);
800 fitresults.AddAt((Int_t)minThFunc->GetParError(1),3);
801 TLine *ly = new TLine((Double_t)min,0.5,(Double_t)fitresults.At(0),0.5); ly->SetLineStyle(6);
803 TLine *lx = new TLine((Double_t)fitresults.At(0),0.,(Double_t)fitresults.At(0),0.5);
811 //_________________________________________________________________________
812 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessMeanMultiplicity() {
813 // process mean multiplicity data
814 if (fScanObj==NULL) {
815 Error("AliITSOnlineSPDscanAnalyzer::ProcessMeanMultiplicity","No data!");
818 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
819 for (UInt_t hs=0; hs<6; hs++) {
820 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
821 // if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
823 if (fMeanMultiplicity[hs][chipNr]!=NULL) {
824 delete fMeanMultiplicity[hs][chipNr];
826 fMeanMultiplicity[hs][chipNr] = new TGraph();
828 Float_t multiplMean=fScanObj->GetAverageMultiplicity(step,hs,chipNr);
829 if (fType==kMINTH || fType==kMEANTH || fType==kDAC || fType==kDELAY) {
830 fMeanMultiplicity[hs][chipNr]->SetPoint(step,((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step),multiplMean);
833 fMeanMultiplicity[hs][chipNr]->SetPoint(step,0,multiplMean);
841 //_________________________________________________________________________
842 TGraph* AliITSOnlineSPDscanAnalyzer::GetMeanMultiplicityG(UInt_t hs, UInt_t chipNr) {
843 // returns mean multiplicity graph
844 if (hs>=6 || chipNr>10) return NULL;
845 if (fMeanMultiplicity[hs][chipNr]==NULL) {
846 if (!ProcessMeanMultiplicity()) {
850 return fMeanMultiplicity[hs][chipNr];
852 //_________________________________________________________________________
853 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessHitEventEfficiency() {
854 // process hit event efficiency
855 if (fScanObj==NULL) {
856 Error("AliITSOnlineSPDscanAnalyzer::ProcessHitEventEfficiency","No data!");
859 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
860 for (UInt_t hs=0; hs<6; hs++) {
861 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
862 // if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
864 if (fHitEventEfficiency[hs][chipNr]!=NULL) {
865 delete fHitEventEfficiency[hs][chipNr];
867 fHitEventEfficiency[hs][chipNr] = new TGraph();
869 Float_t efficiency=fScanObj->GetHitEventsEfficiency(step,hs,chipNr);
870 if (fType==kMINTH || fType==kDAC || fType==kDELAY) {
871 fHitEventEfficiency[hs][chipNr]->SetPoint(step,((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step),efficiency);
872 } else if(fType==kMEANTH){
873 fHitEventEfficiency[hs][chipNr]->SetPoint(step,((AliITSOnlineSPDscanMeanTh*)fScanObj)->GetTPAmp(step,hs),efficiency);
875 fHitEventEfficiency[hs][chipNr]->SetPoint(step,0,efficiency);
883 //_________________________________________________________________________
884 TGraph* AliITSOnlineSPDscanAnalyzer::GetHitEventEfficiencyG(UInt_t hs, UInt_t chipNr) {
885 // returns hit event efficiency graph
886 if (hs>=6 || chipNr>10) return NULL;
887 if (fHitEventEfficiency[hs][chipNr]==NULL) {
888 if (!ProcessHitEventEfficiency()) {
892 return fHitEventEfficiency[hs][chipNr];
894 //_________________________________________________________________________
895 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessNrTriggers() {
896 // process nr of triggers data
897 if (fScanObj==NULL) {
898 Error("AliITSOnlineSPDscanAnalyzer::ProcessNrTriggers","No data!");
901 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
903 if (fTriggers!=NULL) {
906 fTriggers = new TGraph();
908 if (fType==kMINTH || fType==kMEANTH || fType==kDAC || fType==kDELAY) {
909 fTriggers->SetPoint(step,((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step),fScanObj->GetTriggers(step));
912 fTriggers->SetPoint(step,0,fScanObj->GetTriggers(step));
917 //_________________________________________________________________________
918 TGraph* AliITSOnlineSPDscanAnalyzer::GetNrTriggersG() {
919 // returns nr of triggers graph
920 if (fTriggers==NULL) {
921 if (!ProcessNrTriggers()) {
927 //_________________________________________________________________________
928 Bool_t AliITSOnlineSPDscanAnalyzer::GetHalfStavePresent(UInt_t hs) {
929 // returns half stave present info
930 if (hs<6 && fScanObj!=NULL) {
932 for (Int_t chip=0; chip<10; chip++) {
933 chipstatus+=fScanObj->GetChipPresent(hs,chip);
935 if (chipstatus>0) return kTRUE;
939 //_________________________________________________________________________
940 UInt_t AliITSOnlineSPDscanAnalyzer::GetRouterNr() {
941 // returns the router nr of scan obj
942 if (fScanObj!=NULL) return fScanObj->GetRouterNr();
945 //_________________________________________________________________________
946 TH2F* AliITSOnlineSPDscanAnalyzer::GetHitMapTot(UInt_t step) {
947 // creates and returns a pointer to a hitmap histo (half sector style a la spdmood)
948 if (fScanObj==NULL) {
949 Error("AliITSOnlineSPDscanAnalyzer::GetHitMapTot","No data!");
953 if (fType==kMINTH || fType==kMEANTH || fType==kDAC) {
954 histoname = Form("Router %d , DAC %d",GetRouterNr(),((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step));
957 histoname = Form("Router %d ",GetRouterNr());
959 TH2F* fHitMapTot = new TH2F(histoname.Data(),histoname.Data(),32*10,-0.5,32*10-0.5,256*6,-0.5,256*6-0.5);
960 fHitMapTot->SetNdivisions(-10,"X");
961 fHitMapTot->SetNdivisions(-006,"Y");
962 fHitMapTot->SetTickLength(0,"X");
963 fHitMapTot->SetTickLength(0,"Y");
964 fHitMapTot->GetXaxis()->SetLabelColor(gStyle->GetCanvasColor());
965 fHitMapTot->GetYaxis()->SetLabelColor(gStyle->GetCanvasColor());
966 for (UInt_t hs=0; hs<6; hs++) {
967 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
968 for (UInt_t col=0; col<32; col++) {
969 for (UInt_t row=0; row<256; row++) {
970 fHitMapTot->Fill(chipNr*32+col,(5-hs)*256+row,fScanObj->GetHits(step,hs,chipNr,col,row));
977 //_________________________________________________________________________
978 TH2F* AliITSOnlineSPDscanAnalyzer::GetPhysicalHitMapTot(UInt_t step) {
979 // creates and returns a pointer to a hitmap histo (-> 1 equpment opened up)
981 if (fScanObj==NULL) {
982 Error("AliITSOnlineSPDscanAnalyzer::GetHitMapTot","No data!");
986 if (fType==kMINTH || fType==kMEANTH || fType==kDAC) {
987 histoname = Form("Router %d , DAC %d",GetRouterNr(),((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step));
990 histoname = Form("Router %d ",GetRouterNr());
992 TH2F* hPhysicalHitMapTot = new TH2F(histoname.Data(),histoname.Data(),32*10,-0.5,32*10-0.5,256*6,-0.5,256*6-0.5);
993 hPhysicalHitMapTot->SetNdivisions(-10,"X");
994 hPhysicalHitMapTot->SetNdivisions(-006,"Y");
995 hPhysicalHitMapTot->SetTickLength(0,"X");
996 hPhysicalHitMapTot->SetTickLength(0,"Y");
997 hPhysicalHitMapTot->GetXaxis()->SetLabelColor(gStyle->GetCanvasColor());
998 hPhysicalHitMapTot->GetYaxis()->SetLabelColor(gStyle->GetCanvasColor());
999 Int_t correctChip = -1;
1000 for (UInt_t hs=0; hs<6; hs++) {
1001 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
1002 if(GetRouterNr()<10) correctChip = 9-chipNr;
1003 else correctChip=chipNr;
1004 for (UInt_t col=0; col<32; col++) {
1005 for (UInt_t row=0; row<256; row++) {
1006 if(hs<2) hPhysicalHitMapTot->Fill(correctChip*32+col,(5-hs)*256+row,fScanObj->GetHits(step,hs,chipNr,col,row));
1007 else hPhysicalHitMapTot->Fill(correctChip*32+(31-col),(5-hs)*256+row,fScanObj->GetHits(step,hs,chipNr,col,row));
1012 return hPhysicalHitMapTot;
1015 //_________________________________________________________________________
1016 TH2F* AliITSOnlineSPDscanAnalyzer::GetHitMapChip(UInt_t step, UInt_t hs, UInt_t chip) {
1017 // creates and returns a pointer to a hitmap histo (chip style a la spdmood)
1018 if (fScanObj==NULL) {
1019 Error("AliITSOnlineSPDscanAnalyzer::GetHitMapChip","No data!");
1025 histoName = Form("fChipHisto_%d_%d_%d", GetRouterNr(), hs, chip);
1026 histoTitle = Form("Eq ID %d, Half Stave %d, Chip %d ", GetRouterNr(), hs, chip);
1028 TH2F *returnHisto = new TH2F(histoName.Data(), histoTitle.Data(), 32, -0.5, 31.5, 256, -0.5, 255.5);
1029 returnHisto->SetMinimum(0);
1030 for (UInt_t col=0; col<32; col++) {
1031 for (UInt_t row=0; row<256; row++) {
1032 if(hs<2) returnHisto->Fill(31-col,row,fScanObj->GetHits(step,hs,chip,col,row));
1033 else returnHisto->Fill(col,row,fScanObj->GetHits(step,hs,chip,col,row));