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 detector //
22 ///////////////////////////////////////////////////////////////////////////////
29 #include "AliTRDCalDet.h"
30 #include "AliTRDgeometry.h"
31 #include "AliMathBase.h"
32 #include "AliTRDpadPlane.h"
34 ClassImp(AliTRDCalDet)
36 //_____________________________________________________________________________
37 AliTRDCalDet::AliTRDCalDet():TNamed()
40 // AliTRDCalDet default constructor
43 for (Int_t idet = 0; idet < kNdet; idet++) {
49 //_____________________________________________________________________________
50 AliTRDCalDet::AliTRDCalDet(const Text_t *name, const Text_t *title)
54 // AliTRDCalDet constructor
57 for (Int_t idet = 0; idet < kNdet; idet++) {
63 //_____________________________________________________________________________
64 AliTRDCalDet::AliTRDCalDet(const AliTRDCalDet &c):TNamed(c)
67 // AliTRDCalDet copy constructor
70 ((AliTRDCalDet &) c).Copy(*this);
74 ///_____________________________________________________________________________
75 AliTRDCalDet::~AliTRDCalDet()
78 // AliTRDCalDet destructor
83 //_____________________________________________________________________________
84 AliTRDCalDet &AliTRDCalDet::operator=(const AliTRDCalDet &c)
87 // Assignment operator
90 if (this != &c) ((AliTRDCalDet &) c).Copy(*this);
95 //_____________________________________________________________________________
96 void AliTRDCalDet::Copy(TObject &c) const
102 for (Int_t idet = 0; idet < kNdet; idet++) {
103 ((AliTRDCalDet &) c).fData[idet] = fData[idet];
110 //___________________________________________________________________________________
111 Double_t AliTRDCalDet::GetMean(AliTRDCalDet * const outlierDet) const
114 // Calculate the mean
117 if (!outlierDet) return TMath::Mean(kNdet,fData);
118 Double_t *ddata = new Double_t[kNdet];
120 for (Int_t i=0;i<kNdet;i++) {
121 if (!(outlierDet->GetValue(i))) {
122 ddata[nPoints]= fData[nPoints];
126 Double_t mean = TMath::Mean(nPoints,ddata);
131 //_______________________________________________________________________________________
132 Double_t AliTRDCalDet::GetMedian(AliTRDCalDet * const outlierDet) const
135 // Calculate the median
138 if (!outlierDet) return (Double_t) TMath::Median(kNdet,fData);
139 Double_t *ddata = new Double_t[kNdet];
141 for (Int_t i=0;i<kNdet;i++) {
142 if (!(outlierDet->GetValue(i))) {
143 ddata[nPoints]= fData[nPoints];
147 Double_t mean = TMath::Median(nPoints,ddata);
153 //____________________________________________________________________________________________
154 Double_t AliTRDCalDet::GetRMS(AliTRDCalDet * const outlierDet) const
160 if (!outlierDet) return TMath::RMS(kNdet,fData);
161 Double_t *ddata = new Double_t[kNdet];
163 for (Int_t i=0;i<kNdet;i++) {
164 if (!(outlierDet->GetValue(i))) {
165 ddata[nPoints]= fData[nPoints];
169 Double_t mean = TMath::RMS(nPoints,ddata);
174 //______________________________________________________________________________________________
175 Double_t AliTRDCalDet::GetLTM(Double_t *sigma
177 , AliTRDCalDet * const outlierDet)
180 // Calculate LTM mean and sigma
183 Double_t *ddata = new Double_t[kNdet];
184 Double_t mean=0, lsigma=0;
186 for (Int_t i=0;i<kNdet;i++) {
187 if (!outlierDet || !(outlierDet->GetValue(i))) {
188 ddata[nPoints]= fData[nPoints];
192 Int_t hh = TMath::Min(TMath::Nint(fraction *nPoints), Int_t(nPoints));
193 AliMathBase::EvaluateUni(nPoints,ddata, mean, lsigma, hh);
194 if (sigma) *sigma=lsigma;
199 //_________________________________________________________________________
200 TH1F * AliTRDCalDet::MakeHisto1Distribution(Float_t min, Float_t max,Int_t type)
204 // type -1 = user defined range
205 // 0 = nsigma cut nsigma=min
206 // 1 = delta cut around median delta=min
212 Float_t mean = GetMean();
213 Float_t sigma = GetRMS();
214 Float_t nsigma = TMath::Abs(min);
215 min = mean-nsigma*sigma;
216 max = mean+nsigma*sigma;
220 Float_t mean = GetMedian();
227 // LTM mean +- nsigma
230 Float_t mean = GetLTM(&sigma,max);
237 sprintf(name,"%s CalDet 1Distribution",GetTitle());
238 TH1F * his = new TH1F(name,name,100, min,max);
239 for (Int_t idet=0; idet<kNdet; idet++){
240 his->Fill(GetValue(idet));
245 //________________________________________________________________________________
246 TH1F * AliTRDCalDet::MakeHisto1DAsFunctionOfDet(Float_t min, Float_t max,Int_t type)
250 // type -1 = user defined range
251 // 0 = nsigma cut nsigma=min
252 // 1 = delta cut around median delta=min
258 Float_t mean = GetMean();
259 Float_t sigma = GetRMS();
260 Float_t nsigma = TMath::Abs(min);
261 min = mean-nsigma*sigma;
262 max = mean+nsigma*sigma;
266 Float_t mean = GetMedian();
273 Float_t mean = GetLTM(&sigma,max);
282 sprintf(name,"%s CalDet as function of det",GetTitle());
283 TH1F * his = new TH1F(name,name,kNdet, 0, kNdet);
284 for(Int_t det = 0; det< kNdet; det++){
285 his->Fill(det+0.5,GetValue(det));
288 his->SetMaximum(max);
289 his->SetMinimum(min);
293 //_____________________________________________________________________________
294 TH2F *AliTRDCalDet::MakeHisto2DCh(Int_t ch, Float_t min, Float_t max, Int_t type)
299 // type -1 = user defined range
300 // 0 = nsigma cut nsigma=min
301 // 1 = delta cut around median delta=min
304 gStyle->SetPalette(1);
308 Float_t mean = GetMean();
309 Float_t sigma = GetRMS();
310 Float_t nsigma = TMath::Abs(min);
311 min = mean-nsigma*sigma;
312 max = mean+nsigma*sigma;
316 Float_t mean = GetMedian();
323 Float_t mean = GetLTM(&sigma,max);
331 AliTRDgeometry *trdGeo = new AliTRDgeometry();
333 Double_t poslocal[3] = {0.0,0.0,0.0};
334 Double_t posglobal[3] = {0.0,0.0,0.0};
337 sprintf(name,"%s CalDet 2D ch %d",GetTitle(),ch);
338 TH2F * his = new TH2F(name, name, 400,-400.0,400.0,400,-400.0,400.0);
341 Int_t offsetch = 6*ch;
344 for (Int_t isec = 0; isec < kNsect; isec++){
345 for(Int_t ipl = 0; ipl < kNplan; ipl++){
346 Int_t det = offsetch+isec*30+ipl;
347 AliTRDpadPlane *padPlane = trdGeo->GetPadPlane(ipl,ch);
348 for (Int_t icol=0; icol<padPlane->GetNcols(); icol++){
349 poslocal[0] = trdGeo->GetTime0(ipl);
350 poslocal[2] = padPlane->GetRowPos(0);
351 poslocal[1] = padPlane->GetColPos(icol);
352 trdGeo->RotateBack(det,poslocal,posglobal);
353 Int_t binx = 1+TMath::Nint((posglobal[0]+400.0)*0.5);
354 Int_t biny = 1+TMath::Nint((posglobal[1]+400.0)*0.5);
355 his->SetBinContent(binx,biny,fData[det]);
359 his->SetXTitle("x (cm)");
360 his->SetYTitle("y (cm)");
362 his->SetMaximum(max);
363 his->SetMinimum(min);
368 //_____________________________________________________________________________
369 TH2F *AliTRDCalDet::MakeHisto2DSmPl(Int_t sm, Int_t pl, Float_t min, Float_t max, Int_t type)
373 // sm - supermodule number
375 // type -1 = user defined range
376 // 0 = nsigma cut nsigma=min
377 // 1 = delta cut around median delta=min
380 gStyle->SetPalette(1);
384 Float_t mean = GetMean();
385 Float_t sigma = GetRMS();
386 Float_t nsigma = TMath::Abs(min);
387 min = mean-nsigma*sigma;
388 max = mean+nsigma*sigma;
392 Float_t mean = GetMedian();
399 Float_t mean = GetLTM(&sigma,max);
407 AliTRDgeometry *trdGeo = new AliTRDgeometry();
408 AliTRDpadPlane *padPlane0 = trdGeo->GetPadPlane(pl,0);
409 Double_t row0 = padPlane0->GetRow0();
410 Double_t col0 = padPlane0->GetCol0();
413 sprintf(name,"%s CalDet 2D sm %d and pl %d",GetTitle(),sm,pl);
414 TH2F * his = new TH2F( name, name, 5, -TMath::Abs(row0), TMath::Abs(row0)
415 , 4,-2*TMath::Abs(col0),2*TMath::Abs(col0));
418 Int_t offsetsmpl = 30*sm+pl;
421 for (Int_t k = 0; k < kNcham; k++){
422 Int_t det = offsetsmpl+k*6;
423 Int_t kb = kNcham-1-k;
424 his->SetBinContent(kb+1,2,fData[det]);
425 his->SetBinContent(kb+1,3,fData[det]);
427 his->SetXTitle("z (cm)");
428 his->SetYTitle("y (cm)");
430 his->SetMaximum(max);
431 his->SetMinimum(min);
436 //_____________________________________________________________________________
437 void AliTRDCalDet::Add(Float_t c1)
440 // Add constant for all detectors
443 for (Int_t idet = 0; idet < kNdet; idet++) {
448 //_____________________________________________________________________________
449 void AliTRDCalDet::Multiply(Float_t c1)
452 // multiply constant for all detectors
454 for (Int_t idet = 0; idet < kNdet; idet++) {
459 //_____________________________________________________________________________
460 void AliTRDCalDet::Add(const AliTRDCalDet * calDet, Double_t c1)
463 // add caldet channel by channel multiplied by c1
465 for (Int_t idet = 0; idet < kNdet; idet++) {
466 fData[idet] += calDet->GetValue(idet)*c1;
470 //_____________________________________________________________________________
471 void AliTRDCalDet::Multiply(const AliTRDCalDet * calDet)
474 // multiply caldet channel by channel
476 for (Int_t idet = 0; idet < kNdet; idet++) {
477 fData[idet] *= calDet->GetValue(idet);
481 //_____________________________________________________________________________
482 void AliTRDCalDet::Divide(const AliTRDCalDet * calDet)
485 // divide caldet channel by channel
487 Float_t kEpsilon=0.00000000000000001;
488 for (Int_t idet = 0; idet < kNdet; idet++) {
489 if (TMath::Abs(calDet->GetValue(idet))>kEpsilon){
490 fData[idet] /= calDet->GetValue(idet);