1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 ///////////////////////////////////////////////////////////////////////////////
20 // TRD calibration class for parameters which saved per pad //
22 ///////////////////////////////////////////////////////////////////////////////
29 //#include "AliMathBase.h"
31 #include "AliTRDCalPad.h"
32 #include "AliTRDCalROC.h"
33 #include "AliTRDCalDet.h"
34 #include "AliTRDpadPlane.h"
35 #include "AliTRDgeometry.h"
37 ClassImp(AliTRDCalPad)
39 //_____________________________________________________________________________
40 AliTRDCalPad::AliTRDCalPad():TNamed()
43 // AliTRDCalPad default constructor
46 for (Int_t idet = 0; idet < kNdet; idet++) {
52 //_____________________________________________________________________________
53 AliTRDCalPad::AliTRDCalPad(const Text_t *name, const Text_t *title)
57 // AliTRDCalPad constructor
60 for (Int_t isec = 0; isec < kNsect; isec++) {
61 for (Int_t ipla = 0; ipla < kNplan; ipla++) {
62 for (Int_t icha = 0; icha < kNcham; icha++) {
63 Int_t idet = GetDet(ipla,icha,isec);
64 fROC[idet] = new AliTRDCalROC(ipla,icha);
71 //_____________________________________________________________________________
72 AliTRDCalPad::AliTRDCalPad(const AliTRDCalPad &c)
76 // AliTRDCalPad copy constructor
79 for (Int_t idet = 0; idet < kNdet; idet++) {
80 fROC[idet] = new AliTRDCalROC(*((AliTRDCalPad &) c).fROC[idet]);
85 //_____________________________________________________________________________
86 AliTRDCalPad::~AliTRDCalPad()
89 // AliTRDCalPad destructor
92 for (Int_t idet = 0; idet < kNdet; idet++) {
101 //_____________________________________________________________________________
102 AliTRDCalPad &AliTRDCalPad::operator=(const AliTRDCalPad &c)
105 // Assignment operator
108 if (this != &c) ((AliTRDCalPad &) c).Copy(*this);
113 //_____________________________________________________________________________
114 void AliTRDCalPad::Copy(TObject &c) const
120 for (Int_t idet = 0; idet < kNdet; idet++) {
121 if (((AliTRDCalPad &) c).fROC[idet]) {
122 delete ((AliTRDCalPad &) c).fROC[idet];
124 ((AliTRDCalPad &) c).fROC[idet] = new AliTRDCalROC();
126 fROC[idet]->Copy(*((AliTRDCalPad &) c).fROC[idet]);
134 //_____________________________________________________________________________
135 Bool_t AliTRDCalPad::ScaleROCs(const AliTRDCalDet* values)
138 // Scales ROCs of this class with the values from the class <values>
139 // Is used if an AliTRDCalPad object defines local variations of a parameter
140 // defined per detector using a AliTRDCalDet class
145 Bool_t result = kTRUE;
146 for (Int_t idet = 0; idet < kNdet; idet++) {
148 if(!fROC[idet]->Multiply(values->GetValue(idet))) result = kFALSE;
154 //_____________________________________________________________________________
155 void AliTRDCalPad::SetCalROC(Int_t det, AliTRDCalROC* calroc)
158 // Set the AliTRDCalROC to this one
163 for(Int_t icol = 0; icol < calroc->GetNcols(); icol++){
164 for(Int_t irow = 0; irow < calroc->GetNrows(); irow++){
165 fROC[det]->SetValue(icol,irow,calroc->GetValue(icol,irow));
171 //_____________________________________________________________________________
172 Double_t AliTRDCalPad::GetMeanRMS(Double_t &rms, const AliTRDCalDet *calDet, Int_t type)
175 // Calculate mean an RMS of all rocs
176 // If calDet correct the CalROC from the detector coefficient
177 // type == 0 for gain and vdrift
180 Double_t factor = 0.0;
181 if(type == 0) factor = 1.0;
182 Double_t sum = 0, sum2 = 0, n=0, val=0;
183 for (Int_t idet = 0; idet < kNdet; idet++) {
184 if(calDet) factor = calDet->GetValue(idet);
185 AliTRDCalROC *calRoc = fROC[idet];
187 for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
188 for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
189 if(type == 0) val = calRoc->GetValue(icol,irow)*factor;
190 else val = calRoc->GetValue(icol,irow)+factor;
199 Double_t mean = sum*n1;
200 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
204 //_____________________________________________________________________________
205 Double_t AliTRDCalPad::GetMean(const AliTRDCalDet *calDet, Int_t type
206 , AliTRDCalPad* const outlierPad)
209 // return mean of the mean of all ROCs
210 // If calDet correct the CalROC from the detector coefficient
211 // type == 0 for gain and vdrift
214 Double_t factor = 0.0;
215 if(type == 0) factor = 1.0;
218 for (Int_t idet = 0; idet < kNdet; idet++) {
219 if(calDet) factor = calDet->GetValue(idet);
220 AliTRDCalROC *calRoc = fROC[idet];
222 AliTRDCalROC* outlierROC = 0;
223 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
224 if(type == 0) arr[n] = calRoc->GetMean(outlierROC)*factor;
225 else arr[n] = calRoc->GetMean(outlierROC)+factor;
229 return TMath::Mean(n,arr);
232 //_____________________________________________________________________________
233 Double_t AliTRDCalPad::GetRMS(const AliTRDCalDet *calDet, Int_t type
234 , AliTRDCalPad* const outlierPad)
237 // return mean of the RMS of all ROCs
238 // If calDet correct the CalROC from the detector coefficient
239 // type == 0 for gain and vdrift
242 Double_t factor = 0.0;
243 if(type == 0) factor = 1.0;
246 for (Int_t idet = 0; idet < kNdet; idet++) {
247 if(calDet) factor = calDet->GetValue(idet);
248 AliTRDCalROC *calRoc = fROC[idet];
250 AliTRDCalROC* outlierROC = 0;
251 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
252 if(type == 0) arr[n] = calRoc->GetRMS(outlierROC)*factor;
253 else arr[n] = calRoc->GetRMS(outlierROC);
257 return TMath::Mean(n,arr);
260 //_____________________________________________________________________________
261 Double_t AliTRDCalPad::GetMedian(const AliTRDCalDet *calDet, Int_t type
262 , AliTRDCalPad* const outlierPad)
265 // return mean of the median of all ROCs
266 // If AliTRDCalDet, the correct the CalROC from the detector coefficient
267 // type == 0 for gain and vdrift
270 Double_t factor = 0.0;
271 if(type == 0) factor = 1.0;
274 for (Int_t idet = 0; idet < kNdet; idet++) {
275 if(calDet) factor = calDet->GetValue(idet);
276 AliTRDCalROC *calRoc = fROC[idet];
278 AliTRDCalROC* outlierROC = 0;
279 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
280 if(type == 0) arr[n] = calRoc->GetMedian(outlierROC)*factor;
281 else arr[n] = calRoc->GetMedian(outlierROC)+factor;
285 return TMath::Mean(n,arr);
288 //_____________________________________________________________________________
289 Double_t AliTRDCalPad::GetLTM(Double_t *sigma, Double_t fraction
290 , const AliTRDCalDet *calDet, Int_t type
291 , AliTRDCalPad* const outlierPad)
294 // return mean of the LTM and sigma of all ROCs
295 // If calDet correct the CalROC from the detector coefficient
296 // type == 0 for gain and vdrift
299 Double_t factor = 0.0;
300 if(type == 0) factor = 1.0;
301 Double_t arrm[kNdet];
302 Double_t arrs[kNdet];
306 for (Int_t idet = 0; idet < kNdet; idet++) {
307 if(calDet) factor = calDet->GetValue(idet);
308 AliTRDCalROC *calRoc = fROC[idet];
310 if ( sigma ) sTemp=arrs+n;
311 AliTRDCalROC* outlierROC = 0;
312 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
313 if(type == 0) arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC)*factor;
314 else arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC)+factor;
318 if ( sigma ) *sigma = TMath::Mean(n,arrs);
319 return TMath::Mean(n,arrm);
322 //_____________________________________________________________________________
323 TH1F * AliTRDCalPad::MakeHisto1D(const AliTRDCalDet *calDet, Int_t typedet
324 , Float_t min, Float_t max,Int_t type)
328 // type -1 = user defined range
329 // 0 = nsigma cut nsigma=min
330 // If calDet correct the CalROC from the detector coefficient
331 // typedet == 0 for gain and vdrift
332 // typedet == 1 for t0
334 Float_t kEpsilonr = 0.005;
336 Double_t factor = 0.0;
337 if(typedet == 0) factor = 1.0;
342 Float_t mean = GetMean(calDet,typedet);
344 if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
347 sigma = GetMeanRMS(rms,calDet,typedet);
349 Float_t nsigma = TMath::Abs(min);
351 if(sigma < kEpsilonr) sigma = kEpsilonr;
357 Float_t mean = GetMedian(calDet,typedet);
364 // LTM mean +- nsigma
367 Float_t mean = GetLTM(&sigma,max,calDet,typedet);
369 if(sigma < kEpsilonr) sigma = kEpsilonr;
375 snprintf(name,1000,"%s Pad 1D",GetTitle());
376 TH1F * his = new TH1F(name,name,100, min,max);
377 for (Int_t idet = 0; idet < kNdet; idet++) {
378 if(calDet) factor = calDet->GetValue(idet);
380 for (Int_t irow=0; irow<fROC[idet]->GetNrows(); irow++){
381 for (Int_t icol=0; icol<fROC[idet]->GetNcols(); icol++){
382 if(typedet == 0) his->Fill(fROC[idet]->GetValue(irow,icol)*factor);
383 else his->Fill(fROC[idet]->GetValue(irow,icol)+factor);
391 //_____________________________________________________________________________
392 TH2F *AliTRDCalPad::MakeHisto2DSmPl(Int_t sm, Int_t pl, const AliTRDCalDet *calDet
393 , Int_t typedet, Float_t min, Float_t max,Int_t type)
397 // sm - supermodule number
399 // If calDet correct the CalROC from the detector coefficient
400 // typedet == 0 for gain and vdrift
401 // typedet == 1 for t0
403 gStyle->SetPalette(1);
404 Double_t factor = 0.0;
405 if(typedet == 0) factor = 1.0;
407 Float_t kEpsilon = 0.000000000001;
408 Float_t kEpsilonr = 0.005;
410 AliTRDgeometry *trdGeo = new AliTRDgeometry();
415 Float_t mean = GetMean(calDet,typedet);
417 if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
420 sigma = GetMeanRMS(rms,calDet,typedet);
422 Float_t nsigma = TMath::Abs(min);
424 if(sigma < kEpsilonr) sigma = kEpsilonr;
430 Float_t mean = GetMedian(calDet,typedet);
437 // LTM mean +- nsigma
440 Float_t mean = GetLTM(&sigma,max,calDet,typedet);
442 if(sigma < kEpsilonr) sigma = kEpsilonr;
448 AliTRDpadPlane *padPlane0 = trdGeo->GetPadPlane(pl,0);
449 Double_t row0 = padPlane0->GetRow0();
450 Double_t col0 = padPlane0->GetCol0();
453 snprintf(name,1000,"%s Pad 2D sm %d pl %d",GetTitle(),sm,pl);
454 TH2F * his = new TH2F( name, name, 76,-TMath::Abs(row0),TMath::Abs(row0),144,-TMath::Abs(col0),TMath::Abs(col0));
457 Int_t offsetsmpl = 30*sm+pl;
460 for (Int_t k = 0; k < kNcham; k++){
461 Int_t det = offsetsmpl+k*6;
462 if(calDet) factor = calDet->GetValue(det);
464 AliTRDCalROC * calRoc = fROC[det];
465 for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
466 for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
467 if (TMath::Abs(calRoc->GetValue(icol,irow))>kEpsilon){
469 Int_t kb = kNcham-1-k;
470 Int_t krow = calRoc->GetNrows()-1-irow;
471 Int_t kcol = calRoc->GetNcols()-1-icol;
472 if(kb > 2) binz = 16*(kb-1)+12+krow+1;
473 else binz = 16*kb+krow+1;
475 Float_t value = calRoc->GetValue(icol,irow);
476 if(typedet == 0) his->SetBinContent(binz,biny,value*factor);
477 else his->SetBinContent(binz,biny,value+factor);
483 his->SetXTitle("z (cm)");
484 his->SetYTitle("y (cm)");
486 his->SetMaximum(max);
487 his->SetMinimum(min);
492 //_____________________________________________________________________________
493 TH2F *AliTRDCalPad::MakeHisto2DCh(Int_t ch, const AliTRDCalDet *calDet, Int_t typedet
494 , Float_t min, Float_t max,Int_t type)
497 // Make 2D graph mean value in z direction
499 // If calDet correct the CalROC from the detector coefficient
500 // typedet == 0 for gain and vdrift
501 // typedet == 1 for t0
503 gStyle->SetPalette(1);
504 Double_t factor = 0.0;
505 if(typedet == 0) factor = 1.0;
507 Float_t kEpsilonr = 0.005;
512 Float_t mean = GetMean(calDet,typedet);
514 if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
517 sigma = GetMeanRMS(rms,calDet,typedet);
519 Float_t nsigma = TMath::Abs(min);
521 if(sigma < kEpsilonr) sigma = kEpsilonr;
527 Float_t mean = GetMedian(calDet,typedet);
534 // LTM mean +- nsigma
537 Float_t mean = GetLTM(&sigma,max,calDet,typedet);
539 if(sigma < kEpsilonr) sigma = kEpsilonr;
545 AliTRDgeometry *trdGeo = new AliTRDgeometry();
547 Float_t kEpsilon = 0.000000000001;
549 Double_t poslocal[3] = {0.0,0.0,0.0};
550 Double_t posglobal[3] = {0.0,0.0,0.0};
553 snprintf(name,1000,"%s Pad 2D ch %d",GetTitle(),ch);
554 TH2F * his = new TH2F( name, name, 400,-400.0,400.0,400,-400.0,400.0);
557 Int_t offsetch = 6*ch;
560 for (Int_t isec = 0; isec < kNsect; isec++){
561 for(Int_t ipl = 0; ipl < kNplan; ipl++){
562 Int_t det = offsetch+isec*30+ipl;
563 if(calDet) factor = calDet->GetValue(det);
565 AliTRDCalROC * calRoc = fROC[det];
566 AliTRDpadPlane *padPlane = trdGeo->GetPadPlane(ipl,ch);
567 for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
568 poslocal[0] = trdGeo->GetTime0(ipl);
569 poslocal[2] = padPlane->GetRowPos(0);
570 poslocal[1] = padPlane->GetColPos(icol);
571 trdGeo->RotateBack(det,poslocal,posglobal);
572 Int_t binx = 1+TMath::Nint((posglobal[0]+400.0)*0.5);
573 Int_t biny = 1+TMath::Nint((posglobal[1]+400.0)*0.5);
576 for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
577 if (TMath::Abs(calRoc->GetValue(icol,irow))>kEpsilon){
578 value += calRoc->GetValue(icol,irow);
583 if(typedet == 0) his->SetBinContent(binx,biny,value*factor);
584 else his->SetBinContent(binx,biny,value+factor);
589 his->SetXTitle("x (cm)");
590 his->SetYTitle("y (cm)");
592 his->SetMaximum(max);
593 his->SetMinimum(min);
598 //_____________________________________________________________________________
599 Bool_t AliTRDCalPad::Add(Float_t c1)
602 // add constant for all channels of all ROCs
605 Bool_t result = kTRUE;
606 for (Int_t idet = 0; idet < kNdet; idet++) {
608 if(!fROC[idet]->Add(c1)) result = kFALSE;
614 //_____________________________________________________________________________
615 Bool_t AliTRDCalPad::Multiply(Float_t c1)
618 // multiply constant for all channels of all ROCs
620 Bool_t result = kTRUE;
621 for (Int_t idet = 0; idet < kNdet; idet++) {
623 if(!fROC[idet]->Multiply(c1)) result = kFALSE;
629 //_____________________________________________________________________________
630 Bool_t AliTRDCalPad::Add(const AliTRDCalPad * pad, Double_t c1
631 , const AliTRDCalDet * calDet1, const AliTRDCalDet *calDet2
635 // add calpad channel by channel multiplied by c1 - all ROCs
636 // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
637 // then you have calDet1 and the calPad together
638 // calDet2 and pad together
639 // typedet == 0 for gain and vdrift
640 // typedet == 1 for t0
642 Float_t kEpsilon = 0.000000000001;
644 Double_t factor1 = 0.0;
645 Double_t factor2 = 0.0;
650 Bool_t result = kTRUE;
651 for (Int_t idet = 0; idet < kNdet; idet++) {
652 if(calDet1) factor1 = calDet1->GetValue(idet);
653 if(calDet2) factor2 = calDet2->GetValue(idet);
656 if(TMath::Abs(factor1) > kEpsilon){
657 if(!fROC[idet]->Add(pad->GetCalROC(idet),c1*factor2/factor1)) result = kFALSE;
659 else result = kFALSE;
662 AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
663 if(!croc->Add(factor2)) result = kFALSE;
664 if(!fROC[idet]->Add(croc,c1)) result = kFALSE;
671 //_____________________________________________________________________________
672 Bool_t AliTRDCalPad::Multiply(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1
673 , const AliTRDCalDet *calDet2, Int_t type)
676 // multiply calpad channel by channel - all ROCs
677 // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
678 // then you have calDet1 and the calPad together
679 // calDet2 and pad together
680 // typedet == 0 for gain and vdrift
681 // typedet == 1 for t0
683 Float_t kEpsilon = 0.000000000001;
684 Bool_t result = kTRUE;
685 Double_t factor1 = 0.0;
686 Double_t factor2 = 0.0;
691 for (Int_t idet = 0; idet < kNdet; idet++) {
692 if(calDet1) factor1 = calDet1->GetValue(idet);
693 if(calDet2) factor2 = calDet2->GetValue(idet);
696 if(TMath::Abs(factor1) > kEpsilon){
697 AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
698 if(!croc->Multiply(factor2)) result = kFALSE;
699 if(!fROC[idet]->Multiply(croc)) result = kFALSE;
701 else result = kFALSE;
704 AliTRDCalROC *croc2 = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
705 if(!croc2->Add(factor2)) result = kFALSE;
706 if(!fROC[idet]->Add(factor1)) result = kFALSE;
707 if(!fROC[idet]->Multiply(croc2)) result = kFALSE;
708 if(!fROC[idet]->Add(-factor1)) result = kFALSE;
715 //_____________________________________________________________________________
716 Bool_t AliTRDCalPad::Divide(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1
717 , const AliTRDCalDet *calDet2, Int_t type)
720 // divide calpad channel by channel - all ROCs
721 // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
722 // then you have calDet1 and the calPad together
723 // calDet2 and pad together
724 // typedet == 0 for gain and vdrift
725 // typedet == 1 for t0
727 Float_t kEpsilon = 0.000000000001;
728 Bool_t result = kTRUE;
729 Double_t factor1 = 0.0;
730 Double_t factor2 = 0.0;
735 for (Int_t idet = 0; idet < kNdet; idet++) {
736 if(calDet1) factor1 = calDet1->GetValue(idet);
737 if(calDet2) factor2 = calDet2->GetValue(idet);
740 if(TMath::Abs(factor1) > kEpsilon){
741 AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
742 if(!croc->Multiply(factor2)) result = kFALSE;
743 if(!fROC[idet]->Divide(croc)) result = kFALSE;
745 else result = kFALSE;
748 AliTRDCalROC *croc2 = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
749 if(!croc2->Add(factor2)) result = kFALSE;
750 if(!fROC[idet]->Add(factor1)) result = kFALSE;
751 if(!fROC[idet]->Divide(croc2)) result = kFALSE;
752 if(!fROC[idet]->Add(-factor1)) result = kFALSE;