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 // Offline TRD calibration task
24 // R. Bailhache (R.Bailhache@gsi.de)
26 //////////////////////////////////////////////////////////////////////////////////////
30 #include "Riostream.h"
34 #include "TProfile2D.h"
41 #include "TObjArray.h"
47 #include "TIterator.h"
49 #include "AliAnalysisTask.h"
50 #include "AliAnalysisManager.h"
52 #include "AliExternalTrackParam.h"
53 #include "AliESDVertex.h"
54 #include "AliESDEvent.h"
55 #include "AliESDfriend.h"
56 #include "AliCentrality.h"
57 #include "AliESDInputHandler.h"
58 #include "AliESDtrack.h"
59 #include "AliESDfriendTrack.h"
60 #include "AliTRDtrackV1.h"
61 #include "AliTRDseedV1.h"
62 #include "AliTRDcluster.h"
63 #include "AliTRDgeometry.h"
64 #include "AliESDtrackCuts.h"
65 #include "AliESDVertex.h"
66 #include "AliTRDCalDet.h"
68 #include "AliTRDCalibraVector.h"
69 #include "AliTRDCalibraFillHisto.h"
70 #include "AliTRDCalibraVdriftLinearFit.h"
72 #include "AliTRDcalibDB.h"
77 #include "AliTRDCalibTask.h"
80 ClassImp(AliTRDCalibTask)
82 //________________________________________________________________________
83 AliTRDCalibTask::AliTRDCalibTask(const char *name)
84 : AliAnalysisTaskSE(name), fESD(0),
92 fTRDCalibraFillHisto(0),
96 fNbTRDTrackOffline(0),
97 fNbTRDTrackStandalone(0),
101 fNbTimeBinOffline(0),
102 fNbTimeBinStandalone(0),
104 fNbClustersOffline(0),
105 fNbClustersStandalone(0),
107 fNbTrackletsOffline(0),
108 fNbTrackletsStandalone(0),
116 fVdriftLinear(kTRUE),
118 fSelectedTrigger(new TObjArray()),
121 fRequirePrimaryVertex(kFALSE),
124 fMinNbContributors(0),
125 fRangePrimaryVertexZ(9999999.0),
131 fNormalizeNbOfCluster(kFALSE),
135 fOfflineTracks(kFALSE),
136 fStandaloneTracks(kFALSE),
138 fVersionGainUsed(-1),
139 fSubVersionGainUsed(-1),
140 fFirstRunGainLocal(-1),
141 fVersionGainLocalUsed(-1),
142 fSubVersionGainLocalUsed(-1),
144 fVersionVdriftUsed(-1),
145 fSubVersionVdriftUsed(-1),
152 // Default constructor
163 // Define input and output slots here
164 // Input slot #0 works with a TChain
165 DefineInput(0, TChain::Class());
167 // Output slot #0 writes into a TList container
168 DefineOutput(1, TList::Class());
172 //____________________________________________________________________________________
173 AliTRDCalibTask::~AliTRDCalibTask()
176 // AliTRDCalibTask destructor
180 if(fNEvents) delete fNEvents;
181 if(fNEventsInput) delete fNEventsInput;
182 if(fNbTRDTrack) delete fNbTRDTrack;
183 if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;
184 if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;
185 if(fNbTPCTRDtrack) delete fNbTPCTRDtrack;
186 if(fNbGoodTracks) delete fNbGoodTracks;
187 if(fNbTimeBin) delete fNbTimeBin;
188 if(fNbTimeBinOffline) delete fNbTimeBinOffline;
189 if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;
190 if(fNbClusters) delete fNbClusters;
191 if(fNbClustersOffline) delete fNbClustersOffline;
192 if(fNbClustersStandalone) delete fNbClustersStandalone;
193 if(fNbTracklets) delete fNbTracklets;
194 if(fNbTrackletsOffline) delete fNbTrackletsOffline;
195 if(fNbTrackletsStandalone) delete fNbTrackletsStandalone;
196 if(fAbsoluteGain) delete fAbsoluteGain;
197 if(fCH2dSum) delete fCH2dSum;
198 if(fPH2dSum) delete fPH2dSum;
199 if(fCH2dSM) delete fCH2dSM;
200 if(fPH2dSM) delete fPH2dSM;
201 if(fCalDetGain) delete fCalDetGain;
203 if(fSelectedTrigger) {
204 fSelectedTrigger->Delete();
205 delete fSelectedTrigger;
208 delete fEsdTrackCuts;
213 //________________________________________________________________________
214 void AliTRDCalibTask::UserCreateOutputObjects()
219 //cout << "AliTRDCalibTask::CreateOutputObjects() IN" << endl;
221 // Number of time bins
223 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
224 fNbTimeBins = cal->GetNumberOfTimeBinsDCS();
225 if(fNbTimeBins <= 0){
226 AliWarning(Form("No of TimeBins from DB [%d] use default [30]", fNbTimeBins));
231 // instance calibration
232 fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
233 fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
234 fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
235 fTRDCalibraFillHisto->SetCH2dOn(); // choose to calibrate the gain
236 fTRDCalibraFillHisto->SetPH2dOn(); // choose to calibrate the drift velocity
237 fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
238 fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
239 fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
240 for(Int_t k = 0; k < 3; k++){
241 if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
242 fTRDCalibraFillHisto->SetNz(k,fNz[k]); // Mode calibration
243 fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]); // Mode calibration
246 if((fNz[k] == 100) && (fNrphi[k] == 100)) {
247 if(fVector2d) AliInfo("The mode all together is not supported by the vector method");
248 fTRDCalibraFillHisto->SetAllTogether(k);
250 if((fNz[k] == 10) && (fNrphi[k] == 10)) {
251 if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");
252 fTRDCalibraFillHisto->SetPerSuperModule(k);
256 // Variables for how to fill
257 fTRDCalibraFillHisto->SetFillWithZero(fFillZero);
258 fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fNormalizeNbOfCluster);
259 fTRDCalibraFillHisto->SetMaxCluster(fMaxCluster);
260 fTRDCalibraFillHisto->SetNbMaxCluster(fNbMaxCluster);
262 // Init with 30 timebins
263 fTRDCalibraFillHisto->Init2Dhistos(fNbTimeBins); // initialise the histos
264 fTRDCalibraFillHisto->SetNumberClusters(fLow); // At least 11 clusters
265 fTRDCalibraFillHisto->SetNumberClustersf(fHigh); // At least 11 clusters
266 fRelativeScale = fTRDCalibraFillHisto->GetRelativeScale(); // Get the relative scale for the gain
269 if(fDebug > 2) fTRDCalibraFillHisto->SetDebugLevel(1); //debug stuff
272 fListHist = new TList();
273 fListHist->SetOwner();
275 fListHist->Add(fTRDCalibraFillHisto->GetCH2d());
276 fListHist->Add(fTRDCalibraFillHisto->GetPH2d());
277 fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());
279 if(fVdriftLinear) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetVdriftLinearFit());
280 if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
281 fNEvents = new TH1I("NEvents","NEvents", 2, 0, 2);
282 fListHist->Add(fNEvents);
283 fNEventsInput = new TH1I("NEventsInput","NEventsInput", 2, 0, 2);
284 fListHist->Add(fNEventsInput);
286 // absolute gain calibration even without AliESDfriend
288 Double_t minPt = 0.001;
289 Double_t maxPt = 10.0;
291 Double_t *binLimLogPt = new Double_t[nBinsPt+1];
292 Double_t *binLimPt = new Double_t[nBinsPt+1];
293 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 ;
294 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
296 fAbsoluteGain = new TH2F("AbsoluteGain","AbsoluteGain", 200, 0.0, 700.0, nBinsPt, binLimPt);
297 fAbsoluteGain->SetYTitle("Momentum at TRD");
298 fAbsoluteGain->SetXTitle("charge deposit [a.u]");
299 fAbsoluteGain->SetZTitle("counts");
300 fAbsoluteGain->SetStats(0);
301 fAbsoluteGain->Sumw2();
302 fListHist->Add(fAbsoluteGain);
306 /////////////////////////////////////////
308 ///////////////////////////////////////
311 // Standart with AliESDfriend
312 fPH2dSM = new TProfile2D("PH2dSM","Nz10Nrphi10"
313 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
315 fPH2dSM->SetYTitle("Det/pad groups");
316 fPH2dSM->SetXTitle("time [#mus]");
317 fPH2dSM->SetZTitle("<PH> [a.u.]");
318 fPH2dSM->SetStats(0);
320 fCH2dSM = new TH2I("CH2dSM","Nz10Nrphi10",50,0,300,18,0,18);
321 fCH2dSM->SetYTitle("Det/pad groups");
322 fCH2dSM->SetXTitle("charge deposit [a.u]");
323 fCH2dSM->SetZTitle("counts");
324 fCH2dSM->SetStats(0);
327 fPH2dSum = new TProfile2D("PH2dSum","Nz100Nrphi100"
328 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
330 fPH2dSum->SetYTitle("Det/pad groups");
331 fPH2dSum->SetXTitle("time [#mus]");
332 fPH2dSum->SetZTitle("<PH> [a.u.]");
333 fPH2dSum->SetStats(0);
335 fCH2dSum = new TH2I("CH2dSum","Nz100Nrphi100",50,0,300,1,0,1);
336 fCH2dSum->SetYTitle("Det/pad groups");
337 fCH2dSum->SetXTitle("charge deposit [a.u]");
338 fCH2dSum->SetZTitle("counts");
339 fCH2dSum->SetStats(0);
342 fNbGoodTracks = new TH2F("NbGoodTracks","NbGoodTracks",500,0.0,2500.0,200,0.0,100.0);
343 fNbGoodTracks->SetXTitle("Nb of good tracks");
344 fNbGoodTracks->SetYTitle("Centrality");
345 fNbGoodTracks->SetStats(0);
349 fListHist->Add(fPH2dSM);
350 fListHist->Add(fCH2dSM);
351 fListHist->Add(fPH2dSum);
352 fListHist->Add(fCH2dSum);
353 fListHist->Add(fNbGoodTracks);
356 /////////////////////////////////////////
357 // Second debug level
358 ///////////////////////////////////////
361 fNbTRDTrack = new TH1F("TRDTrack","TRDTrack",50,0,50);
362 fNbTRDTrack->Sumw2();
363 fNbTRDTrackOffline = new TH1F("TRDTrackOffline","TRDTrackOffline",50,0,50);
364 fNbTRDTrackOffline->Sumw2();
365 fNbTRDTrackStandalone = new TH1F("TRDTrackStandalone","TRDTrackStandalone",50,0,50);
366 fNbTRDTrackStandalone->Sumw2();
367 fNbTPCTRDtrack = new TH2F("NbTPCTRDtrack","NbTPCTRDtrack",100,0,100,100,0,100);
368 fNbTPCTRDtrack->Sumw2();
370 fNbTimeBin = new TH1F("NbTimeBin","NbTimeBin",35,0,35);
372 fNbTimeBinOffline = new TH1F("NbTimeBinOffline","NbTimeBinOffline",35,0,35);
373 fNbTimeBinOffline->Sumw2();
374 fNbTimeBinStandalone = new TH1F("NbTimeBinStandalone","NbTimeBinStandalone",35,0,35);
375 fNbTimeBinStandalone->Sumw2();
377 fNbClusters = new TH1F("NbClusters","",35,0,35);
378 fNbClusters->Sumw2();
379 fNbClustersOffline = new TH1F("NbClustersOffline","",35,0,35);
380 fNbClustersOffline->Sumw2();
381 fNbClustersStandalone = new TH1F("NbClustersStandalone","",35,0,35);
382 fNbClustersStandalone->Sumw2();
384 fNbTracklets = new TH1F("NbTracklets","NbTracklets",540,0.,540.);
385 fNbTracklets->Sumw2();
386 fNbTrackletsOffline = new TH1F("NbTrackletsOffline","NbTrackletsOffline",540,0.,540.);
387 fNbTrackletsOffline->Sumw2();
388 fNbTrackletsStandalone = new TH1F("NbTrackletsStandalone","NbTrackletsStandalone",540,0.,540.);
389 fNbTrackletsStandalone->Sumw2();
391 fListHist->Add(fNbTRDTrack);
392 fListHist->Add(fNbTRDTrackOffline);
393 fListHist->Add(fNbTRDTrackStandalone);
394 fListHist->Add(fNbTPCTRDtrack);
396 fListHist->Add(fNbTimeBin);
397 fListHist->Add(fNbTimeBinOffline);
398 fListHist->Add(fNbTimeBinStandalone);
399 fListHist->Add(fNbClusters);
400 fListHist->Add(fNbClustersOffline);
401 fListHist->Add(fNbClustersStandalone);
402 fListHist->Add(fNbTracklets);
403 fListHist->Add(fNbTrackletsOffline);
404 fListHist->Add(fNbTrackletsStandalone);
408 delete [] binLimLogPt;
411 PostData(1,fListHist);
413 //cout << "AliTRDCalibTask::UserCreateOutputObjects() OUT" << endl;
417 //________________________________________________________________________
418 void AliTRDCalibTask::UserExec(Option_t *)
421 // Filling of the histos
423 //cout << "AliTRDCalibTask::Exec() IN" << endl;
425 // Init Versions and subversions used
426 if((fFirstRunGain==-1) || (fVersionGainUsed==-1) || (fSubVersionGainUsed==-1) || (fFirstRunGainLocal==-1) || (fVersionGainLocalUsed==-1) || (fSubVersionGainLocalUsed==-1) || (fFirstRunVdrift==-1) || (fVersionVdriftUsed==-1) || (fSubVersionVdriftUsed==-1)) {
427 if(!SetVersionSubversion()) {
428 PostData(1, fListHist);
433 fTRDCalibraFillHisto->SetFirstRunGain(fFirstRunGain); // Gain Used
434 fTRDCalibraFillHisto->SetVersionGainUsed(fVersionGainUsed); // Gain Used
435 fTRDCalibraFillHisto->SetSubVersionGainUsed(fSubVersionGainUsed); // Gain Used
436 fTRDCalibraFillHisto->SetFirstRunGainLocal(fFirstRunGainLocal); // Gain Used
437 fTRDCalibraFillHisto->SetVersionGainLocalUsed(fVersionGainLocalUsed); // Gain Used
438 fTRDCalibraFillHisto->SetSubVersionGainLocalUsed(fSubVersionGainLocalUsed); // Gain Used
439 fTRDCalibraFillHisto->SetFirstRunVdrift(fFirstRunVdrift); // Vdrift Used
440 fTRDCalibraFillHisto->SetVersionVdriftUsed(fVersionVdriftUsed); // Vdrift Used
441 fTRDCalibraFillHisto->SetSubVersionVdriftUsed(fSubVersionVdriftUsed); // Vdrift Used
442 fTRDCalibraFillHisto->InitCalDet();
445 // AliLog::SetGlobalLogLevel(AliLog::kError);
446 // cout << "AliTRDCalibTask::Exec() 1" << endl;
447 fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
449 AliError("ESD Event missing");
450 PostData(1, fListHist);
454 const char* type = fESD->GetBeamType();
457 //printf("Counter %d\n",fCounter);
460 fNEventsInput->Fill(1);
462 //cout << "maxEvent = " << fMaxEvent << endl;
463 //if(fCounter%100==0) cout << "fCounter = " << fCounter << endl;
464 if((fMaxEvent != 0) && (fMaxEvent < fCounter)) {
465 PostData(1, fListHist);
468 //if(fCounter%100==0) cout << "fCounter1 = " << fCounter << endl;
469 //cout << "event = " << fCounter << endl;
471 //printf("Counter %d\n",fCounter);
478 if (strstr(type,"p-p")) {
480 //printf("Will check the triggers\n");
482 Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
483 //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
486 for(Int_t k = 0; k < numberOfTriggerSelected; k++){
487 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
488 const TString tString=obString->GetString();
489 if(fESD->IsTriggerClassFired((const char*)tString)) {
496 for(Int_t k = 0; k < numberOfTriggerSelected; k++){
497 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
498 const TString tString=obString->GetString();
499 if(fESD->IsTriggerClassFired((const char*)tString)) {
505 PostData(1, fListHist);
511 //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
512 //printf("Trigger passed\n");
514 ///////////////////////////////
515 // Require a primary vertex
516 //////////////////////////////
517 if(fRequirePrimaryVertex) {
518 const AliESDVertex* vtxESD = 0x0;
519 if (fVtxTPC) vtxESD = fESD->GetPrimaryVertexTPC() ;
520 else if (fVtxSPD) vtxESD = fESD->GetPrimaryVertexSPD() ;
521 else vtxESD = fESD->GetPrimaryVertexTracks() ;
523 PostData(1, fListHist);
526 Int_t nCtrb = vtxESD->GetNContributors();
527 if(nCtrb < fMinNbContributors) {
528 PostData(1, fListHist);
531 Double_t zPosition = vtxESD->GetZ();
532 if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) {
533 PostData(1, fListHist);
539 //printf("Primary vertex passed\n");
541 //////////////////////////////////////
542 // Requirement on number of good tracks
543 //////////////////////////////////////
544 Int_t nGoodParticles = 0;
545 Double_t nbTracks = fESD->GetNumberOfTracks();
546 for(Int_t itrack = 0; itrack < nbTracks; itrack++) {
547 if(ParticleGood(itrack)) nGoodParticles++;
551 AliCentrality *esdCentrality = fESD->GetCentrality();
552 Float_t centrality = esdCentrality->GetCentralityPercentile("V0M");
553 //Float_t centralityb = esdCentrality->GetCentralityPercentile("CL1");
554 fNbGoodTracks->Fill(nGoodParticles,centrality);
555 //printf("centrality %f, centralityb %f\n",centrality,centralityb);
558 if (strstr(type,"Pb-Pb")) {
559 //printf("Will check the number of good tracks\n");
560 if((nGoodParticles < fMinNbTracks) || (nGoodParticles > fMaxNbTracks)) {
561 PostData(1, fListHist);
569 Int_t nbTrdTracks = 0;
571 Int_t nbTrdTracksStandalone = 0;
573 Int_t nbTrdTracksOffline = 0;
575 Int_t nbtrackTPC = 0;
579 if (nbTracks <= 0.0) {
582 fNbTRDTrack->Fill(nbTrdTracks);
583 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
584 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
586 PostData(1, fListHist);
591 fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
593 AliError("fESDfriend not available");
594 PostData(1, fListHist);
598 if(fESDfriend->TestSkipBit()) {
599 PostData(1, fListHist);
603 //printf("has friends\n");
606 ////////////////////////////////////
607 // Check the number of TPC tracks
608 ///////////////////////////////////
609 //printf("Nb of tracks %f\n",nbTracks);
610 for(Int_t itrk = 0; itrk < nbTracks; itrk++){
612 fkEsdTrack = fESD->GetTrack(itrk);
613 ULong_t status = fkEsdTrack->GetStatus();
614 if(status&(AliESDtrack::kTPCout)) nbtrackTPC++;
615 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
617 nbTrdTracksStandalone++;
619 if((status&(AliESDtrack::kTRDin))) {
621 nbTrdTracksOffline++;
625 if((nbtrackTPC>0) && (nbTrdTracks > (3.0*nbtrackTPC))) pass = kFALSE;
629 fNbTRDTrack->Fill(nbTrdTracks);
630 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
631 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
632 fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
637 PostData(1, fListHist);
642 /////////////////////////////////////
643 // Loop on AliESDtrack
644 ////////////////////////////////////
645 //printf("Nb of tracks %f\n",nbTracks);
646 for(int itrk=0; itrk < nbTracks; ++itrk){
649 fkEsdTrack = fESD->GetTrack(itrk);
650 if(!fkEsdTrack) continue;
651 ULong_t status = fkEsdTrack->GetStatus();
652 if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
654 // Quality cuts on the AliESDtrack
655 if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
656 //printf("Not a good track\n");
660 // First Absolute gain calibration
661 Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
662 Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID();
663 if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
664 for(Int_t iPlane = 0; iPlane < 6; ++iPlane){
665 //Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
666 //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
667 //Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
668 //printf("momentum %f, slide %f\n",momentum,slide);
669 if(fkEsdTrack->GetTRDslice(iPlane) > 0.0)
670 fAbsoluteGain->Fill(fkEsdTrack->GetTRDslice(iPlane)*8.0/100.0,
671 fkEsdTrack->GetTRDmomentum(iPlane));
677 Bool_t standalonetrack = kFALSE;
678 Bool_t offlinetrack = kFALSE;
679 //ULong_t status = fkEsdTrack->GetStatus();
681 fFriendTrack = fESDfriend->GetTrack(itrk);
683 //printf("No friend track %d\n",itrk);
686 //////////////////////////////////////
687 // Loop on calibration objects
688 //////////////////////////////////////
691 while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
692 //printf("Name %s\n",fCalibObject->IsA()->GetName());
693 if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
694 //printf("Find the calibration object\n");
697 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
698 standalonetrack = kTRUE;
700 if((status&(AliESDtrack::kTRDin))) {
701 offlinetrack = kTRUE;
708 else if(fStandaloneTracks){
709 if(!standalonetrack){
714 fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
716 //cout << "good" << endl;
717 fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);
718 //printf("Fill fTRDCalibraFillHisto\n");
721 //////////////////////////////////
723 ////////////////////////////////
727 //printf("Enter debug\n");
729 Int_t nbtracklets = 0;
732 Bool_t standalonetracklet = kFALSE;
733 const AliTRDseedV1 *tracklet = 0x0;
734 //////////////////////////////////////
736 /////////////////////////////////////
738 Double_t phtb[AliTRDseedV1::kNtb];
739 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
741 Float_t normalisation = 6.67;
744 for(Int_t itr = 0; itr < 6; ++itr){
746 if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
747 if(!tracklet->IsOK()) continue;
749 standalonetracklet = kFALSE;
750 if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
753 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
757 //Int_t crossrow = 0;
759 // Check no shared clusters
760 //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
761 // if((fcl = tracklet->GetClusters(icc))) crossrow = 1;
768 for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
770 if(!(fCl = tracklet->GetClusters(ic))) continue;
772 time = fCl->GetPadTime();
773 ch = tracklet->GetdQdl(ic);
774 qcl = TMath::Abs(fCl->GetQ());
775 detector = fCl->GetDetector();
776 // Add the charge if shared cluster
777 if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
778 if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
779 qcl += TMath::Abs(fCl->GetQ());
780 //printf("Add the cluster charge\n");
783 if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
784 if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;
785 else sum += ch/normalisation;
788 fNbTimeBin->Fill(time);
789 if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
790 else fNbTimeBinOffline->Fill(time);
793 sector = AliTRDgeometry::GetSector(detector);
796 fNbTracklets->Fill(detector);
797 if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
798 else fNbTrackletsOffline->Fill(detector);
800 fNbClusters->Fill(nbclusters);
801 if(tracklet->IsStandAlone()) fNbClustersStandalone->Fill(nbclusters);
802 else fNbClustersOffline->Fill(nbclusters);
806 if((nbclusters > fLow) && (nbclusters < fHigh)){
807 if(fRelativeScale > 0.0) sum = sum/fRelativeScale;
808 fCH2dSM->Fill(sum,sector+0.5);
809 fCH2dSum->Fill(sum,0.5);
810 Bool_t checknoise = kTRUE;
811 if(fMaxCluster > 0) {
812 if(phtb[0] > fMaxCluster) checknoise = kFALSE;
813 if(fNbTimeBins > fNbMaxCluster) {
814 for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
815 if(phtb[k] > fMaxCluster) checknoise = kFALSE;
820 for(int ic=0; ic<fNbTimeBins; ic++){
822 fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
823 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
827 fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
828 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
835 } // loop on tracklets
839 }// while calibration objects
840 if(nTRDtrackV1 > 0) {
842 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
843 ++nbTrdTracksStandalone;
845 if((status&(AliESDtrack::kTRDin))) {
846 ++nbTrdTracksOffline;
849 //delete fFriendTrack;
853 fNbTRDTrack->Fill(nbTrdTracks);
854 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
855 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
856 fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
860 PostData(1, fListHist);
861 //cout << "AliTRDCalibTask::Exec() OUT" << endl;
864 //________________________________________________________________________
865 void AliTRDCalibTask::Terminate(Option_t *)
871 if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
875 //_______________________________________________________
876 Bool_t AliTRDCalibTask::Load(const Char_t *filename)
879 // Generic container loader
882 if(!TFile::Open(filename)){
883 //AliWarning(Form("Couldn't open file %s.", filename));
887 if(!(o = (TList*)gFile->Get(GetName()))){
888 //AliWarning("Missing histogram container.");
891 fListHist = (TList*)o->Clone(GetName());
895 //_______________________________________________________
896 Bool_t AliTRDCalibTask::Load(TList *lister)
899 // Generic container loader
902 fListHist = (TList*)lister->Clone(GetName());
905 //________________________________________________________________________
906 void AliTRDCalibTask::Plot()
909 // Plot the histos stored in the TList
912 if(!fListHist) return;
914 /////////////////////////////////////
915 // Take the debug stuff
916 /////////////////////////////////////
918 TH1I *nEvents = (TH1I *) fListHist->FindObject("NEvents");
920 TH2F *absoluteGain = (TH2F *) fListHist->FindObject("AbsoluteGain");
922 TH1F *trdTrack = (TH1F *) fListHist->FindObject("TRDTrack");
923 TH1F *trdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");
924 TH1F *trdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");
926 TH2F *tpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");
928 TH1F *nbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");
929 TH1F *nbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");
930 TH1F *nbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");
932 TH1F *nbClusters = (TH1F *) fListHist->FindObject("NbClusters");
933 TH1F *nbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");
934 TH1F *nbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");
936 TH1F *nbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");
937 TH1F *nbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");
938 TH1F *nbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");
940 /////////////////////////////////////
941 // Take the calibration objects
942 /////////////////////////////////////
944 TH2I *ch2d = (TH2I *) fListHist->FindObject("CH2d");
945 TProfile2D *ph2d = (TProfile2D *) fListHist->FindObject("PH2d");
947 TH2I *ch2dSum = (TH2I *) fListHist->FindObject("CH2dSum");
948 TProfile2D *ph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");
950 TH2I *ch2dSM = (TH2I *) fListHist->FindObject("CH2dSM");
951 TProfile2D *ph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");
953 AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
955 ////////////////////////////////////////////////
956 // Add the AliTRDCalibraVdriftLinearFit
957 ///////////////////////////////////////////////
960 TH2S *histolinearfitsum = 0x0;
963 for(Int_t det = 0; det < 540; det++) {
964 if(linearfit->GetLinearFitterHisto(det)){
965 if(TMath::Abs(first)<0.0001){
966 histolinearfitsum = linearfit->GetLinearFitterHisto(det);
970 if (histolinearfitsum) {
971 histolinearfitsum->Add(linearfit->GetLinearFitterHisto(det));
978 ///////////////////////////
980 //////////////////////////
982 gStyle->SetPalette(1);
983 gStyle->SetOptStat(1111);
984 gStyle->SetOptFit(1111);
985 gStyle->SetPadBorderMode(0);
986 gStyle->SetCanvasColor(10);
987 gStyle->SetPadLeftMargin(0.13);
988 gStyle->SetPadRightMargin(0.13);
990 /////////////////////////
992 ////////////////////////
996 TCanvas *debugEvents = new TCanvas("cNEvents","cNEvents",10,10,510,510);
998 if(nEvents) nEvents->Draw();
1004 TCanvas *debugAbsoluteGain = new TCanvas("cAbsoluteGain","cAbsoluteGain",10,10,510,510);
1005 debugAbsoluteGain->cd(1);
1006 if(absoluteGain) absoluteGain->Draw();
1010 if(trdTrack || tpctrdTrack) {
1012 TCanvas *debugtrdtpcTrack = new TCanvas("TRDtracktpctrdtrack","TRDtracktpctrdtrack",10,10,510,510);
1013 debugtrdtpcTrack->Divide(2,1);
1014 debugtrdtpcTrack->cd(1);
1015 if(trdTrack) trdTrack->Draw();
1016 if(trdTrackOffline) trdTrackOffline->Draw("same");
1017 if(trdTrackStandalone) trdTrackStandalone->Draw("same");
1018 TLegend *leg = new TLegend(0.4,0.6,0.89,0.89);
1019 if(trdTrack) leg->AddEntry(trdTrack,"All","p");
1020 if(trdTrackOffline) leg->AddEntry(trdTrackOffline,"Offline","p");
1021 if(trdTrackStandalone) leg->AddEntry(trdTrackStandalone,"Standalone","p");
1023 debugtrdtpcTrack->cd(2);
1024 if(tpctrdTrack) tpctrdTrack->Draw();
1025 TLine *line = new TLine(0.0,0.0,100.0,100.0);
1030 if(nbTimeBin || nbTracklets || nbClusters) {
1032 TCanvas *debugTracklets = new TCanvas("TRDtimebintrackletcluster","TRDtimebintrackletcluster",10,10,510,510);
1033 debugTracklets->Divide(3,1);
1034 debugTracklets->cd(1);
1035 if(nbTimeBin) nbTimeBin->Draw();
1036 if(nbTimeBinOffline) nbTimeBinOffline->Draw("same");
1037 if(nbTimeBinStandalone) nbTimeBinStandalone->Draw("same");
1038 TLegend *lega = new TLegend(0.4,0.6,0.89,0.89);
1039 if(nbTimeBin) lega->AddEntry(nbTimeBin,"All","p");
1040 if(nbTimeBinOffline) lega->AddEntry(nbTimeBinOffline,"Offline","p");
1041 if(nbTimeBinStandalone) lega->AddEntry(nbTimeBinStandalone,"Standalone","p");
1043 debugTracklets->cd(2);
1044 if(nbTracklets) nbTracklets->Draw();
1045 if(nbTrackletsOffline) nbTrackletsOffline->Draw("same");
1046 if(nbTrackletsStandalone) nbTrackletsStandalone->Draw("same");
1047 TLegend *legb = new TLegend(0.4,0.6,0.89,0.89);
1048 if(nbTracklets) legb->AddEntry(nbTracklets,"All","p");
1049 if(nbTrackletsOffline) legb->AddEntry(nbTrackletsOffline,"Offline","p");
1050 if(nbTrackletsStandalone) legb->AddEntry(nbTrackletsStandalone,"Standalone","p");
1052 debugTracklets->cd(3);
1053 if(nbClusters) nbClusters->Draw();
1054 if(nbClustersOffline) nbClustersOffline->Draw("same");
1055 if(nbClustersStandalone) nbClustersStandalone->Draw("same");
1056 TLegend *legc = new TLegend(0.4,0.6,0.89,0.89);
1057 if(nbClusters) legc->AddEntry(nbClusters,"All","p");
1058 if(nbClustersOffline) legc->AddEntry(nbClustersOffline,"Offline","p");
1059 if(nbClustersStandalone) legc->AddEntry(nbClustersStandalone,"Standalone","p");
1064 if(ch2dSum || ph2dSum || histolinearfitsum) {
1066 TCanvas *debugSum = new TCanvas("SumCalibrationObjects","SumCalibrationObjects",10,10,510,510);
1067 debugSum->Divide(3,1);
1069 if(ch2dSum) ch2dSum->Draw("lego");
1071 if(ph2dSum) ph2dSum->Draw("lego");
1073 if(histolinearfitsum) histolinearfitsum->Draw();
1077 if(ch2dSM || ph2dSM) {
1079 TCanvas *debugSM = new TCanvas("SMCalibrationObjects","SMCalibrationObjects",10,10,510,510);
1080 debugSM->Divide(2,1);
1082 if(ch2dSM) ch2dSM->Draw("lego");
1084 if(ph2dSM) ph2dSM->Draw("lego");
1090 TCanvas *debug = new TCanvas("CalibrationObjects","CalibrationObjects",10,10,510,510);
1093 if(ch2d) ch2d->Draw("lego");
1095 if(ph2d) ph2d->Draw("lego");
1100 //_______________________________________________________________________________________
1101 void AliTRDCalibTask::AddTask(const AliTRDCalibTask * calibTask) {
1107 TList *listcalibTask = calibTask->GetList();
1108 if(!listcalibTask) return;
1110 TH1I *nEvents = (TH1I *) listcalibTask->FindObject("NEvents");
1111 TH1I *nEventsInput = (TH1I *) listcalibTask->FindObject("NEventsInput");
1112 TH2F *absoluteGain = (TH2F *) listcalibTask->FindObject("AbsoluteGain");
1114 TH1F *trdTrack = (TH1F *) listcalibTask->FindObject("TRDTrack");
1115 TH1F *trdTrackOffline = (TH1F *) listcalibTask->FindObject("TRDTrackOffline");
1116 TH1F *trdTrackStandalone = (TH1F *) listcalibTask->FindObject("TRDTrackStandalone");
1118 TH2F *tpctrdTrack = (TH2F *) listcalibTask->FindObject("NbTPCTRDtrack");
1120 TH1F *nbTimeBin = (TH1F *) listcalibTask->FindObject("NbTimeBin");
1121 TH1F *nbTimeBinOffline = (TH1F *) listcalibTask->FindObject("NbTimeBinOffline");
1122 TH1F *nbTimeBinStandalone = (TH1F *) listcalibTask->FindObject("NbTimeBinStandalone");
1124 TH1F *nbClusters = (TH1F *) listcalibTask->FindObject("NbClusters");
1125 TH1F *nbClustersOffline = (TH1F *) listcalibTask->FindObject("NbClustersOffline");
1126 TH1F *nbClustersStandalone = (TH1F *) listcalibTask->FindObject("NbClustersStandalone");
1128 TH1F *nbTracklets = (TH1F *) listcalibTask->FindObject("NbTracklets");
1129 TH1F *nbTrackletsOffline = (TH1F *) listcalibTask->FindObject("NbTrackletsOffline");
1130 TH1F *nbTrackletsStandalone = (TH1F *) listcalibTask->FindObject("NbTrackletsStandalone");
1132 TH2I *ch2d = (TH2I *) listcalibTask->FindObject("CH2d");
1133 TProfile2D *ph2d = (TProfile2D *) listcalibTask->FindObject("PH2d");
1134 TProfile2D *prf2d = (TProfile2D *) listcalibTask->FindObject("PRF2d");
1136 TH2I *ch2dSum = (TH2I *) listcalibTask->FindObject("CH2dSum");
1137 TProfile2D *ph2dSum = (TProfile2D *) listcalibTask->FindObject("PH2dSum");
1139 TH2I *ch2dSM = (TH2I *) listcalibTask->FindObject("CH2dSM");
1140 TProfile2D *ph2dSM = (TProfile2D *) listcalibTask->FindObject("PH2dSM");
1142 AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) listcalibTask->FindObject("AliTRDCalibraVdriftLinearFit");
1143 AliTRDCalibraVector *calibraVector = (AliTRDCalibraVector *) listcalibTask->FindObject("AliTRDCalibraVector");
1147 TH1I *inEventsInput = (TH1I *) fListHist->FindObject("NEventsInput");
1148 TH1I *inEvents = (TH1I *) fListHist->FindObject("NEvents");
1149 TH2F *iabsoluteGain = (TH2F *) fListHist->FindObject("AbsoluteGain");
1151 TH1F *itrdTrack = (TH1F *) fListHist->FindObject("TRDTrack");
1152 TH1F *itrdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");
1153 TH1F *itrdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");
1155 TH2F *itpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");
1157 TH1F *inbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");
1158 TH1F *inbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");
1159 TH1F *inbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");
1161 TH1F *inbClusters = (TH1F *) fListHist->FindObject("NbClusters");
1162 TH1F *inbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");
1163 TH1F *inbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");
1165 TH1F *inbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");
1166 TH1F *inbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");
1167 TH1F *inbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");
1169 TH2I *ich2d = (TH2I *) fListHist->FindObject("CH2d");
1170 TProfile2D *iph2d = (TProfile2D *) fListHist->FindObject("PH2d");
1171 TProfile2D *iprf2d = (TProfile2D *) fListHist->FindObject("PRF2d");
1173 TH2I *ich2dSum = (TH2I *) fListHist->FindObject("CH2dSum");
1174 TProfile2D *iph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");
1176 TH2I *ich2dSM = (TH2I *) fListHist->FindObject("CH2dSM");
1177 TProfile2D *iph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");
1179 AliTRDCalibraVdriftLinearFit *ilinearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
1180 AliTRDCalibraVector *icalibraVector = (AliTRDCalibraVector *) fListHist->FindObject("AliTRDCalibraVector");
1187 inEventsInput->Add(nEventsInput);
1188 //printf("Add Events\n");
1191 //printf("Create new Events\n");
1192 inEventsInput = new TH1I(*nEventsInput);
1193 fListHist->Add(inEventsInput);
1199 inEvents->Add(nEvents);
1200 //printf("Add Events\n");
1203 //printf("Create new Events\n");
1204 inEvents = new TH1I(*nEvents);
1205 fListHist->Add(inEvents);
1210 if(iabsoluteGain) iabsoluteGain->Add(absoluteGain);
1212 iabsoluteGain = new TH2F(*absoluteGain);
1213 fListHist->Add(iabsoluteGain);
1218 if(itrdTrack) itrdTrack->Add(trdTrack);
1220 itrdTrack = new TH1F(*trdTrack);
1221 fListHist->Add(itrdTrack);
1225 if(trdTrackOffline) {
1226 if(itrdTrackOffline) itrdTrackOffline->Add(trdTrackOffline);
1228 itrdTrackOffline = new TH1F(*trdTrackOffline);
1229 fListHist->Add(itrdTrackOffline);
1233 if(trdTrackStandalone) {
1234 if(itrdTrackStandalone) itrdTrackStandalone->Add(trdTrackStandalone);
1236 itrdTrackStandalone = new TH1F(*trdTrackStandalone);
1237 fListHist->Add(itrdTrackStandalone);
1242 if(itpctrdTrack) itpctrdTrack->Add(tpctrdTrack);
1244 itpctrdTrack = new TH2F(*tpctrdTrack);
1245 fListHist->Add(itpctrdTrack);
1250 if(inbTimeBin) inbTimeBin->Add(nbTimeBin);
1252 inbTimeBin = new TH1F(*inbTimeBin);
1253 fListHist->Add(inbTimeBin);
1257 if(nbTimeBinOffline) {
1258 if(inbTimeBinOffline) inbTimeBinOffline->Add(nbTimeBinOffline);
1260 inbTimeBinOffline = new TH1F(*nbTimeBinOffline);
1261 fListHist->Add(inbTimeBinOffline);
1265 if(nbTimeBinStandalone) {
1266 if(inbTimeBinStandalone) inbTimeBinStandalone->Add(nbTimeBinStandalone);
1268 inbTimeBinStandalone = new TH1F(*nbTimeBinStandalone);
1269 fListHist->Add(inbTimeBinStandalone);
1274 if(inbClusters) inbClusters->Add(nbClusters);
1276 inbClusters = new TH1F(*nbClusters);
1277 fListHist->Add(inbClusters);
1281 if(nbClustersOffline) {
1282 if(inbClustersOffline) inbClustersOffline->Add(nbClustersOffline);
1284 inbClustersOffline = new TH1F(*nbClustersOffline);
1285 fListHist->Add(inbClustersOffline);
1289 if(nbClustersStandalone) {
1290 if(inbClustersStandalone) inbClustersStandalone->Add(nbClustersStandalone);
1292 inbClustersStandalone = new TH1F(*nbClustersStandalone);
1293 fListHist->Add(inbClustersStandalone);
1298 if(inbTracklets) inbTracklets->Add(nbTracklets);
1300 inbTracklets = new TH1F(*nbTracklets);
1301 fListHist->Add(inbTracklets);
1305 if(nbTrackletsOffline) {
1306 if(inbTrackletsOffline) inbTrackletsOffline->Add(nbTrackletsOffline);
1308 inbTrackletsOffline = new TH1F(*nbTrackletsOffline);
1309 fListHist->Add(inbTrackletsOffline);
1313 if(nbTrackletsStandalone) {
1314 if(inbTrackletsStandalone) inbTrackletsStandalone->Add(nbTrackletsStandalone);
1316 inbTrackletsStandalone = new TH1F(*nbTrackletsStandalone);
1317 fListHist->Add(inbTrackletsStandalone);
1322 if(ich2d) ich2d->Add(ch2d);
1324 ich2d = new TH2I(*ch2d);
1325 fListHist->Add(ich2d);
1330 if(iph2d) iph2d->Add(ph2d);
1332 iph2d = new TProfile2D(*ph2d);
1333 fListHist->Add(iph2d);
1338 if(iprf2d) iprf2d->Add(prf2d);
1340 iprf2d = new TProfile2D(*prf2d);
1341 fListHist->Add(iprf2d);
1346 if(ich2dSum) ich2dSum->Add(ch2dSum);
1348 ich2dSum = new TH2I(*ch2dSum);
1349 fListHist->Add(ich2dSum);
1354 if(iph2dSum) iph2dSum->Add(ph2dSum);
1356 iph2dSum = new TProfile2D(*ph2dSum);
1357 fListHist->Add(iph2dSum);
1362 if(ich2dSM) ich2dSM->Add(ch2dSM);
1364 ich2dSM = new TH2I(*ch2dSM);
1365 fListHist->Add(ich2dSM);
1370 if(iph2dSM) iph2dSM->Add(ph2dSM);
1372 iph2dSM = new TProfile2D(*ph2dSM);
1373 fListHist->Add(iph2dSM);
1378 if(ilinearfit) ilinearfit->Add(linearfit);
1380 ilinearfit = new AliTRDCalibraVdriftLinearFit(*linearfit);
1381 fListHist->Add(ilinearfit);
1386 if(icalibraVector) icalibraVector->Add(calibraVector);
1388 icalibraVector = new AliTRDCalibraVector(*calibraVector);
1389 fListHist->Add(icalibraVector);
1394 //________________________________________________________________________________
1395 Long64_t AliTRDCalibTask::Merge(TCollection *li) {
1401 TIterator* iter = li->MakeIterator();
1402 AliTRDCalibTask* cal = 0;
1404 while ((cal = (AliTRDCalibTask*)iter->Next())) {
1405 if (!cal->InheritsFrom(AliTRDCalibTask::Class())) {
1406 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
1410 // add histograms here...
1418 //_____________________________________________________
1419 Bool_t AliTRDCalibTask::SetVersionSubversion(){
1421 // Load Chamber Gain factors into the Tender supply
1424 printf("SetVersionSubversion\n");
1426 //find previous entry from the UserInfo
1427 TTree *tree=((TChain*)GetInputData(0))->GetTree();
1429 AliError("Tree not found in ESDhandler");
1433 TList *userInfo=(TList*)tree->GetUserInfo();
1435 AliError("No UserInfo found in tree");
1439 TList *cdbList=(TList*)userInfo->FindObject("cdbList");
1441 AliError("No cdbList found in UserInfo");
1442 if (AliLog::GetGlobalLogLevel()>=AliLog::kError) userInfo->Print();
1446 TIter nextCDB(cdbList);
1448 while ( (os=(TObjString*)nextCDB()) ){
1449 if(os->GetString().Contains("TRD/Calib/ChamberGainFactor")){
1450 // Get Old gain calibration
1451 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1452 fFirstRunGain = id->GetFirstRun();
1453 fVersionGainUsed = id->GetVersion();
1454 fSubVersionGainUsed = id->GetSubVersion();
1455 } else if(os->GetString().Contains("TRD/Calib/ChamberVdrift")){
1456 // Get Old drift velocity calibration
1457 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1458 fFirstRunVdrift = id->GetFirstRun();
1459 fVersionVdriftUsed = id->GetVersion();
1460 fSubVersionVdriftUsed = id->GetSubVersion();
1461 } else if(os->GetString().Contains("TRD/Calib/LocalGainFactor")){
1462 // Get Old drift velocity calibration
1463 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1464 fFirstRunGainLocal = id->GetFirstRun();
1465 fVersionGainLocalUsed = id->GetVersion();
1466 fSubVersionGainLocalUsed = id->GetSubVersion();
1470 //printf("VersionGain %d, SubversionGain %d, VersionLocalGain %d, Subversionlocalgain %d, Versionvdrift %d, Subversionvdrift %d\n",fVersionGainUsed,fSubVersionGainUsed,fVersionGainLocalUsed,fSubVersionGainLocalUsed,fVersionVdriftUsed,fSubVersionVdriftUsed);
1473 if((fFirstRunGain < 0) ||
1474 (fFirstRunGainLocal < 0) ||
1475 (fFirstRunVdrift < 0) ||
1476 (fVersionGainUsed < 0) ||
1477 (fVersionGainLocalUsed < 0) ||
1478 (fSubVersionGainUsed < 0) ||
1479 (fSubVersionGainLocalUsed < 0) ||
1480 (fVersionVdriftUsed < 0) ||
1481 (fSubVersionVdriftUsed < 0)) {
1482 AliError("No recent calibration found");
1488 //_________________________________________________________________________________________________________________________
1489 Bool_t AliTRDCalibTask::ParticleGood(int i) const {
1492 // Definition of good tracks
1496 AliESDtrack *track = fESD->GetTrack(i);
1497 if (!track->IsOn(AliESDtrack::kTPCrefit)) return 0; // TPC refit
1498 if (track->GetTPCNcls() < 90) return 0; // number of TPC clusters
1499 if (fabs(track->Eta())>0.8) return 0; // fiducial pseudorapidity
1501 track->GetImpactParametersTPC(r,z);
1502 if (fabs(z)>2.0) return 0; // impact parameter in z
1503 if (fabs(r)>2.0) return 0; // impact parameter in xy