2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
17 /////////////////////////////////////////////////////////////////////////////////
21 // Task to run the calibration offline.
23 // R. Bailhache (rbailhache@ikf.uni-frankfurt.de, R.Bailhache@gsi.de)
25 //////////////////////////////////////////////////////////////////////////////////
28 #include "Riostream.h"
31 #include "TProfile2D.h"
39 #include "TObjArray.h"
43 #include "TGraphErrors.h"
45 #include "AliTRDrecoTask.h"
46 #include "AliAnalysisManager.h"
48 #include "AliESDInputHandler.h"
49 #include "AliTRDtrackV1.h"
50 #include "AliTRDseedV1.h"
51 #include "AliTRDcluster.h"
52 #include "info/AliTRDtrackInfo.h"
53 #include "AliTRDcalibDB.h"
55 #include "AliTRDCalibraFillHisto.h"
56 #include "AliTRDCalibraFit.h"
57 #include "AliTRDCalibraVdriftLinearFit.h"
58 #include "AliTRDCalibraMode.h"
59 #include "AliTRDCalibraVector.h"
60 #include "./Cal/AliTRDCalPad.h"
61 #include "./Cal/AliTRDCalDet.h"
65 #include "AliTRDcalibration.h"
68 ClassImp(AliTRDcalibration)
70 //________________________________________________________________________
71 AliTRDcalibration::AliTRDcalibration()
72 :AliTRDrecoTask("Calibration", "Calibration on tracks")
76 ,fTRDCalibraFillHisto(0)
78 ,fNbTRDTrackOffline(0)
79 ,fNbTRDTrackStandalone(0)
81 ,fNbTRDTrackletOffline(0)
82 ,fNbTRDTrackletStandalone(0)
84 ,fNbTimeBinOffline(0x0)
85 ,fNbTimeBinStandalone(0x0)
87 ,fNbClustersOffline(0)
88 ,fNbClustersStandalone(0)
102 ,fnormalizeNbOfCluster(kFALSE)
104 ,fOfflineTracks(kFALSE)
105 ,fStandaloneTracks(kFALSE)
106 ,fCompressPerDetector(kFALSE)
109 ,fPostProcess(kFALSE)
115 for(Int_t k = 0; k < 3; k++)
122 //________________________________________________________________________
123 AliTRDcalibration::~AliTRDcalibration()
125 // Default destructor
127 if(fNbTRDTrack) delete fNbTRDTrack;
128 if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;
129 if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;
130 if(fNbTRDTracklet) delete fNbTRDTracklet;
131 if(fNbTRDTrackletOffline) delete fNbTRDTrackletOffline;
132 if(fNbTRDTrackletStandalone) delete fNbTRDTrackletStandalone;
133 if(fNbTimeBin) delete fNbTimeBin;
134 if(fNbTimeBinOffline) delete fNbTimeBinOffline;
135 if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;
136 if(fNbClusters) delete fNbClusters;
137 if(fNbClustersOffline) delete fNbClustersOffline;
138 if(fNbClustersStandalone) delete fNbClustersStandalone;
139 if(fPHSM) delete fPHSM;
140 if(fCHSM) delete fCHSM;
141 if(fPHSum) delete fPHSum;
142 if(fCHSum) delete fCHSum;
143 if(fDetSum) delete fDetSum;
144 if(fDetSumVector) delete fDetSumVector;
145 if(fGraph){fGraph->Delete(); delete fGraph;}
146 if(fArrayCalib){fArrayCalib->Delete(); delete fArrayCalib;}
149 //________________________________________________________________________
150 void AliTRDcalibration::CreateOutputObjects()
152 // Create output objects
154 OpenFile(0, "RECREATE");
156 // Number of time bins
158 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
159 fNbTimeBins = cal->GetNumberOfTimeBinsDCS();
160 if(fNbTimeBins <= 0){
161 AliWarning(Form("No of TimeBins from DB [%d] use default [30]", fNbTimeBins));
166 // instance calibration: what to calibrate
167 fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
168 fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
169 fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
170 fTRDCalibraFillHisto->SetCH2dOn(); // choose to calibrate the gain
171 fTRDCalibraFillHisto->SetPH2dOn(); // choose to calibrate the drift velocity
172 fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
173 fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
174 fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
176 // segmentation (should be per default the max and add at the end)
177 for(Int_t k = 0; k < 3; k++){
178 if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
179 fTRDCalibraFillHisto->SetNz(k,fNz[k]); // Mode calibration
180 fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]); // Mode calibration
183 if((fNz[k] == 100) && (fNrphi[k] == 100)) {
184 if(fVector2d) AliInfo("The mode all together is not supported by the vector method");
185 fTRDCalibraFillHisto->SetAllTogether(k);
187 if((fNz[k] == 10) && (fNrphi[k] == 10)) {
188 if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");
189 fTRDCalibraFillHisto->SetPerSuperModule(k);
195 fTRDCalibraFillHisto->SetDebugLevel(DebugLevel()); //debug stuff
198 fTRDCalibraFillHisto->Init2Dhistos(fNbTimeBins); // initialise the histos
201 fTRDCalibraFillHisto->SetNumberClusters(flow); // At least flow clusters
202 fTRDCalibraFillHisto->SetNumberClustersf(fhigh); // The more fhigh clusters
203 fTRDCalibraFillHisto->SetFillWithZero(ffillZero); // Fill zeros
204 fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fnormalizeNbOfCluster); // For iterations
206 // Add them to the container
207 fContainer = new TObjArray();
209 fContainer->Add(fTRDCalibraFillHisto->GetCH2d()); //TH2I
210 fContainer->Add(fTRDCalibraFillHisto->GetPH2d()); //TProfile2D
211 fContainer->Add(fTRDCalibraFillHisto->GetPRF2d()); //TProfile2D
213 if(fVdriftLinear) fContainer->Add(fTRDCalibraFillHisto->GetVdriftLinearFit()); // Other drift velocity
214 if(fVector2d) fContainer->Add(fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
218 // Init the debug histos
219 fNbTRDTrack = new TH1F("TRDTrack","TRDTrack",500,0,500);
220 fNbTRDTrack->Sumw2();
221 fNbTRDTrackOffline = new TH1F("TRDTrackOffline","TRDTrackOffline",500,0,500);
222 fNbTRDTrackOffline->Sumw2();
223 fNbTRDTrackStandalone = new TH1F("TRDTrackStandalone","TRDTrackStandalone",500,0,500);
224 fNbTRDTrackStandalone->Sumw2();
226 fNbTRDTracklet = new TH1F("TRDTracklet","TRDTracklet",540,0.,540.);
227 fNbTRDTracklet->Sumw2();
228 fNbTRDTrackletOffline = new TH1F("TRDTrackletOffline","TRDTrackletOffline",540,0.,540.);
229 fNbTRDTrackletOffline->Sumw2();
230 fNbTRDTrackletStandalone = new TH1F("TRDTrackletStandalone","TRDTrackletStandalone",540,0.,540.);
231 fNbTRDTrackletStandalone->Sumw2();
233 fNbTimeBin = new TH1F("TimeBin","TimeBin",35,0,35);
235 fNbTimeBinOffline = new TH1F("TimeBinOffline","TimeBinOffline",35,0,35);
236 fNbTimeBinOffline->Sumw2();
237 fNbTimeBinStandalone = new TH1F("TimeBinStandalone","TimeBinStandalone",35,0,35);
238 fNbTimeBinStandalone->Sumw2();
240 fNbClusters = new TH1F("NbClusters","",35,0,35);
241 fNbClusters->Sumw2();
242 fNbClustersOffline = new TH1F("NbClustersOffline","",35,0,35);
243 fNbClustersOffline->Sumw2();
244 fNbClustersStandalone = new TH1F("NbClustersStandalone","",35,0,35);
245 fNbClustersStandalone->Sumw2();
247 fPHSM = new TProfile2D("PH2dSM","Nz10Nrphi10"
248 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
250 fPHSM->SetYTitle("Det/pad groups");
251 fPHSM->SetXTitle("time [#mus]");
252 fPHSM->SetZTitle("<PH> [a.u.]");
255 fCHSM = new TH2I("CH2dSM","Nz10Nrphi10",50,0,300,18,0,18);
256 fCHSM->SetYTitle("Det/pad groups");
257 fCHSM->SetXTitle("charge deposit [a.u]");
258 fCHSM->SetZTitle("counts");
262 fPHSum = new TProfile2D("PH2dSum","Nz100Nrphi100"
263 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
265 fPHSum->SetYTitle("Det/pad groups");
266 fPHSum->SetXTitle("time [#mus]");
267 fPHSum->SetZTitle("<PH> [a.u.]");
270 fCHSum = new TH2I("CH2dSum","Nz100Nrphi100",50,0,300,1,0,1);
271 fCHSum->SetYTitle("Det/pad groups");
272 fCHSum->SetXTitle("charge deposit [a.u]");
273 fCHSum->SetZTitle("counts");
278 fContainer->Add(fNbTRDTrack);
279 fContainer->Add(fNbTRDTrackOffline);
280 fContainer->Add(fNbTRDTrackStandalone);
281 fContainer->Add(fNbTRDTracklet);
282 fContainer->Add(fNbTRDTrackletOffline);
283 fContainer->Add(fNbTRDTrackletStandalone);
284 fContainer->Add(fNbTimeBin);
285 fContainer->Add(fNbTimeBinOffline);
286 fContainer->Add(fNbTimeBinStandalone);
287 fContainer->Add(fNbClusters);
288 fContainer->Add(fNbClustersOffline);
289 fContainer->Add(fNbClustersStandalone);
290 fContainer->Add(fPHSM);
291 fContainer->Add(fCHSM);
292 fContainer->Add(fPHSum);
293 fContainer->Add(fCHSum);
299 //________________________________________________________________________
300 void AliTRDcalibration::Exec(Option_t *)
303 // Execute function where the reference data are filled
309 Int_t nbTrdTracks = 0;
311 Int_t nbTrdTracksStandalone = 0;
313 Int_t nbTrdTracksOffline = 0;
317 // Loop on track in the event
319 //printf("Total of %d\n",fTracks->GetEntriesFast());
320 for(Int_t itrk=0; itrk < fTracks->GetEntriesFast(); itrk++){
322 //printf("itrk %d\n",itrk);
324 fTrackInfo = (AliTRDtrackInfo*)fTracks->UncheckedAt(itrk);
325 ftrdTrack = fTrackInfo->GetTrack();
326 if(!ftrdTrack) continue;
330 fTRDCalibraFillHisto->UpdateHistogramsV1(ftrdTrack);
334 Bool_t standalonetracklet = kFALSE;
335 const AliTRDseedV1 *tracklet = 0x0;
337 // Loop on tracklet in the event
339 for(Int_t itr = 0; itr < 6; itr++){
340 //printf("itr %d\n",itr);
341 if(!(tracklet = ftrdTrack->GetTracklet(itr))) continue;
342 if(!tracklet->IsOK()) continue;
344 if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
345 Int_t nbclusters = 0;
347 Double_t phtb[AliTRDseedV1::kNtb];
348 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
352 Float_t normalisation = 6.67;
356 // Check no shared clusters
357 for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
358 if((fcl = tracklet->GetClusters(icc))) crossrow = 1;
361 for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
362 if(!(fcl = tracklet->GetClusters(ic))) continue;
364 Int_t time = fcl->GetPadTime();
365 Float_t ch = tracklet->GetdQdl(ic);
366 Float_t qcl = TMath::Abs(fcl->GetQ());
367 detector = fcl->GetDetector();
368 if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
369 sum += ch/normalisation;
370 fNbTimeBin->Fill(time);
371 if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
372 else fNbTimeBinOffline->Fill(time);
374 sector = AliTRDgeometry::GetSector(detector);
376 fNbTRDTracklet->Fill(detector);
377 if(tracklet->IsStandAlone()) fNbTRDTrackletStandalone->Fill(detector);
378 else fNbTRDTrackletOffline->Fill(detector);
380 fNbClusters->Fill(nbclusters);
381 if(tracklet->IsStandAlone()) fNbClustersStandalone->Fill(nbclusters);
382 else fNbClustersOffline->Fill(nbclusters);
384 if((nbclusters > flow) && (nbclusters < fhigh)){
385 fCHSM->Fill(sum/20.0,sector+0.5);
386 fCHSum->Fill(sum/20.0,0.5);
387 for(int ic=0; ic<fNbTimeBins; ic++){
389 fPHSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
390 fPHSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
394 fPHSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
395 fPHSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
402 if(standalonetracklet) nbTrdTracksStandalone++;
403 else nbTrdTracksOffline++;
412 fNbTRDTrack->Fill(nbTrdTracks);
413 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
414 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
418 //printf("Nbof tracks %d\n",nbTrdTracks);
421 PostData(0, fContainer);
423 //printf("post container\n");
426 //________________________________________________________________________
427 void AliTRDcalibration::Terminate(Option_t *)
429 // Draw result to the screen
430 // Called once at the end of the query
432 //printf("terminate\n");
434 if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
437 //________________________________________________________
438 Bool_t AliTRDcalibration::GetRefFigure(Int_t ifig)
441 // Draw filled histos
444 gStyle->SetPalette(1);
445 gStyle->SetOptStat(1111);
446 gStyle->SetPadBorderMode(0);
447 gStyle->SetCanvasColor(10);
448 gStyle->SetPadLeftMargin(0.13);
449 gStyle->SetPadRightMargin(0.13);
451 if(!fContainer) return kFALSE;
455 TCanvas *c0 = new TCanvas("c0","c0",10,10,510,510);
456 TLegend *legNbTrack = new TLegend(.7, .7, .98, .98);
457 legNbTrack->SetBorderSize(1);
461 if(!(h = (TH1F *)fContainer->FindObject("TRDTrack"))) break;
462 if(!(ho = (TH1F *)fContainer->FindObject("TRDTrackOffline"))) break;
463 if(!(hs = (TH1F *)fContainer->FindObject("TRDTrackStandalone"))) break;
471 legNbTrack->AddEntry(h, "all", "p");
472 legNbTrack->AddEntry(ho, "offline", "p");
473 legNbTrack->AddEntry(hs, "standalone", "p");
474 legNbTrack->Draw("same");
478 TLegend *legNbTracklet = new TLegend(.7, .7, .98, .98);
479 legNbTracklet->SetBorderSize(1);
483 if(!(h = (TH1F *)fContainer->FindObject("TRDTracklet"))) break;
484 if(!(ho = (TH1F *)fContainer->FindObject("TRDTrackletOffline"))) break;
485 if(!(hs = (TH1F *)fContainer->FindObject("TRDTrackletStandalone"))) break;
489 legNbTracklet->AddEntry(h, "all", "p");
490 legNbTracklet->AddEntry(ho, "offline", "p");
491 legNbTracklet->AddEntry(hs, "standalone", "p");
492 legNbTracklet->Draw("same");
499 TLegend *legNbTimeBin = new TLegend(.7, .7, .98, .98);
500 legNbTimeBin->SetBorderSize(1);
504 if(!(h = (TH1F *)fContainer->FindObject("TimeBin"))) break;
505 if(!(ho = (TH1F *)fContainer->FindObject("TimeBinOffline"))) break;
506 if(!(hs = (TH1F *)fContainer->FindObject("TimeBinStandalone"))) break;
510 legNbTimeBin->AddEntry(h, "all", "p");
511 legNbTimeBin->AddEntry(ho, "offline", "p");
512 legNbTimeBin->AddEntry(hs, "standalone", "p");
513 legNbTimeBin->Draw("same");
520 TLegend *legNbClusters = new TLegend(.7, .7, .98, .98);
521 legNbClusters->SetBorderSize(1);
525 if(!(h = (TH1F *)fContainer->FindObject("NbClusters"))) break;
526 if(!(ho = (TH1F *)fContainer->FindObject("NbClustersOffline"))) break;
527 if(!(hs = (TH1F *)fContainer->FindObject("NbClustersStandalone"))) break;
531 legNbClusters->AddEntry(h, "all", "p");
532 legNbClusters->AddEntry(ho, "offline", "p");
533 legNbClusters->AddEntry(hs, "standalone", "p");
534 legNbClusters->Draw("same");
542 if(!(h = (TProfile2D *)fContainer->FindObject("PH2dSum"))) break;
543 TH1D *projh = h->ProjectionX("projh",1,1,"e");
551 if(!(h = (TH2I *)fContainer->FindObject("CH2dSum"))) break;
552 TH1D *projh = h->ProjectionX("projhh",1,1,"e");
560 AliInfo("Histo was not filled!");
564 if(!(h = (TProfile2D *)fContainer->FindObject("PH2d"))) break;
570 AliInfo("Histo was not filled!");
574 if(!(h = (TH2I *)fContainer->FindObject("CH2d"))) break;
580 AliInfo("Histo was not filled!");
584 if(!(h = (TProfile2D *)fContainer->FindObject("PRF2d"))) break;
590 AliInfo("vector was not filled!");
593 AliTRDCalibraVector *v = 0x0;
594 TGraphErrors *vdet = 0x0;
595 if(!(v = (AliTRDCalibraVector *)fContainer->FindObject("AliTRDCalibraVector"))) break;
596 Int_t detectormax = -1;
598 if(!v->FindTheMaxEntries(1,detectormax,groupmax)) break;
599 if(!(vdet = v->ConvertVectorPHTGraphErrors((Int_t)detectormax,groupmax,"plotPH2dVector"))) break;
600 Int_t nbeentries = 0;
601 TH1F *ko = v->CorrectTheError(vdet,nbeentries);
603 AliInfo(Form("There are %d entries in the detector %d and group %d",nbeentries,detectormax,groupmax));
608 AliInfo("vector was not filled!");
611 AliTRDCalibraVector *v = 0x0;
613 if(!(v = (AliTRDCalibraVector *)fContainer->FindObject("AliTRDCalibraVector"))) break;
614 Int_t detectormax = -1;
616 if(!v->FindTheMaxEntries(0,detectormax,groupmax)) break;
617 if(!(vdet = v->ConvertVectorCHHisto((Int_t)detectormax,groupmax,"plotCH2dVector"))) break;
619 AliInfo(Form("The detectormax and groupmax are %d and %d",detectormax,groupmax));
624 AliInfo("vector was not filled!");
627 AliTRDCalibraVector *v = 0x0;
628 TGraphErrors *vdet = 0x0;
629 if(!(v = (AliTRDCalibraVector *)fContainer->FindObject("AliTRDCalibraVector"))) break;
630 Int_t detectormax = -1;
632 Int_t nbeentries = 0;
633 if(!v->FindTheMaxEntries(2,detectormax,groupmax)) break;
634 if(!(vdet = v->ConvertVectorPRFTGraphErrors((Int_t)detectormax,groupmax,"plotPRF2dVector"))) break;
635 TH1F *ko = v->CorrectTheError(vdet,nbeentries);
637 AliInfo(Form("The detectormax and groupmax are %d and %d",detectormax,groupmax));
642 AliInfo("vdrift linear was not filled!");
645 AliTRDCalibraVdriftLinearFit *h = 0x0;
646 TH2S *hdetector = 0x0;
647 if(!(h = (AliTRDCalibraVdriftLinearFit *)fContainer->FindObject("AliTRDCalibraVdriftLinearFit"))) break;
648 Double_t entries[540];
649 for(Int_t k = 0; k < 540; k++){
652 if(!(hdetector = (TH2S *)h->GetLinearFitterHisto(k,kFALSE))) continue;
653 entries[k] = hdetector->GetEntries();
655 Double_t max = -10.0;
656 Double_t detectormax = -1;
657 for(Int_t k = 0; k < 540; k++){
658 if(entries[k] > max) {
664 if((TMath::Abs(max) <= 0.001) || (detectormax <0.0) || (detectormax >=540.0)) break;
665 if(!(hdetector = (TH2S *)h->GetLinearFitterHisto((Int_t)detectormax,kFALSE))) break;
666 AliInfo(Form("The detector with the maximum of entries is %d",detectormax));
672 if(!PostProcess()) break;
674 TGraph *fgain = (TGraph *) fGraph->At(0);
679 case kVdriftT0Factor:{
681 if(!PostProcess()) break;
683 TCanvas *c = new TCanvas("c","c",10,10,510,510);
685 TGraph *fvd = (TGraph *) fGraph->At(1);
690 TGraph *ft0 = (TGraph *) fGraph->At(2);
697 case kVdriftLorentzAngleFactor:{
699 AliInfo("vdrift linear was not filled!");
703 if(!PostProcess()) break;
705 TCanvas *c = new TCanvas("c","c",10,10,510,510);
707 TGraph *fvdl = (TGraph *) fGraph->At(3);
712 TGraph *flal = (TGraph *) fGraph->At(4);
721 if(!PostProcess()) break;
723 TGraph *fprf = (TGraph *) fGraph->At(5);
733 //________________________________________________________________________
734 Bool_t AliTRDcalibration::PostProcess()
737 // Fit the filled histos
738 // Put the calibration object in fArrayCalib
739 // 0 and 1 AliTRDCalDet and AliTRDCalPad gain
740 // 2 and 3 AliTRDCalDet and AliTRDCalPad vdrift PH
741 // 4 and 5 AliTRDCalDet and AliTRDCalPad timeoffset PH
742 // 6 AliTRDCalPad PRF
743 // 7 and 8 AliTRDCalDet and AliTRDCalPad vdrift second
744 // 9 and 10 AliTRDCalDet and AliTRDCalPad lorentz angle second
748 fArrayCalib = new TObjArray(11);
749 fArrayCalib->SetOwner();
757 fGraph = new TObjArray(6);
765 Bool_t storage[3] = {kFALSE,kFALSE,kFALSE};
767 // Objects for fitting
768 AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
769 calibra->SetDebugLevel(2); // 0 rien, 1 fitvoir, 2 debug files, 3 one detector
773 Printf("ERROR: list not available");
777 if(fHisto2d && fVector2d) AliInfo("We will only look at histos. Set fHisto2d off if you don't want");
778 AliTRDCalibraVector *calibraVector = 0x0;
779 if(fVector2d) calibraVector = (AliTRDCalibraVector *) fContainer->FindObject("CalibraVector");
783 Bool_t pass = kFALSE;
784 AliTRDCalibraVector *vvect = 0x0;
786 TH2I *histogain = (TH2I *) fContainer->FindObject("CH2d");
788 histogain->SetDirectory(0);
789 calibra->SetMinEntries(20);
790 if(fCompressPerDetector){
791 if(AddStatsPerDetector(histogain)) pass = calibra->AnalyseCH(fCHSum);
793 else pass = calibra->AnalyseCH(histogain);
797 if(fVector2d && calibraVector) {
798 calibra->SetMinEntries(20);
799 if(fCompressPerDetector){
800 if(!(vvect = calibraVector->AddStatsPerDetectorCH())) return kFALSE;
801 pass = calibra->AnalyseCH(vvect);
803 else pass = calibra->AnalyseCH(calibraVector);
809 Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(0))
810 + 6* 18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(0));
811 Int_t nbfit = calibra->GetNumberFit(); //Number of fits
812 Int_t nbE = calibra->GetNumberEnt(); //Number of detector mit entries
815 (nbfit >= 0.001*nbE))
817 // create the cal objects
818 calibra->PutMeanValueOtherVectorFit(1,kTRUE);
819 TObjArray object = calibra->GetVectorFit();
820 AliTRDCalDet *objgaindet = calibra->CreateDetObjectGain(&object);
821 TObject *objgainpad = calibra->CreatePadObjectGain();
823 fArrayCalib->AddAt(objgaindet,0);
824 fArrayCalib->AddAt(objgainpad,1);
828 if(FillGraphIndex(&object,graph)){
829 fGraph->AddAt(graph,0);
831 }//if(enough statistics?)
832 calibra->ResetVectorFit();
837 // VDRIFT average pulse height
841 TProfile2D *histodriftvelocity = (TProfile2D *) fContainer->FindObject("PH2d");
842 if(histodriftvelocity) {
843 histodriftvelocity->SetDirectory(0);
844 calibra->SetMinEntries(20*20);
845 if(fCompressPerDetector){
846 if(AddStatsPerDetector(histodriftvelocity,1)) {
847 pass = calibra->AnalysePH(fDetSumVector);
850 else pass = calibra->AnalysePH(histodriftvelocity);
854 if(fVector2d && calibraVector) {
855 calibra->SetMinEntries(20*20);
856 if(fCompressPerDetector){
857 if(!(vvect = calibraVector->AddStatsPerDetectorPH())) return kFALSE;
858 pass = calibra->AnalysePH(vvect);
860 else pass = calibra->AnalysePH(calibraVector);
865 Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(1))
866 + 6* 18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(1));
867 Int_t nbfit = calibra->GetNumberFit();
868 Int_t nbE = calibra->GetNumberEnt();
871 (nbfit >= 0.001*nbE))
873 // create the cal objects
874 calibra->PutMeanValueOtherVectorFit(1,kTRUE);
875 calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
876 TObjArray object = calibra->GetVectorFit();
877 AliTRDCalDet *objdriftvelocitydet = calibra->CreateDetObjectVdrift(&object,kTRUE);
878 TObject *objdriftvelocitypad = calibra->CreatePadObjectVdrift();
879 TObjArray objectt = calibra->GetVectorFit2();
880 AliTRDCalDet *objtime0det = calibra->CreateDetObjectT0(&object,kTRUE);
881 TObject *objtime0pad = calibra->CreatePadObjectT0();
883 fArrayCalib->AddAt(objdriftvelocitydet,2);
884 fArrayCalib->AddAt(objdriftvelocitypad,3);
886 fArrayCalib->AddAt(objtime0det,4);
887 fArrayCalib->AddAt(objtime0pad,5);
890 if(FillGraphIndex(&object,graph)){
891 fGraph->AddAt(graph,1);
893 TGraph *graphh = 0x0;
894 if(FillGraphIndex(&objectt,graphh)){
895 fGraph->AddAt(graphh,2);
897 }//if(enough statistics)
898 calibra->ResetVectorFit();
907 TProfile2D *histoprf = (TProfile2D *) fContainer->FindObject("PRF2d");
909 histoprf->SetDirectory(0);
910 calibra->SetMinEntries(600);
911 if(fCompressPerDetector){
912 if(AddStatsPerDetector(histoprf,2)) pass = calibra->AnalysePRFMarianFit(fDetSumVector);
914 else pass = calibra->AnalysePRFMarianFit(histoprf);
918 if(fVector2d && calibraVector) {
919 calibra->SetMinEntries(600);
920 if(fCompressPerDetector){
921 if(!(vvect =calibraVector->AddStatsPerDetectorPRF())) return kFALSE;
922 pass = calibra->AnalysePRFMarianFit(vvect);
924 else pass = calibra->AnalysePRFMarianFit(calibraVector);
929 Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
930 + 6* 18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
931 Int_t nbfit = calibra->GetNumberFit();
932 Int_t nbE = calibra->GetNumberEnt();
935 (nbfit >= 0.001*nbE)) {
936 TObjArray object = calibra->GetVectorFit();
937 TObject *objPRFpad = calibra->CreatePadObjectPRF(&object);
939 fArrayCalib->AddAt(objPRFpad,6);
942 if(FillGraphIndex(&object,graph)){
943 fGraph->AddAt(graph,5);
946 calibra->ResetVectorFit();
953 AliTRDCalibraVdriftLinearFit *vlinearfit = (AliTRDCalibraVdriftLinearFit *) fContainer->FindObject("LinearVdriftFit");
955 calibra->SetMinEntries(20*20);
956 if(calibra->AnalyseLinearFitters(vlinearfit)) {
958 Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
959 + 6* 18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
960 Int_t nbfit = calibra->GetNumberFit();
961 Int_t nbE = calibra->GetNumberEnt();
964 (nbfit >= 0.001*nbE))
966 // create the cal objects
967 calibra->PutMeanValueOtherVectorFit(1,kTRUE);
968 calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
969 TObjArray object = calibra->GetVectorFit();
970 AliTRDCalDet *objdriftvelocitydet = calibra->CreateDetObjectVdrift(&object,kTRUE);
971 TObject *objdriftvelocitypad = calibra->CreatePadObjectVdrift();
972 TObjArray objectt = calibra->GetVectorFit2();
973 AliTRDCalDet *objtime0det = calibra->CreateDetObjectT0(&object,kTRUE);
974 TObject *objtime0pad = calibra->CreatePadObjectT0();
976 fArrayCalib->AddAt(objdriftvelocitydet,7);
977 fArrayCalib->AddAt(objdriftvelocitypad,8);
979 fArrayCalib->AddAt(objtime0det,9);
980 fArrayCalib->AddAt(objtime0pad,10);
983 if(FillGraphIndex(&object,graph)){
984 fGraph->AddAt(graph,3);
986 TGraph *graphh = 0x0;
987 if(FillGraphIndex(&objectt,graphh)){
988 fGraph->AddAt(graphh,4);
990 }//if(enough statistics)
992 calibra->ResetVectorFit();
996 fPostProcess = kTRUE;
1002 //________________________________________________________________________
1003 Bool_t AliTRDcalibration::FillGraphIndex(const TObjArray *vectora,TGraph *graph) const
1006 // Fill one value (the first one) per detector
1009 Int_t loop = (Int_t) vectora->GetEntriesFast();
1011 AliInfo("The Vector Fit is not complete!");
1017 for (Int_t k = 0; k < loop; k++) {
1018 if(!vectora->At(k)){
1023 x[k] = ((AliTRDCalibraFit::AliTRDFitInfo *) vectora->At(k))->GetDetector();
1024 y[k] = ((Float_t *)((AliTRDCalibraFit::AliTRDFitInfo *) vectora->At(k))->GetCoef())[0];
1027 if(!graph) graph = new TGraph(540,&x[0],&y[0]);
1030 new(graph) TGraph(540,&x[0],&y[0]);
1036 //________________________________________________________________________
1037 Bool_t AliTRDcalibration::AddStatsPerDetector(const TH2I *ch)
1040 // Add statistic per detector
1043 AliTRDCalibraMode calibMode = AliTRDCalibraMode();
1044 const char *name = ch->GetTitle();
1045 if(!SetNzFromTObject(name,0,&calibMode)) return 0x0;
1046 if(!SetNrphiFromTObject(name,0,&calibMode)) return 0x0;
1047 if(((calibMode.GetNz(0) == 100) && (calibMode.GetNrphi(0) == 100)) || ((calibMode.GetNz(0) == 10) && (calibMode.GetNrphi(0) == 10))) return kFALSE;
1049 Int_t nybins = ch->GetNbinsY();// groups number
1050 Int_t nxbins = ch->GetNbinsX();// number of bins X
1051 TAxis *xaxis = ch->GetXaxis();
1052 Double_t lowedge = xaxis->GetBinLowEdge(1);
1053 Double_t upedge = xaxis->GetBinUpEdge(nxbins);
1055 // number per chamber 2
1056 calibMode.ModePadCalibration(2,0);
1057 calibMode.ModePadFragmentation(0,2,0,0);
1058 calibMode.SetDetChamb2(0);
1059 Int_t perChamber2 = (Int_t) calibMode.GetDetChamb2(0);
1061 // number per other chamber
1062 calibMode.ModePadCalibration(0,0);
1063 calibMode.ModePadFragmentation(0,0,0,0);
1064 calibMode.SetDetChamb0(0);
1065 Int_t perChamber0 = (Int_t) calibMode.GetDetChamb0(0);
1067 if(nybins != (6*18*perChamber2+6*4*18*perChamber0)) return kFALSE;
1070 TString nname((const char *)ch->GetName());
1071 nname += "PerDetector";
1072 TString title("Nz");
1076 if(!fCHSum) fCHSum = new TH2I((const char *)nname,(const char *)title
1077 ,nxbins,lowedge,upedge,540,0,540);
1080 new(fCHSum) TH2I((const Char_t *) nname,(const char *)title
1081 ,nxbins,lowedge,upedge,540,0,540);
1083 fCHSum->SetYTitle("Detector number");
1084 fCHSum->SetXTitle("charge deposit [a.u]");
1085 fCHSum->SetZTitle("counts");
1086 fCHSum->SetStats(0);
1091 for(Int_t det = 0; det < 540; det++){
1093 Int_t numberofgroup = 0;
1094 if(AliTRDgeometry::GetStack(det) == 2) numberofgroup = perChamber2;
1095 else numberofgroup = perChamber0;
1096 TH1I *projch = (TH1I *) ch->ProjectionX("projch",counter+1,counter+numberofgroup,(Option_t *)"e");
1097 projch->SetDirectory(0);
1099 for(Int_t nx = 0; nx <= nxbins; nx++) {
1100 fCHSum->SetBinContent(nx,det+1,projch->GetBinContent(nx));
1101 fCHSum->SetBinError(nx,det+1,projch->GetBinError(nx));
1104 counter += numberofgroup;
1113 //_____________________________________________________________________________________________________________________
1114 Bool_t AliTRDcalibration::AddStatsPerDetector(const TProfile2D *ph,Int_t i)
1117 // Add statistic per detector
1120 AliTRDCalibraMode calibMode = AliTRDCalibraMode();
1121 const char *name = ph->GetTitle();
1122 //printf("name %s\n",name);
1123 if(!SetNzFromTObject(name,0,&calibMode)) return kFALSE;
1124 if(!SetNrphiFromTObject(name,0,&calibMode)) return kFALSE;
1125 if(((calibMode.GetNz(0) == 100) && (calibMode.GetNrphi(0) == 100)) || ((calibMode.GetNz(0) == 10) && (calibMode.GetNrphi(0) == 10))) return kFALSE;
1126 //printf("Found mode Mz %d, Nrphi %d\n",calibMode.GetNz(0),calibMode.GetNrphi(0));
1129 Int_t nybins = ph->GetNbinsY();// groups number
1130 Int_t nxbins = ph->GetNbinsX();// number of bins X
1131 TAxis *xaxis = ph->GetXaxis();
1132 Double_t lowedge = xaxis->GetBinLowEdge(1);
1133 Double_t upedge = xaxis->GetBinUpEdge(nxbins);
1135 // number per chamber 2
1136 calibMode.ModePadCalibration(2,0);
1137 calibMode.ModePadFragmentation(0,2,0,0);
1138 calibMode.SetDetChamb2(0);
1139 Int_t perChamber2 = (Int_t) calibMode.GetDetChamb2(0);
1141 // number per other chamber
1142 calibMode.ModePadCalibration(0,0);
1143 calibMode.ModePadFragmentation(0,0,0,0);
1144 calibMode.SetDetChamb0(0);
1145 Int_t perChamber0 = (Int_t) calibMode.GetDetChamb0(0);
1147 if(nybins != (6*18*perChamber2+6*4*18*perChamber0)) return kFALSE;
1150 TString nbname((const char *)ph->GetName());
1151 nbname += "PerDetectorVector";
1152 if(!fDetSumVector) fDetSumVector = new AliTRDCalibraVector();
1154 fDetSumVector->~AliTRDCalibraVector();
1155 new(fDetSumVector) AliTRDCalibraVector();
1158 fDetSumVector->SetTimeMax(nxbins);
1161 fDetSumVector->SetNumberBinPRF(nxbins);
1162 fDetSumVector->SetPRFRange(TMath::Abs(lowedge));
1164 fDetSumVector->SetDetCha0(i,1);
1165 fDetSumVector->SetDetCha2(i,1);
1166 fDetSumVector->SetNzNrphi(i,0,0);
1168 Int_t nbg = GetNumberOfGroupsPRF((const char *)name);
1169 fDetSumVector->SetNbGroupPRF(nbg);
1173 TString nname((const char *)ph->GetName());
1174 nname += "PerDetector";
1175 TString title("Nz");
1179 if(!fDetSum) fDetSum = new TH2D((const char *)nname,(const Char_t *) title
1180 ,nxbins,lowedge,upedge,540,0,540);
1183 new(fDetSum) TH2D((const Char_t *) nname,(const Char_t *) title
1184 ,nxbins,lowedge,upedge,540,0,540);
1186 fDetSum->SetYTitle("Detector number");
1187 fDetSum->SetXTitle(xaxis->GetTitle());
1188 fDetSum->SetStats(0);
1192 for(Int_t det = 0; det < 540; det++){
1194 Int_t numberofgroup = 0;
1195 if(AliTRDgeometry::GetStack(det) == 2) numberofgroup = perChamber2;
1196 else numberofgroup = perChamber0;
1198 for(Int_t nx = 1; nx <= nxbins; nx++) {
1200 Double_t entries = 0.0;
1201 Double_t sumw2 = 0.0;
1202 Double_t sumw = 0.0;
1204 for(Int_t k = counter+1; k <= (counter+numberofgroup); k++){
1205 Int_t binnumber = ph->GetBin(nx,k);
1206 entries += ph->GetBinEntries(binnumber);
1207 sumw2 += (ph->GetBinError(binnumber)*ph->GetBinError(binnumber)+ph->GetBinContent(binnumber)*ph->GetBinContent(binnumber))*ph->GetBinEntries(binnumber);
1208 sumw += ph->GetBinContent(binnumber)*ph->GetBinEntries(binnumber);
1211 Double_t mean = 0.0;
1212 if(entries > 0.0) mean = sumw/entries;
1213 Double_t squaremean = 0.0;
1214 if(entries > 0.0) squaremean = sumw2/entries;
1215 Double_t errorf = squaremean - mean*mean;
1216 Double_t error = 0.0;
1217 if(entries > 0.0) error = TMath::Sqrt(TMath::Abs(errorf)/entries);
1219 fDetSum->SetBinContent(nx,det+1,mean);
1220 fDetSum->SetBinError(nx,det+1,error);
1222 if(i==1) fDetSumVector->FillVectorPH(det,0,nx-1,(Int_t)entries,(Float_t)mean,(Float_t)squaremean);
1223 if(i==2) fDetSumVector->FillVectorPRF(det,0,nx-1,(Int_t)entries,(Float_t)mean,(Float_t)squaremean);
1227 counter += numberofgroup;
1235 //_____________________________________________________________________________
1236 Bool_t AliTRDcalibration::SetNrphiFromTObject(const char *name, Int_t i, AliTRDCalibraMode *calibMode) const
1239 // Set the granularity from object
1242 const Char_t *patternrphi0 = "Nrphi0";
1243 const Char_t *patternrphi1 = "Nrphi1";
1244 const Char_t *patternrphi2 = "Nrphi2";
1245 const Char_t *patternrphi3 = "Nrphi3";
1246 const Char_t *patternrphi4 = "Nrphi4";
1247 const Char_t *patternrphi5 = "Nrphi5";
1248 const Char_t *patternrphi6 = "Nrphi6";
1251 const Char_t *patternrphi10 = "Nrphi10";
1252 const Char_t *patternrphi100 = "Nrphi100";
1253 const Char_t *patternz10 = "Nz10";
1254 const Char_t *patternz100 = "Nz100";
1257 if ((strstr(name,patternrphi100)) && (strstr(name,patternz100))) {
1258 calibMode->SetAllTogether(i);
1261 if ((strstr(name,patternrphi10)) && (strstr(name,patternz10))) {
1262 calibMode->SetPerSuperModule(i);
1266 if (strstr(name,patternrphi0)) {
1267 calibMode->SetNrphi(i ,0);
1270 if (strstr(name,patternrphi1)) {
1271 calibMode->SetNrphi(i, 1);
1274 if (strstr(name,patternrphi2)) {
1275 calibMode->SetNrphi(i, 2);
1278 if (strstr(name,patternrphi3)) {
1279 calibMode->SetNrphi(i, 3);
1282 if (strstr(name,patternrphi4)) {
1283 calibMode->SetNrphi(i, 4);
1286 if (strstr(name,patternrphi5)) {
1287 calibMode->SetNrphi(i, 5);
1290 if (strstr(name,patternrphi6)) {
1291 calibMode->SetNrphi(i, 6);
1295 calibMode->SetNrphi(i ,0);
1299 //_____________________________________________________________________________
1300 Bool_t AliTRDcalibration::SetNzFromTObject(const char *name, Int_t i, AliTRDCalibraMode *calibMode) const
1303 // Set fNz[i] of the AliTRDCalibraFit::Instance()
1304 // corresponding to the given TObject
1308 const Char_t *patternz0 = "Nz0";
1309 const Char_t *patternz1 = "Nz1";
1310 const Char_t *patternz2 = "Nz2";
1311 const Char_t *patternz3 = "Nz3";
1312 const Char_t *patternz4 = "Nz4";
1314 const Char_t *patternrphi10 = "Nrphi10";
1315 const Char_t *patternrphi100 = "Nrphi100";
1316 const Char_t *patternz10 = "Nz10";
1317 const Char_t *patternz100 = "Nz100";
1319 if ((strstr(name,patternrphi100)) && (strstr(name,patternz100))) {
1320 calibMode->SetAllTogether(i);
1323 if ((strstr(name,patternrphi10)) && (strstr(name,patternz10))) {
1324 calibMode->SetPerSuperModule(i);
1327 if (strstr(name,patternz0)) {
1328 calibMode->SetNz(i, 0);
1331 if (strstr(name,patternz1)) {
1332 calibMode->SetNz(i ,1);
1335 if (strstr(name,patternz2)) {
1336 calibMode->SetNz(i ,2);
1339 if (strstr(name,patternz3)) {
1340 calibMode->SetNz(i ,3);
1343 if (strstr(name,patternz4)) {
1344 calibMode->SetNz(i ,4);
1348 calibMode->SetNz(i ,0);
1351 //____________________________________________________________________________________________________
1352 Int_t AliTRDcalibration::GetNumberOfGroupsPRF(const char* nametitle) const
1355 // Get numberofgroupsprf
1359 const Char_t *pattern0 = "Ngp0";
1360 const Char_t *pattern1 = "Ngp1";
1361 const Char_t *pattern2 = "Ngp2";
1362 const Char_t *pattern3 = "Ngp3";
1363 const Char_t *pattern4 = "Ngp4";
1364 const Char_t *pattern5 = "Ngp5";
1365 const Char_t *pattern6 = "Ngp6";
1368 if (strstr(nametitle,pattern0)) {
1371 if (strstr(nametitle,pattern1)) {
1374 if (strstr(nametitle,pattern2)) {
1377 if (strstr(nametitle,pattern3)) {
1380 if (strstr(nametitle,pattern4)) {
1383 if (strstr(nametitle,pattern5)) {
1386 if (strstr(nametitle,pattern6)){