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"
45 Double_t itsSpdErrorf(Double_t *x, Double_t *par){
46 if (par[2]<0) par[2]=0;
47 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.)));
50 //Double_t itsSpdErrorfOrig(Double_t *x, Double_t *par){
51 // return 0.5+0.5*TMath::Erf((x[0]-par[0])/par[1]/sqrt(2.));
53 //_________________________________________________________________________
54 Double_t itsSpdScurveForMeanTh(Double_t *x, Double_t *par){
55 if (par[2]<0) par[2]=0;
56 Double_t val = 1.- par[2]*(1.-TMath::Erf((x[0]-par[0])/par[1]/sqrt(2.)));
57 // 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.)));
61 //_________________________________________________________________________
62 AliITSOnlineSPDscanAnalyzer::AliITSOnlineSPDscanAnalyzer(const Char_t *fileName, AliITSOnlineCalibrationSPDhandler *handler, Bool_t readFromGridFile) :
63 fType(99),fDacId(99),fFileName(fileName),fScanObj(NULL),fHandler(handler),fTriggers(NULL),fTPeff(0),fTPeffHS(NULL),fDeadPixel(0),fDeadPixelHS(NULL),fNoisyPixel(0),fNoisyPixelHS(NULL),
64 fOverWrite(kFALSE),fNoiseThreshold(0.01),fNoiseMinimumEvents(100),
65 fMinNrStepsBeforeIncrease(5),fMinIncreaseFromBaseLine(2),fStepDownDacSafe(5),fMaxBaseLineLevel(10)
68 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
69 for (UInt_t hs=0; hs<6; hs++) {
70 fMeanMultiplicity[hs][chipNr]=NULL;
71 fHitEventEfficiency[hs][chipNr]=NULL;
74 for (UInt_t hs=0; hs<6; hs++) {
76 fDeadPixelChip[hs]=NULL;
77 fNoisyPixelChip[hs]=NULL;
80 for (UInt_t mod=0; mod<240; mod++) {
81 fbModuleScanned[mod]=kFALSE;
84 Init(readFromGridFile);
86 //_________________________________________________________________________
87 AliITSOnlineSPDscanAnalyzer::AliITSOnlineSPDscanAnalyzer(const AliITSOnlineSPDscanAnalyzer& handle) :
88 fType(99),fDacId(99),fFileName("."),fScanObj(NULL),fHandler(NULL),fTriggers(NULL),fTPeff(0),fTPeffHS(NULL),fDeadPixel(0),fDeadPixelHS(NULL),fNoisyPixel(0),fNoisyPixelHS(NULL),
89 fOverWrite(kFALSE),fNoiseThreshold(0.01),fNoiseMinimumEvents(100),
90 fMinNrStepsBeforeIncrease(5),fMinIncreaseFromBaseLine(2),fStepDownDacSafe(5),fMaxBaseLineLevel(10)
92 // copy constructor, only copies the filename and params (not the processed data)
93 fFileName=handle.fFileName;
94 fOverWrite=handle.fOverWrite;
95 fNoiseThreshold=handle.fNoiseThreshold;
96 fNoiseMinimumEvents=handle.fNoiseMinimumEvents;
97 fMinNrStepsBeforeIncrease=handle.fMinNrStepsBeforeIncrease;
98 fMinIncreaseFromBaseLine=handle.fMinIncreaseFromBaseLine;
99 fStepDownDacSafe=handle.fStepDownDacSafe;
100 fMaxBaseLineLevel=handle.fMaxBaseLineLevel;
102 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
103 for (UInt_t hs=0; hs<6; hs++) {
104 fMeanMultiplicity[hs][chipNr]=NULL;
105 fHitEventEfficiency[hs][chipNr]=NULL;
108 for (UInt_t hs=0; hs<6; hs++) {
110 fDeadPixelChip[hs]=NULL;
111 fNoisyPixelChip[hs]=NULL;
114 for (UInt_t mod=0; mod<240; mod++) {
115 fbModuleScanned[mod]=kFALSE;
120 //_________________________________________________________________________
121 AliITSOnlineSPDscanAnalyzer::~AliITSOnlineSPDscanAnalyzer() {
123 for (UInt_t hs=0; hs<6; hs++) {
124 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
125 if (fMeanMultiplicity[hs][chipNr]!=NULL) {
126 delete fMeanMultiplicity[hs][chipNr];
127 fMeanMultiplicity[hs][chipNr]=NULL;
129 if (fHitEventEfficiency[hs][chipNr]!=NULL) {
130 delete fHitEventEfficiency[hs][chipNr];
131 fHitEventEfficiency[hs][chipNr]=NULL;
136 if (fTriggers!=NULL) {
141 DeleteUniformityHistograms();
143 if (fScanObj!=NULL) {
148 //_________________________________________________________________________
149 AliITSOnlineSPDscanAnalyzer& AliITSOnlineSPDscanAnalyzer::operator=(const AliITSOnlineSPDscanAnalyzer& handle) {
150 // assignment operator, only copies the filename and params (not the processed data)
152 for (UInt_t hs=0; hs<6; hs++) {
153 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
154 if (fMeanMultiplicity[hs][chipNr]!=NULL) {
155 delete fMeanMultiplicity[hs][chipNr];
157 if (fHitEventEfficiency[hs][chipNr]!=NULL) {
158 delete fHitEventEfficiency[hs][chipNr];
162 if (fTriggers!=NULL) {
167 DeleteUniformityHistograms();
169 if (fScanObj!=NULL) {
174 fFileName=handle.fFileName;
175 fOverWrite=handle.fOverWrite;
176 fNoiseThreshold=handle.fNoiseThreshold;
177 fNoiseMinimumEvents=handle.fNoiseMinimumEvents;
178 fMinNrStepsBeforeIncrease=handle.fMinNrStepsBeforeIncrease;
179 fMinIncreaseFromBaseLine=handle.fMinIncreaseFromBaseLine;
180 fStepDownDacSafe=handle.fStepDownDacSafe;
181 fMaxBaseLineLevel=handle.fMaxBaseLineLevel;
183 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
184 for (UInt_t hs=0; hs<6; hs++) {
185 fMeanMultiplicity[hs][chipNr]=NULL;
186 fHitEventEfficiency[hs][chipNr]=NULL;
189 for (UInt_t mod=0; mod<240; mod++) {
190 fbModuleScanned[mod]=kFALSE;
202 //_________________________________________________________________________
203 void AliITSOnlineSPDscanAnalyzer::Init(Bool_t readFromGridFile) {
204 // first checks type of container and then initializes container obj
205 if (!readFromGridFile) {
206 FILE* fp0 = fopen(fFileName.Data(), "r");
215 fScanObj = new AliITSOnlineSPDscan(fFileName.Data(),readFromGridFile);
216 fType = fScanObj->GetType();
223 fScanObj = new AliITSOnlineSPDscanSingle(fFileName.Data(),readFromGridFile);
228 fScanObj = new AliITSOnlineSPDscanMultiple(fFileName.Data(),readFromGridFile);
229 fDacId = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacId();
232 fScanObj = new AliITSOnlineSPDscanMeanTh(fFileName.Data(),readFromGridFile);
233 fDacId = ((AliITSOnlineSPDscanMeanTh*)fScanObj)->GetDacId();
236 Error("AliITSOnlineSPDscanAnalyzer::Init","Type %d not defined!",fType);
243 //_________________________________________________________________________
244 void AliITSOnlineSPDscanAnalyzer::SetParam(const Char_t *pname, const Char_t *pval) {
246 TString name = pname;
248 if (name.CompareTo("fOverWrite")==0) {
249 if (val.CompareTo("YES")==0 || val.CompareTo("1")==0) {
252 else fOverWrite = kFALSE;
254 else if (name.CompareTo("fNoiseThreshold")==0) {
255 fNoiseThreshold = val.Atof();
257 else if (name.CompareTo("fNoiseMinimumEvents")==0) {
258 fNoiseMinimumEvents = val.Atoi();
260 else if (name.CompareTo("fMinNrStepsBeforeIncrease")==0) {
261 fMinNrStepsBeforeIncrease = val.Atoi();
263 else if (name.CompareTo("fMinIncreaseFromBaseLine")==0) {
264 fMinIncreaseFromBaseLine = val.Atof();
266 else if (name.CompareTo("fStepDownDacSafe")==0) {
267 fStepDownDacSafe = val.Atoi();
269 else if (name.CompareTo("fMaxBaseLineLevel")==0) {
270 fMaxBaseLineLevel = val.Atof();
273 Error("AliITSOnlineSPDscanAnalyzer::SetParam","Parameter %s in configuration file unknown.",name.Data());
276 //_________________________________________________________________________
277 void AliITSOnlineSPDscanAnalyzer::ReadParamsFromLocation(const Char_t *dirName) {
278 // opens file (default name) in dir dirName and reads parameters from it
279 TString paramsFileName = Form("%s/standal_params.txt",dirName);
281 paramsFile.open(paramsFileName, ifstream::in);
282 if (paramsFile.fail()) {
283 printf("No config file (%s) present. Using default tuning parameters.\n",paramsFileName.Data());
289 paramsFile >> paramN;
290 if (paramsFile.eof()) break;
291 paramsFile >> paramV;
292 SetParam(paramN,paramV);
293 if (paramsFile.eof()) break;
298 //_________________________________________________________________________
299 Bool_t AliITSOnlineSPDscanAnalyzer::IsChipPresent(UInt_t hs, UInt_t chipNr) {
300 // is the chip present?
301 if (fScanObj==NULL) {
302 Warning("AliITSOnlineSPDscanAnalyzer::IsChipPresent","No data!");
305 return fScanObj->GetChipPresent(hs,chipNr);
307 //_________________________________________________________________________
308 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels(/*Char_t *oldcalibDir*/) {
309 // process dead pixel data, for uniformity scan,
310 // NB: This will not be the general way of finding dead pixels.
311 if (fScanObj==NULL) {
312 Warning("AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels","No data!");
315 // should be type kUNIMA
317 Warning("AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels","Dead pixels only for scan type %d.",kUNIMA);
320 // handler should be initialized
321 if (fHandler==NULL) {
322 Error("AliITSOnlineSPDscanAnalyzer::ProcessDeadPixels","Calibration handler is not initialized!");
326 UInt_t routerNr = fScanObj->GetRouterNr();
327 UInt_t rowStart = fScanObj->GetRowStart();
328 UInt_t rowEnd = fScanObj->GetRowEnd();
329 for (UInt_t hs=0; hs<6; hs++) {
330 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
331 if (fScanObj->GetChipPresent(hs,chipNr) && fScanObj->GetAverageMultiplicity(0,hs,chipNr)>0) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
332 if (fOverWrite) {fHandler->ResetDeadForChip(routerNr,hs,chipNr);}
333 for (UInt_t col=0; col<32; col++) {
334 for (UInt_t row=rowStart; row<=rowEnd; row++) {
335 if (col!=1 && col!=9 && col!=17 && col!=25) { //exclude test columns!!!
336 if (fScanObj->GetHits(0,hs,chipNr,col,row)==0) {
337 fHandler->SetDeadPixel(routerNr,hs,chipNr,col,row);
347 //_________________________________________________________________________
348 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessUniformity() {
349 // process uniformity scan data (thanks to Roberta Ferretti for providing this method)
350 if (fScanObj==NULL) {
351 Warning("AliITSOnlineSPDscanAnalyzer::ProcessUniformity","No data!");
354 // should be type kUNIMA
356 Warning("AliITSOnlineSPDscanAnalyzer::ProcessUniformity","Only for scan type %d.",kUNIMA);
360 CreateUniformityHistograms(); // create all histograms that will be filled here
362 // UInt_t routerNr = fScanObj->GetRouterNr();
363 UInt_t rowStart = fScanObj->GetRowStart();
364 UInt_t rowEnd = fScanObj->GetRowEnd();
365 UInt_t nrTriggers = fScanObj->GetTriggers(0)/(rowEnd-rowStart+1);
370 UInt_t numChipsActive=0;
372 for (UInt_t hs=0; hs<6; hs++) {
373 Float_t pixel100hs=0;
376 UInt_t numChipsActiveHS=0;
378 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
379 Float_t pixel100chip=0;
381 Float_t pixelNchip=0;
383 if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
387 for (UInt_t col=0; col<32; col++) {
388 for (UInt_t row=rowStart; row<=rowEnd; row++) {
389 if (col!=1 && col!=9 && col!=17 && col!=25) { //exclude test columns!!!
391 if (fScanObj->GetHits(0,hs,chipNr,col,row)==nrTriggers) {
396 if (fScanObj->GetHits(0,hs,chipNr,col,row)==0) {
401 if (fScanObj->GetHits(0,hs,chipNr,col,row)>nrTriggers) {
410 Float_t tPeffChip=(pixel100chip/(28*(rowEnd-rowStart+1)))*100;
411 fTPeffChip[hs]->Fill(chipNr,tPeffChip);
413 Float_t deadPixelChip=(zerichip/(28*(rowEnd-rowStart+1)))*100;
414 fDeadPixelChip[hs]->Fill(chipNr,deadPixelChip);
416 Float_t noisyPixelChip=(pixelNchip/(28*(rowEnd-rowStart+1)))*100;
417 fNoisyPixelChip[hs]->Fill(chipNr,noisyPixelChip);
421 Float_t tPeffHS=(pixel100hs/(28*numChipsActiveHS*(rowEnd-rowStart+1)))*100;
422 fTPeffHS->Fill(hs,tPeffHS);
424 Float_t deadPixelHS=(zerihs/(28*numChipsActiveHS*(rowEnd-rowStart+1)))*100;
425 fDeadPixelHS->Fill(hs,deadPixelHS);
427 Float_t noisyPixelHS=(pixelNhs/(28*numChipsActiveHS*(rowEnd-rowStart+1)))*100;
428 fNoisyPixelHS->Fill(hs,noisyPixelHS);
431 fTPeff=(pixel100/(28*numChipsActive*(rowEnd-rowStart+1)))*100;
432 fDeadPixel=(zeri/(28*numChipsActive*(rowEnd-rowStart+1)))*100;
433 fNoisyPixel=(pixelN/(28*numChipsActive*(rowEnd-rowStart+1)))*100;
437 //_________________________________________________________________________
438 void AliITSOnlineSPDscanAnalyzer::CreateUniformityHistograms() {
439 // create uniformity histograms to be filled by "ProcessUniformity" method
440 DeleteUniformityHistograms(); // make sure no old histograms are lying around...
441 UInt_t eq = GetRouterNr();
444 label = Form("Ratio of 'Good' Pixels Per HS (eq %d)",eq);
445 fTPeffHS = new TH1F(label.Data(),label.Data(),6,-0.5,5.5);
446 fTPeffHS->SetXTitle("hs");
447 fTPeffHS->SetYTitle("ratio [%]");
448 fTPeffHS->SetFillColor(kBlue);
449 fTPeffHS->SetStats(0);
451 label = Form("Ratio of 'Dead' Pixels Per HS (eq %d)",eq);
452 fDeadPixelHS = new TH1F(label.Data(),label.Data(),6,-0.5,5.5);
453 fDeadPixelHS->SetXTitle("hs");
454 fDeadPixelHS->SetYTitle("ratio [%]");
455 fDeadPixelHS->SetFillColor(kBlue);
456 fDeadPixelHS->SetStats(0);
458 label = Form("Ratio of 'Noisy' Pixels Per HS (eq %d)",eq);
459 fNoisyPixelHS = new TH1F(label.Data(),label.Data(),6,-0.5,5.5);
460 fNoisyPixelHS->SetXTitle("hs");
461 fNoisyPixelHS->SetYTitle("ratio [%]");
462 fNoisyPixelHS->SetFillColor(kBlue);
463 fNoisyPixelHS->SetStats(0);
465 for (UInt_t hs=0; hs<6; hs++) {
466 label = Form("Ratio of 'Good' Pixels Per Chip (eq %d, hs %d)",eq,hs);
467 fTPeffChip[hs] = new TH1F(label.Data(),label.Data(),10,-0.5,9.5);
468 fTPeffChip[hs]->SetXTitle("chip");
469 fTPeffChip[hs]->SetYTitle("ratio [%]");
470 fTPeffChip[hs]->SetFillColor(kBlue);
471 fTPeffChip[hs]->SetStats(0);
473 label = Form("Ratio of 'Dead' Pixels Per Chip (eq %d, hs %d)",eq,hs);
474 fDeadPixelChip[hs] = new TH1F(label.Data(),label.Data(),10,-0.5,9.5);
475 fDeadPixelChip[hs]->SetXTitle("chip");
476 fDeadPixelChip[hs]->SetYTitle("ratio [%]");
477 fDeadPixelChip[hs]->SetFillColor(kBlue);
478 fDeadPixelChip[hs]->SetStats(0);
480 label = Form("Ratio of 'Noisy' Pixels Per Chip (eq %d, hs %d)",eq,hs);
481 fNoisyPixelChip[hs] = new TH1F(label.Data(),label.Data(),10,-0.5,9.5);
482 fNoisyPixelChip[hs]->SetXTitle("chip");
483 fNoisyPixelChip[hs]->SetYTitle("ratio [%]");
484 fNoisyPixelChip[hs]->SetFillColor(kBlue);
485 fNoisyPixelChip[hs]->SetStats(0);
489 //_________________________________________________________________________
490 void AliITSOnlineSPDscanAnalyzer::DeleteUniformityHistograms() {
491 // remove uniformity histograms if they are created
492 if (fTPeffHS!=NULL) {
496 if (fDeadPixelHS!=NULL) {
500 if (fNoisyPixelHS!=NULL) {
501 delete fNoisyPixelHS;
504 for (UInt_t hs=0; hs<6; hs++) {
505 if (fTPeffChip[hs]!=NULL) {
506 delete fTPeffChip[hs];
509 if (fDeadPixelChip[hs]!=NULL) {
510 delete fDeadPixelChip[hs];
511 fDeadPixelChip[hs]=NULL;
513 if (fNoisyPixelChip[hs]!=NULL) {
514 delete fNoisyPixelChip[hs];
515 fNoisyPixelChip[hs]=NULL;
519 //_________________________________________________________________________
520 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels(/*Char_t *oldcalibDir*/) {
521 // process noisy pixel data
522 if (fScanObj==NULL) {
523 Warning("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","No data!");
526 // should be type kNOISE
527 if (fType != kNOISE) {
528 Warning("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","Noisy pixels only for scan type %d.",kNOISE);
531 // handler should be initialized
532 if (fHandler==NULL) {
533 Error("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","Calibration handler is not initialized!");
536 // check if enough statistics
537 if (fScanObj->GetTriggers(0)<fNoiseMinimumEvents) {
538 Warning("AliITSOnlineSPDscanAnalyzer::ProcessNoisyPixels","Process noisy: Too few events.");
542 UInt_t routerNr = fScanObj->GetRouterNr();
543 for (UInt_t hs=0; hs<6; hs++) {
544 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
545 if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
546 if (fOverWrite) {fHandler->ResetNoisyForChip(routerNr,hs,chipNr);}
547 for (UInt_t col=0; col<32; col++) {
548 for (UInt_t row=0; row<256; row++) {
549 if (fScanObj->GetHitsEfficiency(0,hs,chipNr,col,row)>fNoiseThreshold) {
550 fHandler->SetNoisyPixel(routerNr,hs,chipNr,col,row);
559 //_________________________________________________________________________
560 Int_t AliITSOnlineSPDscanAnalyzer::GetDelay(UInt_t hs, UInt_t chipNr) {
562 if (hs>=6 || chipNr>10) return -1;
563 if (fScanObj==NULL) {
564 Warning("AliITSOnlineSPDscanAnalyzer::GetDelay","No data!");
567 // should be type kDELAY or kDAC with id 42 (delay_ctrl)
568 if (fType!=kDELAY && (fType!=kDAC || fDacId!=42)) {
569 Warning("AliITSOnlineSPDscanAnalyzer::GetDelay","Delay only for scan type %d or %d and dac_id 42.",kDELAY,kDAC);
572 if (fMeanMultiplicity[hs][chipNr]==NULL) {
573 if (!ProcessMeanMultiplicity()) {
580 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
583 fMeanMultiplicity[hs][chipNr]->GetPoint(step,thisDac,thisMult);
584 if (thisMult > maxVal) {
591 return ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(maxStep);
598 //_________________________________________________________________________
599 Int_t AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima(UInt_t hs, UInt_t chipNr) {
600 // in case of a uniformity scan, returns the nr of noisy pixels, (here > 200 hits)
601 if (hs>=6 || chipNr>10) return -1;
602 if (fScanObj==NULL) {
603 Error("AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima","No data!");
606 // should be type kUNIMA
607 if (fType != kUNIMA) {
608 Error("AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima","Noisy pixels Unima only for scan type %d.",kUNIMA);
611 if (fScanObj->GetTriggers(0)!=25600) {
612 Error("AliITSOnlineSPDscanAnalyzer::GetNrNoisyUnima","Process noisy unima: Incorrect number of events (!=25600.");
617 if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
618 for (UInt_t col=0; col<32; col++) {
619 for (UInt_t row=0; row<256; row++) {
620 if (fScanObj->GetHits(0,hs,chipNr,col,row)>200) {
631 //_________________________________________________________________________
632 Int_t AliITSOnlineSPDscanAnalyzer::FindLastMinThDac(UInt_t hs, UInt_t chipNr) {
633 // returns dac value where fMinIncreaseFromBaseLine reached
634 if (hs>=6 || chipNr>10) return -1;
635 if (fMeanMultiplicity[hs][chipNr]==NULL) {
636 if (!ProcessMeanMultiplicity()) {
640 Double_t firstVal, dummy1;
641 fMeanMultiplicity[hs][chipNr]->GetPoint(0,dummy1,firstVal);
643 while (step<fScanObj->GetNSteps()-1) {
644 Double_t graphVal, dummy2;
645 fMeanMultiplicity[hs][chipNr]->GetPoint(step+1,dummy2,graphVal);
646 if (graphVal>firstVal+fMinIncreaseFromBaseLine) break;
649 if (step==fScanObj->GetNSteps()-1) return -1;
650 return ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step);
653 Int_t AliITSOnlineSPDscanAnalyzer::FindClosestLowerStep(Float_t dacValueInput) {
654 // returns step closest (lower) to a dacvalue
656 while (step<fScanObj->GetNSteps()-1) {
657 Int_t dacVal = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step+1);
658 if (dacVal>=dacValueInput) break;
663 //_________________________________________________________________________
664 Float_t AliITSOnlineSPDscanAnalyzer::GetCompareLine(UInt_t step, UInt_t hs, UInt_t chipNr, Float_t basePar2) {
665 // returns value to compare mean mult with (when finding min th)
666 if (hs>=6 || chipNr>10) return -1;
667 if (step<fMinNrStepsBeforeIncrease) return -1;
668 Float_t baseLine = basePar2;
669 if (baseLine<0) baseLine=0;
674 for (UInt_t st=1;st<2*step/3;st++) { // skip first point...
675 fMeanMultiplicity[hs][chipNr]->GetPoint(st,d,m);
677 baseS+=(m-baseLine)*(m-baseLine);
679 baseAdd=2*sqrt( baseS/(2*step/3-1) - (baseM/(2*step/3-1))*(baseM/(2*step/3-1)) );
680 baseAdd+=0.03; // magic number
681 if (baseAdd>fMinIncreaseFromBaseLine) baseAdd=fMinIncreaseFromBaseLine;
682 return baseLine + baseAdd;
685 Int_t AliITSOnlineSPDscanAnalyzer::GetMinTh(UInt_t hs, UInt_t chipNr) {
686 // calculates and returns the minimum threshold
687 if (hs>=6 || chipNr>10) return -1;
688 if (fScanObj==NULL) {
689 Error("AliITSOnlineSPDscanAnalyzer::GetMinTh","No data!");
692 // should be type kMINTH or kDAC with id 39 (pre_vth)
693 if (fType!=kMINTH && (fType!=kDAC || fDacId!=39)) {
694 Error("AliITSOnlineSPDscanAnalyzer::GetMinTh","MinTh only for scan type %d OR %d with dac_id 39.",kMINTH,kDAC);
697 if (fMeanMultiplicity[hs][chipNr]==NULL) {
698 if (!ProcessMeanMultiplicity()) {
703 Int_t lastDac = FindLastMinThDac(hs,chipNr);
705 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: Increase of Mean Multiplicity by %1.2f never reached.",hs,chipNr,fMinIncreaseFromBaseLine);
709 Int_t minDac = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(0);
710 TString funcName = Form("Fit minth func HS%d CHIP%d",hs,chipNr);
711 TF1 *minThFunc = new TF1(funcName.Data(),itsSpdErrorf,100,500,3);
712 minThFunc->SetParameter(0,lastDac+10);
713 minThFunc->SetParameter(1,2);
714 minThFunc->SetParameter(2,0);
715 minThFunc->SetParName(0,"Mean");
716 minThFunc->SetParName(1,"Sigma");
717 minThFunc->SetParName(2,"BaseLine");
718 minThFunc->SetLineWidth(1);
719 if (fMeanMultiplicity[hs][chipNr]==NULL) {
720 if (!ProcessMeanMultiplicity()) {
724 fMeanMultiplicity[hs][chipNr]->Fit(funcName,"Q0","",minDac,lastDac);
726 // Double_t mean = fMinThFunc[hs][chipNr]->GetParameter(0);
727 // Double_t sigma = fMinThFunc[hs][chipNr]->GetParameter(1);
728 Double_t baseLine = minThFunc->GetParameter(2);
731 if (baseLine>fMaxBaseLineLevel) {
732 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: BaseLine too large (%1.2f>%1.2f).",hs,chipNr,baseLine,fMaxBaseLineLevel);
735 UInt_t step=FindClosestLowerStep(lastDac);
736 Float_t compareLine = GetCompareLine(step,hs,chipNr,baseLine);
737 if (compareLine==-1) {
738 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: Not enough steps (%d<%d) before increase to get a compare line.",hs,chipNr,step,fMinNrStepsBeforeIncrease);
742 Double_t mult, dummy;
744 while (mult > compareLine && step>0) {
745 fMeanMultiplicity[hs][chipNr]->GetPoint(step,dummy,mult);
748 Int_t minth = ((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step+1)-fStepDownDacSafe;
754 Warning("AliITSOnlineSPDscanAnalyzer::GetMinTh","HS%d, Chip%d: Did not find a point below the compare line (%f).",hs,chipNr,compareLine);
758 //_________________________________________________________________________
759 TArrayI AliITSOnlineSPDscanAnalyzer::GetMeanTh(UInt_t hs, UInt_t chipNr) {
760 // calculates and returns the mean threshold
761 TArrayI fitresults(4);
762 if (hs>=6 || chipNr>10) return fitresults;
763 if (fScanObj==NULL) {
764 Error("AliITSOnlineSPDscanAnalyzer::GetMeanTh","No data!");
767 // should be type kMEANTH or kDAC with id 39
768 if (fType!=kMEANTH && (fType!=kDAC || fDacId!=105)) {
769 Error("AliITSOnlineSPDscanAnalyzer::GetMeanTh","MeanTh only for scan type %d OR %d with dac_id 105.",kMEANTH,kDAC);
772 if (fHitEventEfficiency[hs][chipNr]==NULL) {
773 printf("processing hit efficiency \n");
774 if (!ProcessHitEventEfficiency()) {
775 printf("...not processed!!\n");
780 fHitEventEfficiency[hs][chipNr]->GetPoint(fHitEventEfficiency[hs][chipNr]->GetN()-1,x,y);
782 fHitEventEfficiency[hs][chipNr]->GetPoint(0,x,y);
785 Double_t mean = 0.5*(min+max);
786 TString funcName = Form("Fit meanth func HS%d CHIP%d",hs,chipNr);
787 TF1 *minThFunc = new TF1(funcName.Data(),itsSpdScurveForMeanTh,min,max,3);
788 minThFunc->SetParameter(0,mean);
789 minThFunc->SetParameter(1,264); // 4 (mV) * 66 (el)
790 minThFunc->SetParLimits(1,100,1000); // minimum value is 1 mV (-> 100 in TPAmplitude)
791 minThFunc->SetParameter(2,0.4);
792 minThFunc->SetParName(0,"Mean");
793 minThFunc->SetParName(1,"Sigma");
794 minThFunc->SetParName(2,"Half");
795 minThFunc->SetLineWidth(1);
797 fHitEventEfficiency[hs][chipNr]->Fit(funcName,"Q","",min,max);
799 fitresults.AddAt((Int_t)minThFunc->GetParameter(0),0);
800 fitresults.AddAt((Int_t)minThFunc->GetParError(0),1);
801 fitresults.AddAt((Int_t)minThFunc->GetParameter(1),2);
802 fitresults.AddAt((Int_t)minThFunc->GetParError(1),3);
803 TLine *ly = new TLine((Double_t)min,0.5,(Double_t)fitresults.At(0),0.5); ly->SetLineStyle(6);
805 TLine *lx = new TLine((Double_t)fitresults.At(0),0.,(Double_t)fitresults.At(0),0.5);
813 //_________________________________________________________________________
814 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessMeanMultiplicity() {
815 // process mean multiplicity data
816 if (fScanObj==NULL) {
817 Error("AliITSOnlineSPDscanAnalyzer::ProcessMeanMultiplicity","No data!");
820 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
821 for (UInt_t hs=0; hs<6; hs++) {
822 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
823 // if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
825 if (fMeanMultiplicity[hs][chipNr]!=NULL) {
826 delete fMeanMultiplicity[hs][chipNr];
828 fMeanMultiplicity[hs][chipNr] = new TGraph();
830 Float_t multiplMean=fScanObj->GetAverageMultiplicity(step,hs,chipNr);
831 if (fType==kMINTH || fType==kMEANTH || fType==kDAC || fType==kDELAY) {
832 fMeanMultiplicity[hs][chipNr]->SetPoint(step,((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step),multiplMean);
835 fMeanMultiplicity[hs][chipNr]->SetPoint(step,0,multiplMean);
843 //_________________________________________________________________________
844 TGraph* AliITSOnlineSPDscanAnalyzer::GetMeanMultiplicityG(UInt_t hs, UInt_t chipNr) {
845 // returns mean multiplicity graph
846 if (hs>=6 || chipNr>10) return NULL;
847 if (fMeanMultiplicity[hs][chipNr]==NULL) {
848 if (!ProcessMeanMultiplicity()) {
852 return fMeanMultiplicity[hs][chipNr];
854 //_________________________________________________________________________
855 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessHitEventEfficiency() {
856 // process hit event efficiency
857 if (fScanObj==NULL) {
858 Error("AliITSOnlineSPDscanAnalyzer::ProcessHitEventEfficiency","No data!");
861 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
862 for (UInt_t hs=0; hs<6; hs++) {
863 for (UInt_t chipNr=0; chipNr<11; chipNr++) {
864 // if (fScanObj->GetChipPresent(hs,chipNr)) { // check the status of the chippresent parameter in the mood header!!!!!!!!!!!!!!!!!!!!!!!!!!!!
866 if (fHitEventEfficiency[hs][chipNr]!=NULL) {
867 delete fHitEventEfficiency[hs][chipNr];
869 fHitEventEfficiency[hs][chipNr] = new TGraph();
871 Float_t efficiency=fScanObj->GetHitEventsEfficiency(step,hs,chipNr);
872 if (fType==kMINTH || fType==kDAC || fType==kDELAY) {
873 fHitEventEfficiency[hs][chipNr]->SetPoint(step,((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step),efficiency);
874 } else if(fType==kMEANTH){
875 fHitEventEfficiency[hs][chipNr]->SetPoint(step,((AliITSOnlineSPDscanMeanTh*)fScanObj)->GetTPAmp(step,hs),efficiency);
877 fHitEventEfficiency[hs][chipNr]->SetPoint(step,0,efficiency);
885 //_________________________________________________________________________
886 TGraph* AliITSOnlineSPDscanAnalyzer::GetHitEventEfficiencyG(UInt_t hs, UInt_t chipNr) {
887 // returns hit event efficiency graph
888 if (hs>=6 || chipNr>10) return NULL;
889 if (fHitEventEfficiency[hs][chipNr]==NULL) {
890 if (!ProcessHitEventEfficiency()) {
894 return fHitEventEfficiency[hs][chipNr];
896 //_________________________________________________________________________
897 Bool_t AliITSOnlineSPDscanAnalyzer::ProcessNrTriggers() {
898 // process nr of triggers data
899 if (fScanObj==NULL) {
900 Error("AliITSOnlineSPDscanAnalyzer::ProcessNrTriggers","No data!");
903 for (UInt_t step=0; step<fScanObj->GetNSteps(); step++) {
905 if (fTriggers!=NULL) {
908 fTriggers = new TGraph();
910 if (fType==kMINTH || fType==kMEANTH || fType==kDAC || fType==kDELAY) {
911 fTriggers->SetPoint(step,((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step),fScanObj->GetTriggers(step));
914 fTriggers->SetPoint(step,0,fScanObj->GetTriggers(step));
919 //_________________________________________________________________________
920 TGraph* AliITSOnlineSPDscanAnalyzer::GetNrTriggersG() {
921 // returns nr of triggers graph
922 if (fTriggers==NULL) {
923 if (!ProcessNrTriggers()) {
929 //_________________________________________________________________________
930 Bool_t AliITSOnlineSPDscanAnalyzer::GetHalfStavePresent(UInt_t hs) {
931 // returns half stave present info
932 if (hs<6 && fScanObj!=NULL) {
934 for (Int_t chip=0; chip<10; chip++) {
935 chipstatus+=fScanObj->GetChipPresent(hs,chip);
937 if (chipstatus>0) return kTRUE;
941 //_________________________________________________________________________
942 UInt_t AliITSOnlineSPDscanAnalyzer::GetRouterNr() {
943 // returns the router nr of scan obj
944 if (fScanObj!=NULL) return fScanObj->GetRouterNr();
947 //_________________________________________________________________________
948 TH2F* AliITSOnlineSPDscanAnalyzer::GetHitMapTot(UInt_t step) {
949 // creates and returns a pointer to a hitmap histo (half sector style a la spdmood)
950 if (fScanObj==NULL) {
951 Error("AliITSOnlineSPDscanAnalyzer::GetHitMapTot","No data!");
955 if (fType==kMINTH || fType==kMEANTH || fType==kDAC) {
956 histoname = Form("Router %d , DAC %d",GetRouterNr(),((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step));
959 histoname = Form("Router %d ",GetRouterNr());
961 TH2F* fHitMapTot = new TH2F(histoname.Data(),histoname.Data(),32*10,-0.5,32*10-0.5,256*6,-0.5,256*6-0.5);
962 fHitMapTot->SetNdivisions(-10,"X");
963 fHitMapTot->SetNdivisions(-006,"Y");
964 fHitMapTot->SetTickLength(0,"X");
965 fHitMapTot->SetTickLength(0,"Y");
966 fHitMapTot->GetXaxis()->SetLabelColor(gStyle->GetCanvasColor());
967 fHitMapTot->GetYaxis()->SetLabelColor(gStyle->GetCanvasColor());
968 for (UInt_t hs=0; hs<6; hs++) {
969 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
970 for (UInt_t col=0; col<32; col++) {
971 for (UInt_t row=0; row<256; row++) {
972 fHitMapTot->Fill(chipNr*32+col,(5-hs)*256+row,fScanObj->GetHits(step,hs,chipNr,col,row));
979 //_________________________________________________________________________
980 TH2F* AliITSOnlineSPDscanAnalyzer::GetPhysicalHitMapTot(UInt_t step) {
981 // creates and returns a pointer to a hitmap histo (-> 1 equpment opened up)
983 if (fScanObj==NULL) {
984 Error("AliITSOnlineSPDscanAnalyzer::GetHitMapTot","No data!");
988 if (fType==kMINTH || fType==kMEANTH || fType==kDAC) {
989 histoname = Form("Router %d , DAC %d",GetRouterNr(),((AliITSOnlineSPDscanMultiple*)fScanObj)->GetDacValue(step));
992 histoname = Form("Router %d ",GetRouterNr());
994 TH2F* hPhysicalHitMapTot = new TH2F(histoname.Data(),histoname.Data(),32*10,-0.5,32*10-0.5,256*6,-0.5,256*6-0.5);
995 hPhysicalHitMapTot->SetNdivisions(-10,"X");
996 hPhysicalHitMapTot->SetNdivisions(-006,"Y");
997 hPhysicalHitMapTot->SetTickLength(0,"X");
998 hPhysicalHitMapTot->SetTickLength(0,"Y");
999 hPhysicalHitMapTot->GetXaxis()->SetLabelColor(gStyle->GetCanvasColor());
1000 hPhysicalHitMapTot->GetYaxis()->SetLabelColor(gStyle->GetCanvasColor());
1001 Int_t correctChip = -1;
1002 for (UInt_t hs=0; hs<6; hs++) {
1003 for (UInt_t chipNr=0; chipNr<10; chipNr++) {
1004 if(GetRouterNr()<10) correctChip = 9-chipNr;
1005 else correctChip=chipNr;
1006 for (UInt_t col=0; col<32; col++) {
1007 for (UInt_t row=0; row<256; row++) {
1008 if(hs<2) hPhysicalHitMapTot->Fill(correctChip*32+col,(5-hs)*256+row,fScanObj->GetHits(step,hs,chipNr,col,row));
1009 else hPhysicalHitMapTot->Fill(correctChip*32+(31-col),(5-hs)*256+row,fScanObj->GetHits(step,hs,chipNr,col,row));
1014 return hPhysicalHitMapTot;
1017 //_________________________________________________________________________
1018 TH2F* AliITSOnlineSPDscanAnalyzer::GetHitMapChip(UInt_t step, UInt_t hs, UInt_t chip) {
1019 // creates and returns a pointer to a hitmap histo (chip style a la spdmood)
1020 if (fScanObj==NULL) {
1021 Error("AliITSOnlineSPDscanAnalyzer::GetHitMapChip","No data!");
1027 histoName = Form("fChipHisto_%d_%d_%d", GetRouterNr(), hs, chip);
1028 histoTitle = Form("Eq ID %d, Half Stave %d, Chip %d ", GetRouterNr(), hs, chip);
1030 TH2F *returnHisto = new TH2F(histoName.Data(), histoTitle.Data(), 32, -0.5, 31.5, 256, -0.5, 255.5);
1031 returnHisto->SetMinimum(0);
1032 for (UInt_t col=0; col<32; col++) {
1033 for (UInt_t row=0; row<256; row++) {
1034 if(hs<2) returnHisto->Fill(31-col,row,fScanObj->GetHits(step,hs,chip,col,row));
1035 else returnHisto->Fill(col,row,fScanObj->GetHits(step,hs,chip,col,row));