Use noise information in cluster finder
[u/mrichter/AliRoot.git] / TRD / Cal / AliTRDCalPad.cxx
CommitLineData
7754cd1f 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
15
16/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
20// TRD calibration class for parameters which saved per pad //
21// //
22///////////////////////////////////////////////////////////////////////////////
23
4c87f26e 24#include <TMath.h>
25#include <TH2F.h>
26#include <TH1F.h>
27#include <TStyle.h>
28
7754cd1f 29#include "AliTRDCalPad.h"
30#include "AliTRDCalROC.h"
31#include "AliTRDCalDet.h"
4c87f26e 32#include "AliTRDpadPlane.h"
33#include "AliMathBase.h"
34#include "AliTRDgeometry.h"
7754cd1f 35
36ClassImp(AliTRDCalPad)
37
38//_____________________________________________________________________________
39AliTRDCalPad::AliTRDCalPad():TNamed()
40{
41 //
42 // AliTRDCalPad default constructor
43 //
44
45 for (Int_t idet = 0; idet < kNdet; idet++) {
46 fROC[idet] = 0;
47 }
48
49}
50
51//_____________________________________________________________________________
52AliTRDCalPad::AliTRDCalPad(const Text_t *name, const Text_t *title)
53 :TNamed(name,title)
54{
55 //
56 // AliTRDCalPad constructor
57 //
58
59 for (Int_t isec = 0; isec < kNsect; isec++) {
60 for (Int_t ipla = 0; ipla < kNplan; ipla++) {
61 for (Int_t icha = 0; icha < kNcham; icha++) {
62 Int_t idet = GetDet(ipla,icha,isec);
63 fROC[idet] = new AliTRDCalROC(ipla,icha);
64 }
65 }
66 }
67
68}
69
7754cd1f 70//_____________________________________________________________________________
4c87f26e 71AliTRDCalPad::AliTRDCalPad(const AliTRDCalPad &c)
72 :TNamed(c)
7754cd1f 73{
74 //
75 // AliTRDCalPad copy constructor
76 //
77
4c87f26e 78 for (Int_t idet = 0; idet < kNdet; idet++) {
79 fROC[idet] = new AliTRDCalROC(*((AliTRDCalPad &) c).fROC[idet]);
80 }
7754cd1f 81
82}
83
4c87f26e 84//_____________________________________________________________________________
7754cd1f 85AliTRDCalPad::~AliTRDCalPad()
86{
87 //
88 // AliTRDCalPad destructor
89 //
90
91 for (Int_t idet = 0; idet < kNdet; idet++) {
92 if (fROC[idet]) {
93 delete fROC[idet];
94 fROC[idet] = 0;
95 }
96 }
97
98}
99
100//_____________________________________________________________________________
101AliTRDCalPad &AliTRDCalPad::operator=(const AliTRDCalPad &c)
102{
103 //
104 // Assignment operator
105 //
106
107 if (this != &c) ((AliTRDCalPad &) c).Copy(*this);
108 return *this;
109
110}
111
112//_____________________________________________________________________________
113void AliTRDCalPad::Copy(TObject &c) const
114{
115 //
116 // Copy function
117 //
118
119 for (Int_t idet = 0; idet < kNdet; idet++) {
ea2a3ed8 120 if (((AliTRDCalPad &) c).fROC[idet]) {
121 delete ((AliTRDCalPad &) c).fROC[idet];
122 }
123 ((AliTRDCalPad &) c).fROC[idet] = new AliTRDCalROC();
7754cd1f 124 if (fROC[idet]) {
125 fROC[idet]->Copy(*((AliTRDCalPad &) c).fROC[idet]);
126 }
127 }
128
129 TObject::Copy(c);
2745a409 130
7754cd1f 131}
132
133//_____________________________________________________________________________
4c87f26e 134Bool_t AliTRDCalPad::ScaleROCs(const AliTRDCalDet* values)
7754cd1f 135{
136 //
137 // Scales ROCs of this class with the values from the class <values>
138 // Is used if an AliTRDCalPad object defines local variations of a parameter
139 // defined per detector using a AliTRDCalDet class
140 //
141
142 if (!values)
4c87f26e 143 return kFALSE;
144 Bool_t result = kTRUE;
7754cd1f 145 for (Int_t idet = 0; idet < kNdet; idet++) {
146 if (fROC[idet]) {
4c87f26e 147 if(!fROC[idet]->Multiply(values->GetValue(idet))) result = kFALSE;
148 }
149 }
150 return result;
151}
152
153//_____________________________________________________________________________
154Double_t AliTRDCalPad::GetMeanRMS(Double_t &rms, const AliTRDCalDet *calDet, Int_t type)
155{
156 //
157 // Calculate mean an RMS of all rocs
158 // If calDet correct the CalROC from the detector coefficient
159 // type == 0 for gain and vdrift
160 // type == 1 for t0
161 //
162 Double_t factor = 0.0;
163 if(type == 0) factor = 1.0;
164 Double_t sum = 0, sum2 = 0, n=0, val=0;
165 for (Int_t idet = 0; idet < kNdet; idet++) {
166 if(calDet) factor = calDet->GetValue(idet);
167 AliTRDCalROC *calRoc = fROC[idet];
168 if ( calRoc ){
169 for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
170 for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
171 if(type == 0) val = calRoc->GetValue(icol,irow)*factor;
172 else val = calRoc->GetValue(icol,irow)+factor;
173 sum+=val;
174 sum2+=val*val;
175 n++;
176 }
177 }
178 }
179 }
180 Double_t n1 = 1./n;
181 Double_t mean = sum*n1;
182 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
183 return mean;
184}
185
186//_____________________________________________________________________________
187Double_t AliTRDCalPad::GetMean(const AliTRDCalDet *calDet, Int_t type, AliTRDCalPad* outlierPad)
188{
189 //
190 // return mean of the mean of all ROCs
191 // If calDet correct the CalROC from the detector coefficient
192 // type == 0 for gain and vdrift
193 // type == 1 for t0
194 //
195 Double_t factor = 0.0;
196 if(type == 0) factor = 1.0;
197 Double_t arr[kNdet];
198 Int_t n=0;
199 for (Int_t idet = 0; idet < kNdet; idet++) {
200 if(calDet) factor = calDet->GetValue(idet);
201 AliTRDCalROC *calRoc = fROC[idet];
202 if ( calRoc ){
203 AliTRDCalROC* outlierROC = 0;
204 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
205 if(type == 0) arr[n] = calRoc->GetMean(outlierROC)*factor;
206 else arr[n] = calRoc->GetMean(outlierROC)+factor;
207 n++;
208 }
209 }
210 return TMath::Mean(n,arr);
211}
212
213//_____________________________________________________________________________
214Double_t AliTRDCalPad::GetRMS(const AliTRDCalDet *calDet, Int_t type, AliTRDCalPad* outlierPad)
215{
216 //
217 // return mean of the RMS of all ROCs
218 // If calDet correct the CalROC from the detector coefficient
219 // type == 0 for gain and vdrift
220 // type == 1 for t0
221 //
222 Double_t factor = 0.0;
223 if(type == 0) factor = 1.0;
224 Double_t arr[kNdet];
225 Int_t n=0;
226 for (Int_t idet = 0; idet < kNdet; idet++) {
227 if(calDet) factor = calDet->GetValue(idet);
228 AliTRDCalROC *calRoc = fROC[idet];
229 if ( calRoc ){
230 AliTRDCalROC* outlierROC = 0;
231 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
232 if(type == 0) arr[n] = calRoc->GetRMS(outlierROC)*factor;
233 else arr[n] = calRoc->GetRMS(outlierROC);
234 n++;
235 }
236 }
237 return TMath::Mean(n,arr);
238}
239
240//_____________________________________________________________________________
241Double_t AliTRDCalPad::GetMedian(const AliTRDCalDet *calDet, Int_t type
242 , AliTRDCalPad* outlierPad)
243{
244 //
245 // return mean of the median of all ROCs
246 // If AliTRDCalDet, the correct the CalROC from the detector coefficient
247 // type == 0 for gain and vdrift
248 // type == 1 for t0
249 //
250 Double_t factor = 0.0;
251 if(type == 0) factor = 1.0;
252 Double_t arr[kNdet];
253 Int_t n=0;
254 for (Int_t idet = 0; idet < kNdet; idet++) {
255 if(calDet) factor = calDet->GetValue(idet);
256 AliTRDCalROC *calRoc = fROC[idet];
257 if ( calRoc ){
258 AliTRDCalROC* outlierROC = 0;
259 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
260 if(type == 0) arr[n] = calRoc->GetMedian(outlierROC)*factor;
261 else arr[n] = calRoc->GetMedian(outlierROC)+factor;
262 n++;
263 }
264 }
265 return TMath::Mean(n,arr);
266}
267
268//_____________________________________________________________________________
269Double_t AliTRDCalPad::GetLTM(Double_t *sigma, Double_t fraction
270 , const AliTRDCalDet *calDet, Int_t type
271 , AliTRDCalPad* outlierPad)
272{
273 //
274 // return mean of the LTM and sigma of all ROCs
275 // If calDet correct the CalROC from the detector coefficient
276 // type == 0 for gain and vdrift
277 // type == 1 for t0
278 //
279 Double_t factor = 0.0;
280 if(type == 0) factor = 1.0;
281 Double_t arrm[kNdet];
282 Double_t arrs[kNdet];
283 Double_t *sTemp=0x0;
284 Int_t n=0;
285
286 for (Int_t idet = 0; idet < kNdet; idet++) {
287 if(calDet) factor = calDet->GetValue(idet);
288 AliTRDCalROC *calRoc = fROC[idet];
289 if ( calRoc ){
290 if ( sigma ) sTemp=arrs+n;
291 AliTRDCalROC* outlierROC = 0;
292 if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
293 if(type == 0) arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC)*factor;
294 else arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC)+factor;
295 n++;
296 }
297 }
298 if ( sigma ) *sigma = TMath::Mean(n,arrs);
299 return TMath::Mean(n,arrm);
300}
301
302//_____________________________________________________________________________
303TH1F * AliTRDCalPad::MakeHisto1D(const AliTRDCalDet *calDet, Int_t typedet
304 , Float_t min, Float_t max,Int_t type)
305{
306 //
307 // make 1D histo
308 // type -1 = user defined range
309 // 0 = nsigma cut nsigma=min
310 // If calDet correct the CalROC from the detector coefficient
311 // typedet == 0 for gain and vdrift
312 // typedet == 1 for t0
313 //
314 Float_t kEpsilonr = 0.005;
315
316 Double_t factor = 0.0;
317 if(typedet == 0) factor = 1.0;
318
319 if (type>=0){
320 if (type==0){
321 // nsigma range
322 Float_t mean = GetMean(calDet,typedet);
323 Float_t sigma = 0.0;
324 if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
325 else {
326 Double_t rms = 0.0;
327 sigma = GetMeanRMS(rms,calDet,typedet);
328 }
329 Float_t nsigma = TMath::Abs(min);
330 sigma *= nsigma;
331 if(sigma < kEpsilonr) sigma = kEpsilonr;
332 min = mean-sigma;
333 max = mean+sigma;
334 }
335 if (type==1){
336 // fixed range
337 Float_t mean = GetMedian(calDet,typedet);
338 Float_t delta = min;
339 min = mean-delta;
340 max = mean+delta;
341 }
342 if (type==2){
343 //
344 // LTM mean +- nsigma
345 //
346 Double_t sigma;
347 Float_t mean = GetLTM(&sigma,max,calDet,typedet);
348 sigma*=min;
349 if(sigma < kEpsilonr) sigma = kEpsilonr;
350 min = mean-sigma;
351 max = mean+sigma;
352 }
353 }
354 char name[1000];
355 sprintf(name,"%s Pad 1D",GetTitle());
356 TH1F * his = new TH1F(name,name,100, min,max);
357 for (Int_t idet = 0; idet < kNdet; idet++) {
358 if(calDet) factor = calDet->GetValue(idet);
359 if (fROC[idet]){
360 for (Int_t irow=0; irow<fROC[idet]->GetNrows(); irow++){
361 for (Int_t icol=0; icol<fROC[idet]->GetNcols(); icol++){
362 if(typedet == 0) his->Fill(fROC[idet]->GetValue(irow,icol)*factor);
363 else his->Fill(fROC[idet]->GetValue(irow,icol)+factor);
364 }
365 }
366 }
367 }
368 return his;
369}
370
371//_____________________________________________________________________________
372TH2F *AliTRDCalPad::MakeHisto2DSmPl(Int_t sm, Int_t pl, const AliTRDCalDet *calDet
373 , Int_t typedet, Float_t min, Float_t max,Int_t type)
374{
375 //
376 // Make 2D graph
377 // sm - supermodule number
378 // pl - plane number
379 // If calDet correct the CalROC from the detector coefficient
380 // typedet == 0 for gain and vdrift
381 // typedet == 1 for t0
382 //
383 gStyle->SetPalette(1);
384 Double_t factor = 0.0;
385 if(typedet == 0) factor = 1.0;
386
387 Float_t kEpsilon = 0.000000000001;
388 Float_t kEpsilonr = 0.005;
389
015cd5ba 390 AliTRDgeometry *trdGeo = new AliTRDgeometry();
391
4c87f26e 392 if (type>=0){
393 if (type==0){
394 // nsigma range
395 Float_t mean = GetMean(calDet,typedet);
396 Float_t sigma = 0.0;
397 if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
398 else {
399 Double_t rms = 0.0;
400 sigma = GetMeanRMS(rms,calDet,typedet);
401 }
402 Float_t nsigma = TMath::Abs(min);
403 sigma *= nsigma;
404 if(sigma < kEpsilonr) sigma = kEpsilonr;
405 min = mean-sigma;
406 max = mean+sigma;
407 }
408 if (type==1){
409 // fixed range
410 Float_t mean = GetMedian(calDet,typedet);
411 Float_t delta = min;
412 min = mean-delta;
413 max = mean+delta;
414 }
415 if (type==2){
416 //
417 // LTM mean +- nsigma
418 //
419 Double_t sigma;
420 Float_t mean = GetLTM(&sigma,max,calDet,typedet);
421 sigma*=min;
422 if(sigma < kEpsilonr) sigma = kEpsilonr;
423 min = mean-sigma;
424 max = mean+sigma;
425 }
426 }
427
015cd5ba 428 AliTRDpadPlane *padPlane0 = trdGeo->GetPadPlane(pl,0);
4c87f26e 429 Double_t row0 = padPlane0->GetRow0();
430 Double_t col0 = padPlane0->GetCol0();
431
432 char name[1000];
433 sprintf(name,"%s Pad 2D sm %d pl %d",GetTitle(),sm,pl);
434 TH2F * his = new TH2F( name, name, 76,-TMath::Abs(row0),TMath::Abs(row0),144,-TMath::Abs(col0),TMath::Abs(col0));
435
436 // Where we begin
437 Int_t offsetsmpl = 30*sm+pl;
438
439
440 for (Int_t k = 0; k < kNcham; k++){
441 Int_t det = offsetsmpl+k*6;
442 if(calDet) factor = calDet->GetValue(det);
443 if (fROC[det]){
444 AliTRDCalROC * calRoc = fROC[det];
445 for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
446 for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
447 if (TMath::Abs(calRoc->GetValue(icol,irow))>kEpsilon){
448 Int_t binz = 0;
449 Int_t kb = kNcham-1-k;
413153cb 450 Int_t krow = calRoc->GetNrows()-1-irow;
451 Int_t kcol = calRoc->GetNcols()-1-icol;
452 if(kb > 2) binz = 16*(kb-1)+12+krow+1;
453 else binz = 16*kb+krow+1;
454 Int_t biny = kcol+1;
4c87f26e 455 Float_t value = calRoc->GetValue(icol,irow);
456 if(typedet == 0) his->SetBinContent(binz,biny,value*factor);
457 else his->SetBinContent(binz,biny,value+factor);
458 }
459 }
460 }
461 }
462 }
463 his->SetXTitle("z (cm)");
464 his->SetYTitle("y (cm)");
465 his->SetStats(0);
466 his->SetMaximum(max);
467 his->SetMinimum(min);
015cd5ba 468 delete trdGeo;
4c87f26e 469 return his;
470}
471
472//_____________________________________________________________________________
473TH2F *AliTRDCalPad::MakeHisto2DCh(Int_t ch, const AliTRDCalDet *calDet, Int_t typedet
474 , Float_t min, Float_t max,Int_t type)
475{
476 //
477 // Make 2D graph mean value in z direction
478 // ch - chamber
479 // If calDet correct the CalROC from the detector coefficient
480 // typedet == 0 for gain and vdrift
481 // typedet == 1 for t0
482 //
483 gStyle->SetPalette(1);
484 Double_t factor = 0.0;
485 if(typedet == 0) factor = 1.0;
486
487 Float_t kEpsilonr = 0.005;
488
489 if (type>=0){
490 if (type==0){
491 // nsigma range
492 Float_t mean = GetMean(calDet,typedet);
493 Float_t sigma = 0.0;
494 if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
495 else {
496 Double_t rms = 0.0;
497 sigma = GetMeanRMS(rms,calDet,typedet);
498 }
499 Float_t nsigma = TMath::Abs(min);
500 sigma *= nsigma;
501 if(sigma < kEpsilonr) sigma = kEpsilonr;
502 min = mean-sigma;
503 max = mean+sigma;
504 }
505 if (type==1){
506 // fixed range
507 Float_t mean = GetMedian(calDet,typedet);
508 Float_t delta = min;
509 min = mean-delta;
510 max = mean+delta;
511 }
512 if (type==2){
513 //
514 // LTM mean +- nsigma
515 //
516 Double_t sigma;
517 Float_t mean = GetLTM(&sigma,max,calDet,typedet);
518 sigma*=min;
519 if(sigma < kEpsilonr) sigma = kEpsilonr;
520 min = mean-sigma;
521 max = mean+sigma;
522 }
523 }
524
015cd5ba 525 AliTRDgeometry *trdGeo = new AliTRDgeometry();
4c87f26e 526
527 Float_t kEpsilon = 0.000000000001;
528
529 Double_t poslocal[3] = {0.0,0.0,0.0};
530 Double_t posglobal[3] = {0.0,0.0,0.0};
531
532 char name[1000];
533 sprintf(name,"%s Pad 2D ch %d",GetTitle(),ch);
534 TH2F * his = new TH2F( name, name, 400,-400.0,400.0,400,-400.0,400.0);
535
536 // Where we begin
537 Int_t offsetch = 6*ch;
538
539
540 for (Int_t isec = 0; isec < kNsect; isec++){
541 for(Int_t ipl = 0; ipl < kNplan; ipl++){
542 Int_t det = offsetch+isec*30+ipl;
543 if(calDet) factor = calDet->GetValue(det);
544 if (fROC[det]){
545 AliTRDCalROC * calRoc = fROC[det];
015cd5ba 546 AliTRDpadPlane *padPlane = trdGeo->GetPadPlane(ipl,ch);
4c87f26e 547 for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
015cd5ba 548 poslocal[0] = trdGeo->GetTime0(ipl);
4c87f26e 549 poslocal[2] = padPlane->GetRowPos(0);
550 poslocal[1] = padPlane->GetColPos(icol);
015cd5ba 551 trdGeo->RotateBack(det,poslocal,posglobal);
4c87f26e 552 Int_t binx = 1+TMath::Nint((posglobal[0]+400.0)*0.5);
553 Int_t biny = 1+TMath::Nint((posglobal[1]+400.0)*0.5);
554 Float_t value = 0.0;
555 Int_t nb = 0;
556 for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
557 if (TMath::Abs(calRoc->GetValue(icol,irow))>kEpsilon){
558 value += calRoc->GetValue(icol,irow);
559 nb++;
560 }
561 }
562 value = value/nb;
563 if(typedet == 0) his->SetBinContent(binx,biny,value*factor);
564 else his->SetBinContent(binx,biny,value+factor);
565 }
566 }
567 }
568 }
569 his->SetXTitle("x (cm)");
570 his->SetYTitle("y (cm)");
571 his->SetStats(0);
572 his->SetMaximum(max);
573 his->SetMinimum(min);
015cd5ba 574 delete trdGeo;
4c87f26e 575 return his;
576}
577
578//_____________________________________________________________________________
579Bool_t AliTRDCalPad::Add(Float_t c1)
580{
581 //
582 // add constant for all channels of all ROCs
583 //
584
585 Bool_t result = kTRUE;
586 for (Int_t idet = 0; idet < kNdet; idet++) {
587 if (fROC[idet]){
588 if(!fROC[idet]->Add(c1)) result = kFALSE;
7754cd1f 589 }
590 }
4c87f26e 591 return result;
592}
593
594//_____________________________________________________________________________
595Bool_t AliTRDCalPad::Multiply(Float_t c1)
596{
597 //
598 // multiply constant for all channels of all ROCs
599 //
600 Bool_t result = kTRUE;
601 for (Int_t idet = 0; idet < kNdet; idet++) {
602 if (fROC[idet]){
603 if(!fROC[idet]->Multiply(c1)) result = kFALSE;
604 }
605 }
606 return result;
607}
2745a409 608
4c87f26e 609//_____________________________________________________________________________
610Bool_t AliTRDCalPad::Add(const AliTRDCalPad * pad, Double_t c1
611 , const AliTRDCalDet * calDet1, const AliTRDCalDet *calDet2
612 , Int_t type)
613{
614 //
615 // add calpad channel by channel multiplied by c1 - all ROCs
616 // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
617 // then you have calDet1 and the calPad together
618 // calDet2 and pad together
619 // typedet == 0 for gain and vdrift
620 // typedet == 1 for t0
621 //
622 Float_t kEpsilon = 0.000000000001;
623
624 Double_t factor1 = 0.0;
625 Double_t factor2 = 0.0;
626 if(type == 0) {
627 factor1 = 1.0;
628 factor2 = 1.0;
629 }
630 Bool_t result = kTRUE;
631 for (Int_t idet = 0; idet < kNdet; idet++) {
632 if(calDet1) factor1 = calDet1->GetValue(idet);
633 if(calDet2) factor2 = calDet2->GetValue(idet);
634 if (fROC[idet]){
635 if(type == 0){
636 if(TMath::Abs(factor1) > kEpsilon){
637 if(!fROC[idet]->Add(pad->GetCalROC(idet),c1*factor2/factor1)) result = kFALSE;
638 }
639 else result = kFALSE;
640 }
641 else{
642 AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
643 if(!croc->Add(factor2)) result = kFALSE;
644 if(!fROC[idet]->Add(croc,c1)) result = kFALSE;
645 }
646 }
647 }
648 return result;
7754cd1f 649}
650
4c87f26e 651//_____________________________________________________________________________
652Bool_t AliTRDCalPad::Multiply(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1
653 , const AliTRDCalDet *calDet2, Int_t type)
654{
655 //
656 // multiply calpad channel by channel - all ROCs
657 // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
658 // then you have calDet1 and the calPad together
659 // calDet2 and pad together
660 // typedet == 0 for gain and vdrift
661 // typedet == 1 for t0
662 //
663 Float_t kEpsilon = 0.000000000001;
664 Bool_t result = kTRUE;
665 Double_t factor1 = 0.0;
666 Double_t factor2 = 0.0;
667 if(type == 0){
668 factor1 = 1.0;
669 factor2 = 1.0;
670 }
671 for (Int_t idet = 0; idet < kNdet; idet++) {
672 if(calDet1) factor1 = calDet1->GetValue(idet);
673 if(calDet2) factor2 = calDet2->GetValue(idet);
674 if (fROC[idet]){
675 if(type == 0){
676 if(TMath::Abs(factor1) > kEpsilon){
677 AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
678 if(!croc->Multiply(factor2)) result = kFALSE;
679 if(!fROC[idet]->Multiply(croc)) result = kFALSE;
680 }
681 else result = kFALSE;
682 }
683 else{
684 AliTRDCalROC *croc2 = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
685 if(!croc2->Add(factor2)) result = kFALSE;
686 if(!fROC[idet]->Add(factor1)) result = kFALSE;
687 if(!fROC[idet]->Multiply(croc2)) result = kFALSE;
688 if(!fROC[idet]->Add(-factor1)) result = kFALSE;
689 }
690 }
691 }
692 return result;
693}
694
695//_____________________________________________________________________________
696Bool_t AliTRDCalPad::Divide(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1
697 , const AliTRDCalDet *calDet2, Int_t type)
698{
699 //
700 // divide calpad channel by channel - all ROCs
701 // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
702 // then you have calDet1 and the calPad together
703 // calDet2 and pad together
704 // typedet == 0 for gain and vdrift
705 // typedet == 1 for t0
706 //
707 Float_t kEpsilon = 0.000000000001;
708 Bool_t result = kTRUE;
709 Double_t factor1 = 0.0;
710 Double_t factor2 = 0.0;
711 if(type == 0){
712 factor1 = 1.0;
713 factor2 = 1.0;
714 }
715 for (Int_t idet = 0; idet < kNdet; idet++) {
716 if(calDet1) factor1 = calDet1->GetValue(idet);
717 if(calDet2) factor2 = calDet2->GetValue(idet);
718 if (fROC[idet]){
719 if(type == 0){
720 if(TMath::Abs(factor1) > kEpsilon){
721 AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
722 if(!croc->Multiply(factor2)) result = kFALSE;
723 if(!fROC[idet]->Divide(croc)) result = kFALSE;
724 }
725 else result = kFALSE;
726 }
727 else{
728 AliTRDCalROC *croc2 = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
729 if(!croc2->Add(factor2)) result = kFALSE;
730 if(!fROC[idet]->Add(factor1)) result = kFALSE;
731 if(!fROC[idet]->Divide(croc2)) result = kFALSE;
732 if(!fROC[idet]->Add(-factor1)) result = kFALSE;
733 }
734 }
735 }
736 return result;
737}