]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - TRD/AliTRDCalibTask.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
xxx
[u/mrichter/AliRoot.git] / TRD / AliTRDCalibTask.cxx
... / ...
CommitLineData
1
2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
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 **************************************************************************/
16
17/////////////////////////////////////////////////////////////////////////////////
18//
19// AliTRDCalibTask
20//
21// Offline TRD calibration task
22//
23// Author:
24// R. Bailhache (R.Bailhache@gsi.de)
25//
26//////////////////////////////////////////////////////////////////////////////////////
27
28#include <iostream>
29using namespace std;
30#include "Riostream.h"
31#include "TChain.h"
32#include "TTree.h"
33#include "TFile.h"
34#include "TProfile2D.h"
35#include "TH2I.h"
36#include "TH1F.h"
37#include "TH2F.h"
38#include "TH2S.h"
39#include "TList.h"
40#include "TMath.h"
41#include "TObject.h"
42#include "TObjArray.h"
43#include "TString.h"
44#include "TCanvas.h"
45#include "TLegend.h"
46#include "TStyle.h"
47#include "TLine.h"
48#include "TIterator.h"
49#include "TLinearFitter.h"
50#include "TVectorD.h"
51
52#include "AliAnalysisTask.h"
53#include "AliAnalysisManager.h"
54
55#include "AliExternalTrackParam.h"
56#include "AliESDVertex.h"
57#include "AliESDEvent.h"
58#include "AliESDfriend.h"
59#include "AliCentrality.h"
60#include "AliESDInputHandler.h"
61#include "AliESDtrack.h"
62#include "AliESDfriendTrack.h"
63#include "AliTRDtrackV1.h"
64#include "AliTRDseedV1.h"
65#include "AliTRDcluster.h"
66#include "AliTRDgeometry.h"
67#include "AliESDtrackCuts.h"
68#include "AliESDVertex.h"
69#include "AliTRDCalDet.h"
70
71#include "AliTRDCalibraVector.h"
72#include "AliTRDCalibraFillHisto.h"
73#include "AliTRDCalibraVdriftLinearFit.h"
74#include "AliTRDCalibraExbAltFit.h"
75
76#include "AliTRDcalibDB.h"
77#include "AliCDBId.h"
78#include "AliLog.h"
79
80#include "AliTRDCalibChamberStatus.h"
81#include "AliTRDdEdxBaseUtils.h"
82#include "AliTRDdEdxCalibHistArray.h"
83#include "AliTRDdEdxCalibUtils.h"
84
85#include "AliTRDCalibTask.h"
86
87
88ClassImp(AliTRDCalibTask)
89
90//________________________________________________________________________
91AliTRDCalibTask::AliTRDCalibTask(const char *name)
92: AliAnalysisTaskSE(name), fESD(0),
93 fkEsdTrack(0),
94 fFriendTrack(0),
95 fCalibObject(0),
96 fTrdTrack(0),
97 fCl(0),
98 fListHist(0),
99 fTRDCalibraFillHisto(0),
100 fTRDChamberStatus(0),
101 fNEvents(0),
102 fNEventsInput(0),
103 fNbTRDTrack(0),
104 fNbTRDTrackOffline(0),
105 fNbTRDTrackStandalone(0),
106 fNbTPCTRDtrack(0),
107 fNbGoodTracks(0),
108 fNbTimeBin(0),
109 fNbTimeBinOffline(0),
110 fNbTimeBinStandalone(0),
111 fNbClusters(0),
112 fNbClustersOffline(0),
113 fNbClustersStandalone(0),
114 fNbTracklets(0),
115 fNbTrackletsOffline(0),
116 fNbTrackletsStandalone(0),
117 fAbsoluteGain(0),
118 fCH2dSum(0),
119 fPH2dSum(0),
120 fCH2dSM(0),
121 fPH2dSM(0),
122 fCH2dTest(0),
123 fPH2dTest(0),
124 fLinearVdriftTest(0),
125 fOnInstance(kTRUE),
126 fHisto2d(kTRUE),
127 fVector2d(kFALSE),
128 fVdriftLinear(kTRUE),
129 fExbAlt(kFALSE),
130 fDebugLevelTRDCalibraFillHisto(0),
131 fNbTimeBins(0),
132 fNumberBinCharge(100),
133 fRangeCharge(150),
134 fVdBindx(32),
135 fVdBindy(70),
136 fVdRangex(0.8),
137 fVdRangey(1.4),
138 fSelectTrigger(kTRUE),
139 fSelectedTrigger(new TObjArray()),
140 fRejected(kTRUE),
141 fEsdTrackCuts(0),
142 fRequirePrimaryVertex(kFALSE),
143 fVtxTPC(kFALSE),
144 fVtxSPD(kFALSE),
145 fMinNbContributors(0),
146 fRangePrimaryVertexZ(9999999.0),
147 fRejectPileUpWithSPD(kFALSE),
148 fMinNbTracks(9),
149 fMaxNbTracks(999999999),
150 fCutWithVdriftCalib(kFALSE),
151 fMinNbTRDtracklets(0),
152 fMinTRDMomentum(0),
153 fScaleGainWithTPCSignal(kFALSE),
154 fLow(0),
155 fHigh(30),
156 fFillZero(kFALSE),
157 fNormalizeNbOfCluster(kFALSE),
158 fRelativeScale(0.0),
159 fMaxCluster(100.0),
160 fNbMaxCluster(2),
161 fOfflineTracks(kFALSE),
162 fStandaloneTracks(kFALSE),
163 fFirstRunGain(-1),
164 fVersionGainUsed(-1),
165 fSubVersionGainUsed(-1),
166 fFirstRunGainLocal(-1),
167 fVersionGainLocalUsed(-1),
168 fSubVersionGainLocalUsed(-1),
169 fFirstRunVdrift(-1),
170 fVersionVdriftUsed(-1),
171 fSubVersionVdriftUsed(-1),
172 fFirstRunExB(-1),
173 fVersionExBUsed(-1),
174 fSubVersionExBUsed(-1),
175 fCalDetGain(0x0),
176 fMaxEvent(0),
177 fCounter(0),
178 fPHQon(kTRUE)
179{
180 //
181 // Default constructor
182 //
183
184 fNz[0] = 0;
185 fNz[1] = 0;
186 fNz[2] = 0;
187
188 fNrphi[0] = 0;
189 fNrphi[1] = 0;
190 fNrphi[2] = 0;
191
192
193 // Define input and output slots here
194 // Input slot #0 works with a TChain
195 DefineInput(0, TChain::Class());
196
197 // Output slot #0 writes into a TList container
198 DefineOutput(1, TList::Class());
199
200
201}
202//____________________________________________________________________________________
203AliTRDCalibTask::~AliTRDCalibTask()
204{
205 //
206 // AliTRDCalibTask destructor
207 //
208
209 // Pointeur
210 if(fNEvents) delete fNEvents;
211 if(fNEventsInput) delete fNEventsInput;
212 if(fNbTRDTrack) delete fNbTRDTrack;
213 if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;
214 if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;
215 if(fNbTPCTRDtrack) delete fNbTPCTRDtrack;
216 if(fNbGoodTracks) delete fNbGoodTracks;
217 if(fNbTimeBin) delete fNbTimeBin;
218 if(fNbTimeBinOffline) delete fNbTimeBinOffline;
219 if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;
220 if(fNbClusters) delete fNbClusters;
221 if(fNbClustersOffline) delete fNbClustersOffline;
222 if(fNbClustersStandalone) delete fNbClustersStandalone;
223 if(fNbTracklets) delete fNbTracklets;
224 if(fNbTrackletsOffline) delete fNbTrackletsOffline;
225 if(fNbTrackletsStandalone) delete fNbTrackletsStandalone;
226 if(fAbsoluteGain) delete fAbsoluteGain;
227 if(fCH2dSum) delete fCH2dSum;
228 if(fPH2dSum) delete fPH2dSum;
229 if(fCH2dSM) delete fCH2dSM;
230 if(fPH2dSM) delete fPH2dSM;
231 if(fCH2dTest) delete fCH2dTest;
232 if(fPH2dTest) delete fPH2dTest;
233 if(fLinearVdriftTest) delete fLinearVdriftTest;
234 if(IsPHQon()){
235 AliTRDdEdxCalibUtils::DeleteHistArray();
236 }
237
238 if(fCalDetGain) delete fCalDetGain;
239
240 if(fSelectedTrigger) {
241 fSelectedTrigger->Delete();
242 delete fSelectedTrigger;
243 }
244 if(fEsdTrackCuts) {
245 delete fEsdTrackCuts;
246 }
247
248 if(fTRDChamberStatus) delete fTRDChamberStatus;
249
250}
251
252//________________________________________________________________________
253void AliTRDCalibTask::UserCreateOutputObjects()
254{
255 //
256 // Create the histos
257 //
258 //cout << "AliTRDCalibTask::CreateOutputObjects() IN" << endl;
259
260 // Number of time bins
261 if(fNbTimeBins==0) {
262 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
263 fNbTimeBins = cal->GetNumberOfTimeBinsDCS();
264 if(fNbTimeBins <= 0){
265 AliWarning(Form("No of TimeBins from DB [%d] use default [30]", fNbTimeBins));
266 fNbTimeBins = 30;
267 }
268 }
269
270 // output list
271 fListHist = new TList();
272 fListHist->SetOwner();
273
274 // init chamber status
275 fTRDChamberStatus = new AliTRDCalibChamberStatus();
276 fTRDChamberStatus->Init();
277 fListHist->Add(fTRDChamberStatus->GetSparseI());
278
279 // instance calibration
280 fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
281 if(fOnInstance) {
282 fTRDCalibraFillHisto->SetNumberBinCharge(fNumberBinCharge); // set number of bin of the charge distribution
283 fTRDCalibraFillHisto->SetCutWithVdriftCalib(fCutWithVdriftCalib); // cut vdrift
284 fTRDCalibraFillHisto->SetMinNbTRDtracklets(fMinNbTRDtracklets); // min number of TRD tracklets
285 fTRDCalibraFillHisto->SetMinTRDMomentum(fMinTRDMomentum); // min TRD momentum
286 fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
287 fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
288 fTRDCalibraFillHisto->SetCH2dOn(); // choose to calibrate the gain
289 fTRDCalibraFillHisto->SetPH2dOn(); // choose to calibrate the drift velocity
290 fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
291 fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
292 fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
293 fTRDCalibraFillHisto->SetExbAltFitOn(fExbAlt); // Alternative method for exb
294 fTRDCalibraFillHisto->SetScaleWithTPCSignal(fScaleGainWithTPCSignal); // Scale Gain with TPC signal
295 for(Int_t k = 0; k < 3; k++){
296 if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
297 fTRDCalibraFillHisto->SetNz(k,fNz[k]); // Mode calibration
298 fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]); // Mode calibration
299 }
300 else {
301 if((fNz[k] == 100) && (fNrphi[k] == 100)) {
302 if(fVector2d) AliInfo("The mode all together is not supported by the vector method");
303 fTRDCalibraFillHisto->SetAllTogether(k);
304 }
305 if((fNz[k] == 10) && (fNrphi[k] == 10)) {
306 if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");
307 fTRDCalibraFillHisto->SetPerSuperModule(k);
308 }
309 }
310 }
311 // Variables for how to fill
312 fTRDCalibraFillHisto->SetFillWithZero(fFillZero);
313 fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fNormalizeNbOfCluster);
314 fTRDCalibraFillHisto->SetMaxCluster(fMaxCluster);
315 fTRDCalibraFillHisto->SetNbMaxCluster(fNbMaxCluster);
316
317 // Init with 30 timebins
318 fTRDCalibraFillHisto->Init2Dhistos(fNbTimeBins); // initialise the histos
319 fTRDCalibraFillHisto->SetNumberClusters(fLow); // At least 11 clusters
320 fTRDCalibraFillHisto->SetNumberClustersf(fHigh); // At least 11 clusters
321
322 // For testing only
323 fTRDCalibraFillHisto->SetDebugLevel(fDebugLevelTRDCalibraFillHisto); //debug stuff
324
325 if(fHisto2d) {
326 fListHist->Add(fTRDCalibraFillHisto->GetCH2d());
327 fListHist->Add(fTRDCalibraFillHisto->GetPH2d());
328 fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());
329 }
330 if(fVdriftLinear) {
331 AliTRDCalibraVdriftLinearFit *fvdl = fTRDCalibraFillHisto->GetVdriftLinearFit();
332 if(fvdl){
333 fvdl->SetNbBindx(fVdBindx);
334 fvdl->SetNbBindy(fVdBindy);
335 fvdl->SetRangedx(fVdRangex);
336 fvdl->SetRangedy(fVdRangey);
337 }
338 fListHist->Add((TObject *)fTRDCalibraFillHisto->GetVdriftLinearFit());
339 }
340 if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
341 if(fExbAlt) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetExbAltFit());
342 }
343 fRelativeScale = fTRDCalibraFillHisto->GetRelativeScale(); // Get the relative scale for the gain
344
345 fNEvents = new TH1I(Form("NEvents_%s",(const char*)fName),"NEvents", 2, 0, 2);
346 fListHist->Add(fNEvents);
347 fNEventsInput = new TH1I(Form("NEventsInput_%s",(const char*)fName),"NEventsInput", 2, 0, 2);
348 fListHist->Add(fNEventsInput);
349
350 // absolute gain calibration even without AliESDfriend
351 Int_t nBinsPt = 25;
352 Double_t minPt = 0.001;
353 Double_t maxPt = 10.0;
354
355 Double_t *binLimLogPt = new Double_t[nBinsPt+1];
356 Double_t *binLimPt = new Double_t[nBinsPt+1];
357 for(Int_t i=0; i<=nBinsPt; i++) binLimLogPt[i]=(Double_t)TMath::Log10(minPt) + (TMath::Log10(maxPt)-TMath::Log10(minPt))/nBinsPt*(Double_t)i ;
358 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
359
360 fAbsoluteGain = new TH2F(Form("AbsoluteGain_%s",(const char*)fName),"AbsoluteGain", 200, 0.0, 700.0, nBinsPt, binLimPt);
361 fAbsoluteGain->SetYTitle("Momentum at TRD");
362 fAbsoluteGain->SetXTitle("charge deposit [a.u]");
363 fAbsoluteGain->SetZTitle("counts");
364 fAbsoluteGain->SetStats(0);
365 fAbsoluteGain->Sumw2();
366 fListHist->Add(fAbsoluteGain);
367
368 if(IsPHQon()){
369 printf("\n AliTRDCalibTask PHQ is on!! \n\n");
370 AliTRDdEdxBaseUtils::PrintControl();
371 AliTRDdEdxCalibUtils::IniHistArray(fListHist, kTRUE);
372 }
373 else{
374 printf("\n AliTRDCalibTask PHQ is off!! \n\n");
375 }
376
377 /////////////////////////////////////////
378 // First debug level
379 ///////////////////////////////////////
380 if(fDebug > 0) {
381
382 fLinearVdriftTest = new TH2S(Form("LFDV0testversion_%s",(const char*)fName),"LFDV0testversion",36,-0.9,0.9,48,-1.2,1.2);
383 fLinearVdriftTest->SetXTitle("tan(phi_{track})");
384 fLinearVdriftTest->SetYTitle("dy/dt");
385 fLinearVdriftTest->SetZTitle("Number of tracklets");
386 fLinearVdriftTest->SetStats(0);
387 fLinearVdriftTest->SetDirectory(0);
388
389 // Standart with AliESDfriend
390 fPH2dTest = new TProfile2D(Form("PH2dTest_%s",(const char*)fName),"Nz0Nrphi0"
391 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
392 ,540,0,540);
393 fPH2dTest->SetYTitle("Det/pad groups");
394 fPH2dTest->SetXTitle("time [#mus]");
395 fPH2dTest->SetZTitle("<PH> [a.u.]");
396 fPH2dTest->SetStats(0);
397 //
398 fCH2dTest = new TH2I(Form("CH2dTest_%s",(const char*)fName),"Nz0Nrphi0",50,0,300,540,0,540);
399 fCH2dTest->SetYTitle("Det/pad groups");
400 fCH2dTest->SetXTitle("charge deposit [a.u]");
401 fCH2dTest->SetZTitle("counts");
402 fCH2dTest->SetStats(0);
403 fCH2dTest->Sumw2();
404
405 //
406 fPH2dSM = new TProfile2D(Form("PH2dSM_%s",(const char*)fName),"Nz10Nrphi10"
407 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
408 ,18,0,18);
409 fPH2dSM->SetYTitle("Det/pad groups");
410 fPH2dSM->SetXTitle("time [#mus]");
411 fPH2dSM->SetZTitle("<PH> [a.u.]");
412 fPH2dSM->SetStats(0);
413 //
414 fCH2dSM = new TH2I(Form("CH2dSM_%s",(const char*)fName),"Nz10Nrphi10",50,0,300,18,0,18);
415 fCH2dSM->SetYTitle("Det/pad groups");
416 fCH2dSM->SetXTitle("charge deposit [a.u]");
417 fCH2dSM->SetZTitle("counts");
418 fCH2dSM->SetStats(0);
419 fCH2dSM->Sumw2();
420 //
421 fPH2dSum = new TProfile2D(Form("PH2dSum_%s",(const char*)fName),"Nz100Nrphi100"
422 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
423 ,1,0,1);
424 fPH2dSum->SetYTitle("Det/pad groups");
425 fPH2dSum->SetXTitle("time [#mus]");
426 fPH2dSum->SetZTitle("<PH> [a.u.]");
427 fPH2dSum->SetStats(0);
428 //
429 fCH2dSum = new TH2I(Form("CH2dSum_%s",(const char*)fName),"Nz100Nrphi100",50,0,300,1,0,1);
430 fCH2dSum->SetYTitle("Det/pad groups");
431 fCH2dSum->SetXTitle("charge deposit [a.u]");
432 fCH2dSum->SetZTitle("counts");
433 fCH2dSum->SetStats(0);
434 fCH2dSum->Sumw2();
435
436
437 // Add them
438 fListHist->Add(fLinearVdriftTest);
439 fListHist->Add(fPH2dTest);
440 fListHist->Add(fCH2dTest);
441 fListHist->Add(fPH2dSM);
442 fListHist->Add(fCH2dSM);
443 fListHist->Add(fPH2dSum);
444 fListHist->Add(fCH2dSum);
445
446 }
447
448 /////////////////////////////////////////
449 // Second debug level
450 ///////////////////////////////////////
451 if(fDebug > 1) {
452
453 fNbGoodTracks = new TH2F(Form("NbGoodTracks_%s",(const char*)fName),"NbGoodTracks",500,0.0,2500.0,200,0.0,100.0);
454 fNbGoodTracks->SetXTitle("Nb of good tracks");
455 fNbGoodTracks->SetYTitle("Centrality");
456 fNbGoodTracks->SetStats(0);
457
458 fNbTRDTrack = new TH1F(Form("TRDTrack_%s",(const char*)fName),"TRDTrack",50,0,50);
459 fNbTRDTrack->Sumw2();
460 fNbTRDTrackOffline = new TH1F(Form("TRDTrackOffline_%s",(const char*)fName),"TRDTrackOffline",50,0,50);
461 fNbTRDTrackOffline->Sumw2();
462 fNbTRDTrackStandalone = new TH1F(Form("TRDTrackStandalone_%s",(const char*)fName),"TRDTrackStandalone",50,0,50);
463 fNbTRDTrackStandalone->Sumw2();
464 fNbTPCTRDtrack = new TH2F(Form("NbTPCTRDtrack_%s",(const char*)fName),"NbTPCTRDtrack",100,0,100,100,0,100);
465 fNbTPCTRDtrack->Sumw2();
466 //
467 fNbTimeBin = new TH1F(Form("NbTimeBin_%s",(const char*)fName),"NbTimeBin",35,0,35);
468 fNbTimeBin->Sumw2();
469 fNbTimeBinOffline = new TH1F(Form("NbTimeBinOffline_%s",(const char*)fName),"NbTimeBinOffline",35,0,35);
470 fNbTimeBinOffline->Sumw2();
471 fNbTimeBinStandalone = new TH1F(Form("NbTimeBinStandalone_%s",(const char*)fName),"NbTimeBinStandalone",35,0,35);
472 fNbTimeBinStandalone->Sumw2();
473 //
474 fNbClusters = new TH1F(Form("NbClusters_%s",(const char*)fName),"",35,0,35);
475 fNbClusters->Sumw2();
476 fNbClustersOffline = new TH1F(Form("NbClustersOffline_%s",(const char*)fName),"",35,0,35);
477 fNbClustersOffline->Sumw2();
478 fNbClustersStandalone = new TH1F(Form("NbClustersStandalone_%s",(const char*)fName),"",35,0,35);
479 fNbClustersStandalone->Sumw2();
480 //
481 fNbTracklets = new TH1F(Form("NbTracklets_%s",(const char*)fName),"NbTracklets",540,0.,540.);
482 fNbTracklets->Sumw2();
483 fNbTrackletsOffline = new TH1F(Form("NbTrackletsOffline_%s",(const char*)fName),"NbTrackletsOffline",540,0.,540.);
484 fNbTrackletsOffline->Sumw2();
485 fNbTrackletsStandalone = new TH1F(Form("NbTrackletsStandalone_%s",(const char*)fName),"NbTrackletsStandalone",540,0.,540.);
486 fNbTrackletsStandalone->Sumw2();
487
488 fListHist->Add(fNbGoodTracks);
489
490 fListHist->Add(fNbTRDTrack);
491 fListHist->Add(fNbTRDTrackOffline);
492 fListHist->Add(fNbTRDTrackStandalone);
493 fListHist->Add(fNbTPCTRDtrack);
494
495 fListHist->Add(fNbTimeBin);
496 fListHist->Add(fNbTimeBinOffline);
497 fListHist->Add(fNbTimeBinStandalone);
498 fListHist->Add(fNbClusters);
499 fListHist->Add(fNbClustersOffline);
500 fListHist->Add(fNbClustersStandalone);
501 fListHist->Add(fNbTracklets);
502 fListHist->Add(fNbTrackletsOffline);
503 fListHist->Add(fNbTrackletsStandalone);
504
505 }
506
507 delete [] binLimLogPt;
508 delete [] binLimPt;
509
510 PostData(1,fListHist);
511
512 //cout << "AliTRDCalibTask::UserCreateOutputObjects() OUT" << endl;
513
514}
515
516//________________________________________________________________________
517void AliTRDCalibTask::UserExec(Option_t *)
518{
519 //
520 // Filling of the histos
521 //
522 //cout << "AliTRDCalibTask::Exec() IN" << endl;
523
524 // Init Versions and subversions used
525 if((fFirstRunGain==-1) || (fVersionGainUsed==-1) || (fSubVersionGainUsed==-1) || (fFirstRunGainLocal==-1) || (fVersionGainLocalUsed==-1) || (fSubVersionGainLocalUsed==-1) || (fFirstRunVdrift==-1) || (fVersionVdriftUsed==-1) || (fSubVersionVdriftUsed==-1)) {
526 if(!SetVersionSubversion()) {
527 PostData(1, fListHist);
528 return;
529 }
530 }
531 if(fCounter==0) {
532 if(fOnInstance) {
533 fTRDCalibraFillHisto->SetFirstRunGain(fFirstRunGain); // Gain Used
534 fTRDCalibraFillHisto->SetVersionGainUsed(fVersionGainUsed); // Gain Used
535 fTRDCalibraFillHisto->SetSubVersionGainUsed(fSubVersionGainUsed); // Gain Used
536 fTRDCalibraFillHisto->SetFirstRunGainLocal(fFirstRunGainLocal); // Gain Used
537 fTRDCalibraFillHisto->SetVersionGainLocalUsed(fVersionGainLocalUsed); // Gain Used
538 fTRDCalibraFillHisto->SetSubVersionGainLocalUsed(fSubVersionGainLocalUsed); // Gain Used
539 fTRDCalibraFillHisto->SetFirstRunVdrift(fFirstRunVdrift); // Vdrift Used
540 fTRDCalibraFillHisto->SetVersionVdriftUsed(fVersionVdriftUsed); // Vdrift Used
541 fTRDCalibraFillHisto->SetSubVersionVdriftUsed(fSubVersionVdriftUsed); // Vdrift Used
542 if((fFirstRunExB != -1) && (fVersionExBUsed != -1) && (fSubVersionExBUsed != -1)){
543 fTRDCalibraFillHisto->SetFirstRunExB(fFirstRunExB); // ExB Used
544 fTRDCalibraFillHisto->SetVersionExBUsed(fVersionExBUsed); // ExB Used
545 fTRDCalibraFillHisto->SetSubVersionExBUsed(fSubVersionExBUsed); // ExB Used
546 }
547 fTRDCalibraFillHisto->InitCalDet();
548 }
549 if(fDebug > 1){
550 // title CH2dTest
551 TString name("Ver");
552 name += fVersionGainUsed;
553 name += "Subver";
554 name += fSubVersionGainUsed;
555 name += "FirstRun";
556 name += fFirstRunGain;
557 name += "Nz0Nrphi0";
558 fCH2dTest->SetTitle(name);
559 // title PH2dTest
560 TString namee("Ver");
561 namee += fVersionVdriftUsed;
562 namee += "Subver";
563 namee += fSubVersionVdriftUsed;
564 namee += "FirstRun";
565 namee += fFirstRunVdrift;
566 namee += "Nz0Nrphi0";
567 fPH2dTest->SetTitle(namee);
568 }
569 }
570
571 // AliLog::SetGlobalLogLevel(AliLog::kError);
572 // cout << "AliTRDCalibTask::Exec() 1" << endl;
573 fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
574 if(!fESD){
575 AliError("ESD Event missing");
576 PostData(1, fListHist);
577 return;
578 }
579
580 const char* type = fESD->GetBeamType();
581
582
583 //printf("Counter %d\n",fCounter);
584
585 fCounter++;
586 fNEventsInput->Fill(1);
587
588 //cout << "maxEvent = " << fMaxEvent << endl;
589 //if(fCounter%100==0) cout << "fCounter = " << fCounter << endl;
590 if((fMaxEvent != 0) && (fMaxEvent < fCounter)) {
591 PostData(1, fListHist);
592 return;
593 }
594 //if(fCounter%100==0) cout << "fCounter1 = " << fCounter << endl;
595 //cout << "event = " << fCounter << endl;
596
597 //printf("Counter %d\n",fCounter);
598
599 ///////////////////
600 // Check trigger
601 ///////////////////
602 Bool_t pass = kTRUE;
603
604 if (fSelectTrigger) {
605
606 //printf("Will check the triggers\n");
607
608 Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
609 //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
610 if(fRejected) {
611 pass = kTRUE;
612 for(Int_t k = 0; k < numberOfTriggerSelected; k++){
613 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
614 const TString tString=obString->GetString();
615 if(fESD->IsTriggerClassFired((const char*)tString)) {
616 pass = kFALSE;
617 }
618 }
619 }
620 else {
621 pass = kFALSE;
622 for(Int_t k = 0; k < numberOfTriggerSelected; k++){
623 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
624 const TString tString=obString->GetString();
625 if(fESD->IsTriggerClassFired((const char*)tString)) {
626 pass = kTRUE;
627 }
628 }
629 }
630 if(!pass) {
631 PostData(1, fListHist);
632 return;
633 }
634
635 }
636
637 //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
638 //printf("Trigger passed\n");
639
640 ///////////////////////////////
641 // Require a primary vertex
642 //////////////////////////////
643 if(fRequirePrimaryVertex) {
644 const AliESDVertex* vtxESD = 0x0;
645 if (fVtxTPC) vtxESD = fESD->GetPrimaryVertexTPC() ;
646 else if (fVtxSPD) vtxESD = fESD->GetPrimaryVertexSPD() ;
647 else vtxESD = fESD->GetPrimaryVertexTracks() ;
648 if(!vtxESD){
649 PostData(1, fListHist);
650 return;
651 }
652 Int_t nCtrb = vtxESD->GetNContributors();
653 if(nCtrb < fMinNbContributors) {
654 PostData(1, fListHist);
655 return;
656 }
657 Double_t zPosition = vtxESD->GetZ();
658 if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) {
659 PostData(1, fListHist);
660 return;
661 }
662
663 }
664
665 //printf("Primary vertex passed\n");
666
667 //////////////////////////////////
668 // Reject pile-up with SPD
669 //////////////////////////////////
670
671 if(fRejectPileUpWithSPD) {
672 if(fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5)){
673 //printf("test\n");
674 PostData(1, fListHist);
675 return;
676 }
677 }
678
679 //////////////////////////////////////
680 // Requirement on number of good tracks
681 //////////////////////////////////////
682 Int_t nGoodParticles = 0;
683 Double_t nbTracks = fESD->GetNumberOfTracks();
684 for(Int_t itrack = 0; itrack < nbTracks; itrack++) {
685 if(ParticleGood(itrack)) nGoodParticles++;
686 }
687 if(fDebug > 1) {
688 // Centrality
689 AliCentrality *esdCentrality = fESD->GetCentrality();
690 Float_t centrality = esdCentrality->GetCentralityPercentile("V0M");
691 //Float_t centralityb = esdCentrality->GetCentralityPercentile("CL1");
692 fNbGoodTracks->Fill(nGoodParticles,centrality);
693 //printf("centrality %f, centralityb %f\n",centrality,centralityb);
694 }
695
696 //printf("Beam type %s\n",(const char*)type);
697 if (!strstr(type,"p-p")) {
698 //if (strstr(type,"A-A")) {
699 //printf("Will check the number of good tracks %f %f %f\n",nGoodParticles,fMinNbTracks,fMaxNbTracks);
700 if((nGoodParticles < fMinNbTracks) || (nGoodParticles > fMaxNbTracks)) {
701 PostData(1, fListHist);
702 return;
703 }
704 }
705
706 fNEvents->Fill(1);
707
708 // In total
709 Int_t nbTrdTracks = 0;
710 // standalone
711 Int_t nbTrdTracksStandalone = 0;
712 // offline
713 Int_t nbTrdTracksOffline = 0;
714 // TPC
715 Int_t nbtrackTPC = 0;
716
717
718
719 if (nbTracks <= 0.0) {
720
721 if(fDebug > 1) {
722 fNbTRDTrack->Fill(nbTrdTracks);
723 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
724 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
725 }
726 PostData(1, fListHist);
727 return;
728 }
729
730
731 fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
732 if(!fESDfriend){
733 AliError("fESDfriend not available");
734 PostData(1, fListHist);
735 return;
736 }
737
738 if(fESDfriend->TestSkipBit()) {
739 PostData(1, fListHist);
740 return;
741 }
742
743 //printf("has friends\n");
744
745 /////////////////////////////////////
746 // Loop on AliESDtrack
747 ////////////////////////////////////
748 //printf("Nb of tracks %f\n",nbTracks);
749 Int_t nbTracksfriends = fESDfriend->GetNumberOfTracks();
750 for(int itrk=0; itrk < nbTracksfriends; ++itrk){
751
752 // Get ESD track
753 fkEsdTrack = fESD->GetTrack(itrk);
754 if(!fkEsdTrack) continue;
755 ULong_t status = fkEsdTrack->GetStatus();
756 if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
757
758 // Fix suggested by Alex
759 fFriendTrack = fESDfriend->GetTrack(itrk);
760 //printf("itrk %d\n",itrk);
761 //fFriendTrack = (fESDfriend->GetNumberOfTracks()>itrk)?fESDfriend->GetTrack(itrk):NULL;
762 if(!fFriendTrack) {
763 //printf("No friend track %d\n",itrk);
764 continue;
765 }
766
767 // Other cuts
768 fTrdTrack = 0x0;
769 Bool_t good = kTRUE;
770 Bool_t standalonetrack = kFALSE;
771 Bool_t offlinetrack = kFALSE;
772 //ULong_t status = fkEsdTrack->GetStatus();
773
774 //////////////////////////////////////
775 // Loop on calibration objects
776 //////////////////////////////////////
777 Int_t icalib=0;
778 Int_t nTRDtrackV1=0;
779 while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
780 //printf("Name %s\n",fCalibObject->IsA()->GetName());
781 if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
782 //printf("Find the calibration object\n");
783 ++nTRDtrackV1;
784
785 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
786 standalonetrack = kTRUE;
787 }
788 if((status&(AliESDtrack::kTRDin))) {
789 offlinetrack = kTRUE;
790 }
791 if(fOfflineTracks){
792 if(!offlinetrack){
793 good = kFALSE;
794 }
795 }
796 else if(fStandaloneTracks){
797 if(!standalonetrack){
798 good = kFALSE;
799 }
800 }
801
802 fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
803 // process chamberstatus
804 fTRDChamberStatus->ProcessTrack(fTrdTrack);
805 }
806
807 // Quality cuts on the AliESDtrack
808 if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
809 //printf("Not a good track\n");
810 continue;
811 }
812
813 // First Absolute gain calibration
814 Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
815 Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID();
816 //printf("Number of trd tracklets %d and PID trd tracklets %d\n",trdNTracklets,trdNTrackletsPID);
817 if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
818 for(Int_t iPlane = 0; iPlane < 6; ++iPlane){
819 //Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
820 //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
821 //Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
822 //printf("momentum %f, slide %f\n",momentum,slide);
823 if(fkEsdTrack->GetTRDslice(iPlane) > 0.0)
824 fAbsoluteGain->Fill(fkEsdTrack->GetTRDslice(iPlane)*8.0/100.0,
825 fkEsdTrack->GetTRDmomentum(iPlane));
826 }
827 }
828
829
830 if(!fTrdTrack) continue;
831
832 if(good && fOnInstance) {
833 //cout << "good" << endl;
834 fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack,fkEsdTrack);
835 //printf("Fill fTRDCalibraFillHisto\n");
836 }
837
838 if(IsPHQon()){
839 const Double_t mag = AliTRDdEdxBaseUtils::IsExBOn() ? fESD->GetMagneticField() : -1;
840 const Int_t charge = AliTRDdEdxBaseUtils::IsExBOn() ? fkEsdTrack->Charge() : -1;
841 const Double_t toTPCscale = AliTRDdEdxCalibUtils::GetCalibTPCscale(fkEsdTrack->GetTPCncls(), fkEsdTrack->GetTPCsignal());
842 if(toTPCscale>0){
843 AliTRDdEdxCalibUtils::FillHist(fTrdTrack, 0, mag, charge, toTPCscale);
844 }
845 }
846
847 //////////////////////////////////
848 // Debug
849 ////////////////////////////////
850
851 if(fDebug > 0) {
852
853 //printf("Enter debug\n");
854
855 Int_t nbtracklets = 0;
856
857 // Check some stuff
858 Bool_t standalonetracklet = kFALSE;
859 const AliTRDseedV1 *tracklet = 0x0;
860 //////////////////////////////////////
861 // Loop tracklets
862 /////////////////////////////////////
863 Int_t nbclusters=0;
864 Double_t phtb[AliTRDseedV1::kNtb];
865 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
866 Double_t sum = 0.0;
867 Float_t normalisation = 1.13;
868 Int_t detector = 0;
869 Int_t sector = 0;
870 for(Int_t itr = 0; itr < 6; ++itr){
871
872 if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
873 if(!tracklet->IsOK()) continue;
874 ++nbtracklets;
875 standalonetracklet = kFALSE;
876 if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
877
878 nbclusters = 0;
879 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
880 sum = 0.0;
881 detector = 0;
882 sector = 0;
883 //Int_t crossrow = 0;
884
885 // Check no shared clusters
886 //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
887 // if((fcl = tracklet->GetClusters(icc))) crossrow = 1;
888 // }
889
890 // Loop on clusters
891 Int_t time = 0;
892 Float_t ch = 0;
893 Float_t qcl = 0;
894 for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
895
896 if(!(fCl = tracklet->GetClusters(ic))) continue;
897 ++nbclusters;
898 time = fCl->GetPadTime();
899 //ch = tracklet->GetdQdl(ic);
900 ch = tracklet->GetQperTB(ic);
901 qcl = TMath::Abs(fCl->GetQ());
902 detector = fCl->GetDetector();
903 // Add the charge if shared cluster
904 if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
905 if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
906 qcl += TMath::Abs(fCl->GetQ());
907 //printf("Add the cluster charge\n");
908 }
909 }
910 if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
911 if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;
912 else sum += ch/normalisation;
913
914 if(fDebug > 1) {
915 fNbTimeBin->Fill(time);
916 if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
917 else fNbTimeBinOffline->Fill(time);
918 }
919 }
920 sector = AliTRDgeometry::GetSector(detector);
921
922 if(fDebug > 1) {
923 fNbTracklets->Fill(detector);
924 if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
925 else fNbTrackletsOffline->Fill(detector);
926
927 fNbClusters->Fill(nbclusters);
928 if(tracklet->IsStandAlone()) fNbClustersStandalone->Fill(nbclusters);
929 else fNbClustersOffline->Fill(nbclusters);
930 }
931
932 if((nbclusters > fLow) && (nbclusters < fHigh)){
933 if(fRelativeScale > 0.0) sum = sum/fRelativeScale;
934 fCH2dTest->Fill(sum,detector+0.5);
935 fCH2dSM->Fill(sum,sector+0.5);
936 fCH2dSum->Fill(sum,0.5);
937 Bool_t checknoise = kTRUE;
938 if(fMaxCluster > 0) {
939 if(phtb[0] > fMaxCluster) checknoise = kFALSE;
940 if(fNbTimeBins > fNbMaxCluster) {
941 for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
942 if(phtb[k] > fMaxCluster) checknoise = kFALSE;
943 }
944 }
945 }
946 if(checknoise) {
947 for(int ic=0; ic<fNbTimeBins; ic++){
948 if(fFillZero) {
949 fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
950 fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
951 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
952 }
953 else {
954 if(phtb[ic] > 0.0) {
955 fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
956 fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
957 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
958 }
959 }
960 }
961 }
962 }
963 if(detector == 0) FindP1TrackPHtrackletV1Test(tracklet,nbclusters);
964
965 } // loop on tracklets
966
967
968 } // debug
969
970 if(nTRDtrackV1 > 0) {
971 ++nbTrdTracks;
972 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
973 ++nbTrdTracksStandalone;
974 }
975 if((status&(AliESDtrack::kTRDin))) {
976 ++nbTrdTracksOffline;
977 }
978 }
979 //delete fFriendTrack;
980 } // loop ESD track
981
982 if(fDebug > 1) {
983 fNbTRDTrack->Fill(nbTrdTracks);
984 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
985 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
986 fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
987 }
988
989 // Post output data
990 PostData(1, fListHist);
991 //cout << "AliTRDCalibTask::Exec() OUT" << endl;
992}
993
994//________________________________________________________________________
995void AliTRDCalibTask::Terminate(Option_t *)
996{
997 //
998 // Terminate
999 //
1000
1001 if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
1002
1003
1004}
1005//_______________________________________________________
1006Bool_t AliTRDCalibTask::Load(const Char_t *filename)
1007{
1008 //
1009 // Generic container loader
1010 //
1011
1012 if(!TFile::Open(filename)){
1013 //AliWarning(Form("Couldn't open file %s.", filename));
1014 return kFALSE;
1015 }
1016 TList *o = 0x0;
1017 if(!(o = (TList*)gFile->Get(GetName()))){
1018 //AliWarning("Missing histogram container.");
1019 return kFALSE;
1020 }
1021 fListHist = (TList*)o->Clone(GetName());
1022 gFile->Close();
1023 return kTRUE;
1024}
1025//_______________________________________________________
1026Bool_t AliTRDCalibTask::Load(TList *lister)
1027{
1028 //
1029 // Generic container loader
1030 //
1031
1032 fListHist = (TList*)lister->Clone(GetName());
1033 return kTRUE;
1034}
1035//_______________________________________________________________________________________
1036void AliTRDCalibTask::AddTask(const AliTRDCalibTask * calibTask) {
1037
1038 //
1039 // Add stats
1040 //
1041
1042 TList *listcalibTask = calibTask->GetList();
1043 if(!listcalibTask) return;
1044
1045 THnSparseI *histoEntries = (THnSparseI *) listcalibTask->FindObject("NumberOfEntries");
1046
1047 TH1I *nEvents = (TH1I *) listcalibTask->FindObject(Form("NEvents_%s",(const char*)calibTask->GetName()));
1048 TH1I *nEventsInput = (TH1I *) listcalibTask->FindObject(Form("NEventsInput_%s",(const char*)calibTask->GetName()));
1049 TH2F *absoluteGain = (TH2F *) listcalibTask->FindObject(Form("AbsoluteGain_%s",(const char*)calibTask->GetName()));
1050
1051 TH1F *trdTrack = (TH1F *) listcalibTask->FindObject(Form("TRDTrack_%s",(const char*)calibTask->GetName()));
1052 TH1F *trdTrackOffline = (TH1F *) listcalibTask->FindObject(Form("TRDTrackOffline_%s",(const char*)calibTask->GetName()));
1053 TH1F *trdTrackStandalone = (TH1F *) listcalibTask->FindObject(Form("TRDTrackStandalone_%s",(const char*)calibTask->GetName()));
1054
1055 TH2F *tpctrdTrack = (TH2F *) listcalibTask->FindObject(Form("NbTPCTRDtrack_%s",(const char*)calibTask->GetName()));
1056
1057 TH1F *nbTimeBin = (TH1F *) listcalibTask->FindObject(Form("NbTimeBin_%s",(const char*)calibTask->GetName()));
1058 TH1F *nbTimeBinOffline = (TH1F *) listcalibTask->FindObject(Form("NbTimeBinOffline_%s",(const char*)calibTask->GetName()));
1059 TH1F *nbTimeBinStandalone = (TH1F *) listcalibTask->FindObject(Form("NbTimeBinStandalone_%s",(const char*)calibTask->GetName()));
1060
1061 TH1F *nbClusters = (TH1F *) listcalibTask->FindObject(Form("NbClusters_%s",(const char*)calibTask->GetName()));
1062 TH1F *nbClustersOffline = (TH1F *) listcalibTask->FindObject(Form("NbClustersOffline_%s",(const char*)calibTask->GetName()));
1063 TH1F *nbClustersStandalone = (TH1F *) listcalibTask->FindObject(Form("NbClustersStandalone_%s",(const char*)calibTask->GetName()));
1064
1065 TH1F *nbTracklets = (TH1F *) listcalibTask->FindObject(Form("NbTracklets_%s",(const char*)calibTask->GetName()));
1066 TH1F *nbTrackletsOffline = (TH1F *) listcalibTask->FindObject(Form("NbTrackletsOffline_%s",(const char*)calibTask->GetName()));
1067 TH1F *nbTrackletsStandalone = (TH1F *) listcalibTask->FindObject(Form("NbTrackletsStandalone_%s",(const char*)calibTask->GetName()));
1068
1069 TH2I *ch2d = (TH2I *) listcalibTask->FindObject("CH2d");
1070 TProfile2D *ph2d = (TProfile2D *) listcalibTask->FindObject("PH2d");
1071 TProfile2D *prf2d = (TProfile2D *) listcalibTask->FindObject("PRF2d");
1072
1073 TH2I *ch2dSum = (TH2I *) listcalibTask->FindObject(Form("CH2dSum_%s",(const char*)calibTask->GetName()));
1074 TProfile2D *ph2dSum = (TProfile2D *) listcalibTask->FindObject(Form("PH2dSum_%s",(const char*)calibTask->GetName()));
1075
1076 TH2I *ch2dSM = (TH2I *) listcalibTask->FindObject(Form("CH2dSM_%s",(const char*)calibTask->GetName()));
1077 TProfile2D *ph2dSM = (TProfile2D *) listcalibTask->FindObject(Form("PH2dSM_%s",(const char*)calibTask->GetName()));
1078
1079 AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) listcalibTask->FindObject("AliTRDCalibraVdriftLinearFit");
1080 AliTRDCalibraExbAltFit *exbaltfit = (AliTRDCalibraExbAltFit *) listcalibTask->FindObject("AliTRDCalibraExbAltFit");
1081 AliTRDCalibraVector *calibraVector = (AliTRDCalibraVector *) listcalibTask->FindObject("AliTRDCalibraVector");
1082
1083 //
1084
1085 THnSparseI *inhistoEntries = (THnSparseI *) fListHist->FindObject("NumberOfEntries");
1086
1087 TH1I *inEventsInput = (TH1I *) fListHist->FindObject(Form("NEventsInput_%s",(const char*)fName));
1088 TH1I *inEvents = (TH1I *) fListHist->FindObject(Form("NEvents_%s",(const char*)fName));
1089 TH2F *iabsoluteGain = (TH2F *) fListHist->FindObject(Form("AbsoluteGain_%s",(const char*)fName));
1090
1091 TH1F *itrdTrack = (TH1F *) fListHist->FindObject(Form("TRDTrack_%s",(const char*)fName));
1092 TH1F *itrdTrackOffline = (TH1F *) fListHist->FindObject(Form("TRDTrackOffline_%s",(const char*)fName));
1093 TH1F *itrdTrackStandalone = (TH1F *) fListHist->FindObject(Form("TRDTrackStandalone_%s",(const char*)fName));
1094
1095 TH2F *itpctrdTrack = (TH2F *) fListHist->FindObject(Form("NbTPCTRDtrack_%s",(const char*)fName));
1096
1097 TH1F *inbTimeBin = (TH1F *) fListHist->FindObject(Form("NbTimeBin_%s",(const char*)fName));
1098 TH1F *inbTimeBinOffline = (TH1F *) fListHist->FindObject(Form("NbTimeBinOffline_%s",(const char*)fName));
1099 TH1F *inbTimeBinStandalone = (TH1F *) fListHist->FindObject(Form("NbTimeBinStandalone_%s",(const char*)fName));
1100
1101 TH1F *inbClusters = (TH1F *) fListHist->FindObject(Form("NbClusters_%s",(const char*)fName));
1102 TH1F *inbClustersOffline = (TH1F *) fListHist->FindObject(Form("NbClustersOffline_%s",(const char*)fName));
1103 TH1F *inbClustersStandalone = (TH1F *) fListHist->FindObject(Form("NbClustersStandalone_%s",(const char*)fName));
1104
1105 TH1F *inbTracklets = (TH1F *) fListHist->FindObject(Form("NbTracklets_%s",(const char*)fName));
1106 TH1F *inbTrackletsOffline = (TH1F *) fListHist->FindObject(Form("NbTrackletsOffline_%s",(const char*)fName));
1107 TH1F *inbTrackletsStandalone = (TH1F *) fListHist->FindObject(Form("NbTrackletsStandalone_%s",(const char*)fName));
1108
1109 TH2I *ich2d = (TH2I *) fListHist->FindObject("CH2d");
1110 TProfile2D *iph2d = (TProfile2D *) fListHist->FindObject("PH2d");
1111 TProfile2D *iprf2d = (TProfile2D *) fListHist->FindObject("PRF2d");
1112
1113 TH2I *ich2dSum = (TH2I *) fListHist->FindObject(Form("CH2dSum_%s",(const char*)fName));
1114 TProfile2D *iph2dSum = (TProfile2D *) fListHist->FindObject(Form("PH2dSum_%s",(const char*)fName));
1115
1116 TH2I *ich2dSM = (TH2I *) fListHist->FindObject(Form("CH2dSM_%s",(const char*)fName));
1117 TProfile2D *iph2dSM = (TProfile2D *) fListHist->FindObject(Form("PH2dSM_%s",(const char*)fName));
1118
1119 AliTRDCalibraVdriftLinearFit *ilinearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
1120 AliTRDCalibraExbAltFit *iexbaltfit = (AliTRDCalibraExbAltFit *) fListHist->FindObject("AliTRDCalibraExbAltFit");
1121 AliTRDCalibraVector *icalibraVector = (AliTRDCalibraVector *) fListHist->FindObject("AliTRDCalibraVector");
1122
1123
1124 // Add
1125
1126 if(histoEntries) {
1127 if(inhistoEntries) {
1128 inhistoEntries->Add(histoEntries);
1129 //printf("Add Events\n");
1130 }
1131 else {
1132 //printf("Create new Events\n");
1133 inhistoEntries = (THnSparseI *) histoEntries->Clone();
1134 fListHist->Add(inhistoEntries);
1135 }
1136 }
1137
1138 if(nEventsInput) {
1139 if(inEventsInput) {
1140 inEventsInput->Add(nEventsInput);
1141 //printf("Add Events\n");
1142 }
1143 else {
1144 //printf("Create new Events\n");
1145 inEventsInput = new TH1I(*nEventsInput);
1146 fListHist->Add(inEventsInput);
1147 }
1148 }
1149
1150 if(nEvents) {
1151 if(inEvents) {
1152 inEvents->Add(nEvents);
1153 //printf("Add Events\n");
1154 }
1155 else {
1156 //printf("Create new Events\n");
1157 inEvents = new TH1I(*nEvents);
1158 fListHist->Add(inEvents);
1159 }
1160 }
1161
1162 if(absoluteGain) {
1163 if(iabsoluteGain) iabsoluteGain->Add(absoluteGain);
1164 else {
1165 iabsoluteGain = new TH2F(*absoluteGain);
1166 fListHist->Add(iabsoluteGain);
1167 }
1168 }
1169
1170 if(trdTrack) {
1171 if(itrdTrack) itrdTrack->Add(trdTrack);
1172 else {
1173 itrdTrack = new TH1F(*trdTrack);
1174 fListHist->Add(itrdTrack);
1175 }
1176 }
1177
1178 if(trdTrackOffline) {
1179 if(itrdTrackOffline) itrdTrackOffline->Add(trdTrackOffline);
1180 else {
1181 itrdTrackOffline = new TH1F(*trdTrackOffline);
1182 fListHist->Add(itrdTrackOffline);
1183 }
1184 }
1185
1186 if(trdTrackStandalone) {
1187 if(itrdTrackStandalone) itrdTrackStandalone->Add(trdTrackStandalone);
1188 else {
1189 itrdTrackStandalone = new TH1F(*trdTrackStandalone);
1190 fListHist->Add(itrdTrackStandalone);
1191 }
1192 }
1193
1194 if(tpctrdTrack) {
1195 if(itpctrdTrack) itpctrdTrack->Add(tpctrdTrack);
1196 else {
1197 itpctrdTrack = new TH2F(*tpctrdTrack);
1198 fListHist->Add(itpctrdTrack);
1199 }
1200 }
1201
1202 if(nbTimeBin) {
1203 if(inbTimeBin) inbTimeBin->Add(nbTimeBin);
1204 else {
1205 inbTimeBin = new TH1F(*inbTimeBin);
1206 fListHist->Add(inbTimeBin);
1207 }
1208 }
1209
1210 if(nbTimeBinOffline) {
1211 if(inbTimeBinOffline) inbTimeBinOffline->Add(nbTimeBinOffline);
1212 else {
1213 inbTimeBinOffline = new TH1F(*nbTimeBinOffline);
1214 fListHist->Add(inbTimeBinOffline);
1215 }
1216 }
1217
1218 if(nbTimeBinStandalone) {
1219 if(inbTimeBinStandalone) inbTimeBinStandalone->Add(nbTimeBinStandalone);
1220 else {
1221 inbTimeBinStandalone = new TH1F(*nbTimeBinStandalone);
1222 fListHist->Add(inbTimeBinStandalone);
1223 }
1224 }
1225
1226 if(nbClusters) {
1227 if(inbClusters) inbClusters->Add(nbClusters);
1228 else {
1229 inbClusters = new TH1F(*nbClusters);
1230 fListHist->Add(inbClusters);
1231 }
1232 }
1233
1234 if(nbClustersOffline) {
1235 if(inbClustersOffline) inbClustersOffline->Add(nbClustersOffline);
1236 else {
1237 inbClustersOffline = new TH1F(*nbClustersOffline);
1238 fListHist->Add(inbClustersOffline);
1239 }
1240 }
1241
1242 if(nbClustersStandalone) {
1243 if(inbClustersStandalone) inbClustersStandalone->Add(nbClustersStandalone);
1244 else {
1245 inbClustersStandalone = new TH1F(*nbClustersStandalone);
1246 fListHist->Add(inbClustersStandalone);
1247 }
1248 }
1249
1250 if(nbTracklets) {
1251 if(inbTracklets) inbTracklets->Add(nbTracklets);
1252 else {
1253 inbTracklets = new TH1F(*nbTracklets);
1254 fListHist->Add(inbTracklets);
1255 }
1256 }
1257
1258 if(nbTrackletsOffline) {
1259 if(inbTrackletsOffline) inbTrackletsOffline->Add(nbTrackletsOffline);
1260 else {
1261 inbTrackletsOffline = new TH1F(*nbTrackletsOffline);
1262 fListHist->Add(inbTrackletsOffline);
1263 }
1264 }
1265
1266 if(nbTrackletsStandalone) {
1267 if(inbTrackletsStandalone) inbTrackletsStandalone->Add(nbTrackletsStandalone);
1268 else {
1269 inbTrackletsStandalone = new TH1F(*nbTrackletsStandalone);
1270 fListHist->Add(inbTrackletsStandalone);
1271 }
1272 }
1273
1274 if(ch2d) {
1275 if(ich2d) ich2d->Add(ch2d);
1276 else {
1277 ich2d = new TH2I(*ch2d);
1278 fListHist->Add(ich2d);
1279 }
1280 }
1281
1282 if(ph2d) {
1283 if(iph2d) iph2d->Add(ph2d);
1284 else {
1285 iph2d = new TProfile2D(*ph2d);
1286 fListHist->Add(iph2d);
1287 }
1288 }
1289
1290 if(prf2d) {
1291 if(iprf2d) iprf2d->Add(prf2d);
1292 else {
1293 iprf2d = new TProfile2D(*prf2d);
1294 fListHist->Add(iprf2d);
1295 }
1296 }
1297
1298 if(ch2dSum) {
1299 if(ich2dSum) ich2dSum->Add(ch2dSum);
1300 else {
1301 ich2dSum = new TH2I(*ch2dSum);
1302 fListHist->Add(ich2dSum);
1303 }
1304 }
1305
1306 if(ph2dSum) {
1307 if(iph2dSum) iph2dSum->Add(ph2dSum);
1308 else {
1309 iph2dSum = new TProfile2D(*ph2dSum);
1310 fListHist->Add(iph2dSum);
1311 }
1312 }
1313
1314 if(ch2dSM) {
1315 if(ich2dSM) ich2dSM->Add(ch2dSM);
1316 else {
1317 ich2dSM = new TH2I(*ch2dSM);
1318 fListHist->Add(ich2dSM);
1319 }
1320 }
1321
1322 if(ph2dSM) {
1323 if(iph2dSM) iph2dSM->Add(ph2dSM);
1324 else {
1325 iph2dSM = new TProfile2D(*ph2dSM);
1326 fListHist->Add(iph2dSM);
1327 }
1328 }
1329
1330 if(linearfit) {
1331 if(ilinearfit) ilinearfit->Add(linearfit);
1332 else {
1333 ilinearfit = new AliTRDCalibraVdriftLinearFit(*linearfit);
1334 fListHist->Add(ilinearfit);
1335 }
1336 }
1337
1338 if(exbaltfit) {
1339 if(iexbaltfit) iexbaltfit->Add(exbaltfit);
1340 else {
1341 iexbaltfit = new AliTRDCalibraExbAltFit(*exbaltfit);
1342 fListHist->Add(iexbaltfit);
1343 }
1344 }
1345
1346 if(calibraVector) {
1347 if(icalibraVector) icalibraVector->Add(calibraVector);
1348 else {
1349 icalibraVector = new AliTRDCalibraVector(*calibraVector);
1350 fListHist->Add(icalibraVector);
1351 }
1352 }
1353
1354}
1355//________________________________________________________________________________
1356Long64_t AliTRDCalibTask::Merge(TCollection *li) {
1357
1358 //
1359 // merge component
1360 //
1361
1362 TIterator* iter = li->MakeIterator();
1363 AliTRDCalibTask* cal = 0;
1364
1365 while ((cal = (AliTRDCalibTask*)iter->Next())) {
1366 if (!cal->InheritsFrom(AliTRDCalibTask::Class())) {
1367 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
1368 return -1;
1369 }
1370
1371 // add histograms here...
1372 this->AddTask(cal);
1373
1374 }
1375
1376 return 0;
1377
1378}
1379//_____________________________________________________
1380Bool_t AliTRDCalibTask::SetVersionSubversion(){
1381 //
1382 // Load Chamber Gain factors into the Tender supply
1383 //
1384
1385 //printf("SetVersionSubversion\n");
1386
1387 //find previous entry from the UserInfo
1388 TTree *tree=((TChain*)GetInputData(0))->GetTree();
1389 if (!tree) {
1390 AliError("Tree not found in ESDhandler");
1391 return kFALSE;
1392 }
1393
1394 TList *userInfo=(TList*)tree->GetUserInfo();
1395 if (!userInfo) {
1396 AliError("No UserInfo found in tree");
1397 return kFALSE;
1398 }
1399
1400 TList *cdbList=(TList*)userInfo->FindObject("cdbList");
1401 if (!cdbList) {
1402 AliError("No cdbList found in UserInfo");
1403 if (AliLog::GetGlobalLogLevel()>=AliLog::kError) userInfo->Print();
1404 return kFALSE;
1405 }
1406
1407 TIter nextCDB(cdbList);
1408 TObjString *os=0x0;
1409 while ( (os=(TObjString*)nextCDB()) ){
1410 if(os->GetString().Contains("TRD/Calib/ChamberGainFactor")){
1411 // Get Old gain calibration
1412 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1413 fFirstRunGain = id->GetFirstRun();
1414 fVersionGainUsed = id->GetVersion();
1415 fSubVersionGainUsed = id->GetSubVersion();
1416 } else if(os->GetString().Contains("TRD/Calib/ChamberVdrift")){
1417 // Get Old drift velocity calibration
1418 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1419 fFirstRunVdrift = id->GetFirstRun();
1420 fVersionVdriftUsed = id->GetVersion();
1421 fSubVersionVdriftUsed = id->GetSubVersion();
1422 } else if(os->GetString().Contains("TRD/Calib/LocalGainFactor")){
1423 // Get Old drift velocity calibration
1424 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1425 fFirstRunGainLocal = id->GetFirstRun();
1426 fVersionGainLocalUsed = id->GetVersion();
1427 fSubVersionGainLocalUsed = id->GetSubVersion();
1428 } else if((os->GetString().Contains("TRD/Calib/ChamberExB")) && (!os->GetString().Contains("TRD/Calib/ChamberExBAlt"))){
1429 // Get Old drift velocity calibration
1430 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1431 fFirstRunExB = id->GetFirstRun();
1432 fVersionExBUsed = id->GetVersion();
1433 fSubVersionExBUsed = id->GetSubVersion();
1434 //printf("Version %d and subversion %d\n",fVersionExBUsed,fSubVersionExBUsed);
1435 }
1436 }
1437
1438 //printf("VersionGain %d, SubversionGain %d, VersionLocalGain %d, Subversionlocalgain %d, Versionvdrift %d, Subversionvdrift %d\n",fVersionGainUsed,fSubVersionGainUsed,fVersionGainLocalUsed,fSubVersionGainLocalUsed,fVersionVdriftUsed,fSubVersionVdriftUsed);
1439
1440 // Check
1441 if((fFirstRunGain < 0) ||
1442 (fFirstRunGainLocal < 0) ||
1443 (fFirstRunVdrift < 0) ||
1444 (fVersionGainUsed < 0) ||
1445 (fVersionGainLocalUsed < 0) ||
1446 (fSubVersionGainUsed < 0) ||
1447 (fSubVersionGainLocalUsed < 0) ||
1448 (fVersionVdriftUsed < 0) ||
1449 (fSubVersionVdriftUsed < 0)) {
1450 AliError("No recent calibration found");
1451 return kFALSE;
1452 }
1453 else return kTRUE;
1454
1455}
1456//_________________________________________________________________________________________________________________________
1457Bool_t AliTRDCalibTask::ParticleGood(int i) const {
1458
1459 //
1460 // Definition of good tracks
1461 //
1462
1463
1464 AliESDtrack *track = fESD->GetTrack(i);
1465 if (!track->IsOn(AliESDtrack::kTPCrefit)) return 0; // TPC refit
1466 if (track->GetTPCNcls() < 90) return 0; // number of TPC clusters
1467 if (fabs(track->Eta())>0.8) return 0; // fiducial pseudorapidity
1468 Float_t r,z;
1469 track->GetImpactParametersTPC(r,z);
1470 if (fabs(z)>2.0) return 0; // impact parameter in z
1471 if (fabs(r)>2.0) return 0; // impact parameter in xy
1472 if (r==0) return 0;
1473 return 1;
1474
1475
1476}
1477//______________________________________________________________________________________________________________________
1478Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1479{
1480 //
1481 // Drift velocity calibration:
1482 // Fit the clusters with a straight line
1483 // From the slope find the drift velocity
1484 //
1485
1486 ////////////////////////////////////////////////
1487 //Number of points: if less than 3 return kFALSE
1488 /////////////////////////////////////////////////
1489 if(nbclusters <= 2) return kFALSE;
1490
1491 ////////////
1492 //Variables
1493 ////////////
1494 // results of the linear fit
1495 Double_t dydt = 0.0; // dydt tracklet after straight line fit
1496 Double_t errorpar = 0.0; // error after straight line fit on dy/dt
1497 Double_t pointError = 0.0; // error after straight line fit
1498 // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
1499 Int_t crossrow = 0; // if it crosses a pad row
1500 Int_t rowp = -1; // if it crosses a pad row
1501 Float_t tnt = tracklet->GetTilt(); // tan tiltingangle
1502 TLinearFitter linearFitterTracklet(2,"pol1");
1503 linearFitterTracklet.StoreData(kTRUE);
1504
1505
1506 ///////////////////////////////////////////
1507 // Take the parameters of the track
1508 //////////////////////////////////////////
1509 // take now the snp, tnp and tgl from the track
1510 Double_t snp = tracklet->GetSnp(); // sin dy/dx at the end of the chamber
1511 Double_t tnp = 0.0; // dy/dx at the end of the chamber
1512 if( TMath::Abs(snp) < 1.){
1513 tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1514 }
1515 Double_t tgl = tracklet->GetTgl(); // dz/dl
1516 Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp); // dz/dx calculated from dz/dl
1517 // at the entrance
1518 //Double_t tnp = tracklet->GetYref(1); // dy/dx at the entrance of the chamber
1519 //Double_t tgl = tracklet->GetZref(1); // dz/dl at the entrance of the chamber
1520 //Double_t dzdx = tgl; //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1521 // at the end with correction due to linear fit
1522 //Double_t tnp = tracklet->GetYfit(1); // dy/dx at the end of the chamber after fit correction
1523 //Double_t tgl = tracklet->GetZfit(1); // dz/dl at the end of the chamber after fit correction
1524
1525
1526 ////////////////////////////
1527 // loop over the clusters
1528 ////////////////////////////
1529 Int_t nbli = 0;
1530 AliTRDcluster *cl = 0x0;
1531 //////////////////////////////
1532 // Check no shared clusters
1533 //////////////////////////////
1534 for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
1535 cl = tracklet->GetClusters(icc);
1536 if(cl) crossrow = 1;
1537 }
1538 //////////////////////////////////
1539 // Loop clusters
1540 //////////////////////////////////
1541 for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1542 if(!(cl = tracklet->GetClusters(ic))) continue;
1543 //if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
1544
1545 Double_t ycluster = cl->GetY();
1546 Int_t time = cl->GetPadTime();
1547 Double_t timeis = time/10.0;
1548 //See if cross two pad rows
1549 Int_t row = cl->GetPadRow();
1550 if(rowp==-1) rowp = row;
1551 if(row != rowp) crossrow = 1;
1552
1553 linearFitterTracklet.AddPoint(&timeis,ycluster,1);
1554 nbli++;
1555
1556
1557 }
1558
1559 ////////////////////////////////////
1560 // Do the straight line fit now
1561 ///////////////////////////////////
1562 if(nbli <= 2){
1563 linearFitterTracklet.ClearPoints();
1564 return kFALSE;
1565 }
1566 TVectorD pars;
1567 linearFitterTracklet.Eval();
1568 linearFitterTracklet.GetParameters(pars);
1569 pointError = TMath::Sqrt(linearFitterTracklet.GetChisquare()/(nbli-2));
1570 errorpar = linearFitterTracklet.GetParError(1)*pointError;
1571 dydt = pars[1];
1572 //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",linearFitterTracklet->GetChisquare(),nbli,pointError,linearFitterTracklet->GetParError(1),errorpar);
1573 linearFitterTracklet.ClearPoints();
1574
1575 /////////////////////////
1576 // Cuts quality
1577 ////////////////////////
1578
1579 if(nbclusters < fLow) return kFALSE;
1580 if(nbclusters > fHigh) return kFALSE;
1581 if(pointError >= 0.3) return kFALSE;
1582 if(crossrow == 1) return kTRUE;
1583
1584 ///////////////////////
1585 // Fill
1586 //////////////////////
1587
1588 if(fDebug > 0){
1589 //Add to the linear fitter of the detector
1590 if( TMath::Abs(snp) < 1.){
1591 Double_t x = tnp-dzdx*tnt;
1592 //if(!fLinearVdriftTest) printf("Not there\n");
1593 Double_t nbentries = fLinearVdriftTest->GetEntries();
1594 if(nbentries < (5.0*32767)) fLinearVdriftTest->Fill(x,dydt);
1595 }
1596 }
1597
1598 return kTRUE;
1599}
1600