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"
42 #include "TObjArray.h"
48 #include "TIterator.h"
49 #include "TLinearFitter.h"
52 #include "AliAnalysisTask.h"
53 #include "AliAnalysisManager.h"
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"
71 #include "AliTRDCalibraVector.h"
72 #include "AliTRDCalibraFillHisto.h"
73 #include "AliTRDCalibraVdriftLinearFit.h"
75 #include "AliTRDcalibDB.h"
80 #include "AliTRDCalibTask.h"
83 ClassImp(AliTRDCalibTask)
85 //________________________________________________________________________
86 AliTRDCalibTask::AliTRDCalibTask(const char *name)
87 : AliAnalysisTaskSE(name), fESD(0),
95 fTRDCalibraFillHisto(0),
99 fNbTRDTrackOffline(0),
100 fNbTRDTrackStandalone(0),
104 fNbTimeBinOffline(0),
105 fNbTimeBinStandalone(0),
107 fNbClustersOffline(0),
108 fNbClustersStandalone(0),
110 fNbTrackletsOffline(0),
111 fNbTrackletsStandalone(0),
119 fLinearVdriftTest(0),
123 fVdriftLinear(kTRUE),
125 fSelectedTrigger(new TObjArray()),
128 fRequirePrimaryVertex(kFALSE),
131 fMinNbContributors(0),
132 fRangePrimaryVertexZ(9999999.0),
138 fNormalizeNbOfCluster(kFALSE),
142 fOfflineTracks(kFALSE),
143 fStandaloneTracks(kFALSE),
145 fVersionGainUsed(-1),
146 fSubVersionGainUsed(-1),
147 fFirstRunGainLocal(-1),
148 fVersionGainLocalUsed(-1),
149 fSubVersionGainLocalUsed(-1),
151 fVersionVdriftUsed(-1),
152 fSubVersionVdriftUsed(-1),
159 // Default constructor
170 // Define input and output slots here
171 // Input slot #0 works with a TChain
172 DefineInput(0, TChain::Class());
174 // Output slot #0 writes into a TList container
175 DefineOutput(1, TList::Class());
179 //____________________________________________________________________________________
180 AliTRDCalibTask::~AliTRDCalibTask()
183 // AliTRDCalibTask destructor
187 if(fNEvents) delete fNEvents;
188 if(fNEventsInput) delete fNEventsInput;
189 if(fNbTRDTrack) delete fNbTRDTrack;
190 if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;
191 if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;
192 if(fNbTPCTRDtrack) delete fNbTPCTRDtrack;
193 if(fNbGoodTracks) delete fNbGoodTracks;
194 if(fNbTimeBin) delete fNbTimeBin;
195 if(fNbTimeBinOffline) delete fNbTimeBinOffline;
196 if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;
197 if(fNbClusters) delete fNbClusters;
198 if(fNbClustersOffline) delete fNbClustersOffline;
199 if(fNbClustersStandalone) delete fNbClustersStandalone;
200 if(fNbTracklets) delete fNbTracklets;
201 if(fNbTrackletsOffline) delete fNbTrackletsOffline;
202 if(fNbTrackletsStandalone) delete fNbTrackletsStandalone;
203 if(fAbsoluteGain) delete fAbsoluteGain;
204 if(fCH2dSum) delete fCH2dSum;
205 if(fPH2dSum) delete fPH2dSum;
206 if(fCH2dSM) delete fCH2dSM;
207 if(fPH2dSM) delete fPH2dSM;
208 if(fCH2dTest) delete fCH2dTest;
209 if(fPH2dTest) delete fPH2dTest;
210 if(fLinearVdriftTest) delete fLinearVdriftTest;
211 if(fCalDetGain) delete fCalDetGain;
213 if(fSelectedTrigger) {
214 fSelectedTrigger->Delete();
215 delete fSelectedTrigger;
218 delete fEsdTrackCuts;
223 //________________________________________________________________________
224 void AliTRDCalibTask::UserCreateOutputObjects()
229 //cout << "AliTRDCalibTask::CreateOutputObjects() IN" << endl;
231 // Number of time bins
233 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
234 fNbTimeBins = cal->GetNumberOfTimeBinsDCS();
235 if(fNbTimeBins <= 0){
236 AliWarning(Form("No of TimeBins from DB [%d] use default [30]", fNbTimeBins));
242 fListHist = new TList();
243 fListHist->SetOwner();
245 // instance calibration
246 fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
248 fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
249 fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
250 fTRDCalibraFillHisto->SetCH2dOn(); // choose to calibrate the gain
251 fTRDCalibraFillHisto->SetPH2dOn(); // choose to calibrate the drift velocity
252 fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
253 fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
254 fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
255 for(Int_t k = 0; k < 3; k++){
256 if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
257 fTRDCalibraFillHisto->SetNz(k,fNz[k]); // Mode calibration
258 fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]); // Mode calibration
261 if((fNz[k] == 100) && (fNrphi[k] == 100)) {
262 if(fVector2d) AliInfo("The mode all together is not supported by the vector method");
263 fTRDCalibraFillHisto->SetAllTogether(k);
265 if((fNz[k] == 10) && (fNrphi[k] == 10)) {
266 if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");
267 fTRDCalibraFillHisto->SetPerSuperModule(k);
271 // Variables for how to fill
272 fTRDCalibraFillHisto->SetFillWithZero(fFillZero);
273 fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fNormalizeNbOfCluster);
274 fTRDCalibraFillHisto->SetMaxCluster(fMaxCluster);
275 fTRDCalibraFillHisto->SetNbMaxCluster(fNbMaxCluster);
277 // Init with 30 timebins
278 fTRDCalibraFillHisto->Init2Dhistos(fNbTimeBins); // initialise the histos
279 fTRDCalibraFillHisto->SetNumberClusters(fLow); // At least 11 clusters
280 fTRDCalibraFillHisto->SetNumberClustersf(fHigh); // At least 11 clusters
283 if(fDebug > 2) fTRDCalibraFillHisto->SetDebugLevel(1); //debug stuff
286 fListHist->Add(fTRDCalibraFillHisto->GetCH2d());
287 fListHist->Add(fTRDCalibraFillHisto->GetPH2d());
288 fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());
290 if(fVdriftLinear) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetVdriftLinearFit());
291 if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
293 fRelativeScale = fTRDCalibraFillHisto->GetRelativeScale(); // Get the relative scale for the gain
295 fNEvents = new TH1I(Form("NEvents_%s",(const char*)fName),"NEvents", 2, 0, 2);
296 fListHist->Add(fNEvents);
297 fNEventsInput = new TH1I(Form("NEventsInput_%s",(const char*)fName),"NEventsInput", 2, 0, 2);
298 fListHist->Add(fNEventsInput);
300 // absolute gain calibration even without AliESDfriend
302 Double_t minPt = 0.001;
303 Double_t maxPt = 10.0;
305 Double_t *binLimLogPt = new Double_t[nBinsPt+1];
306 Double_t *binLimPt = new Double_t[nBinsPt+1];
307 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 ;
308 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
310 fAbsoluteGain = new TH2F(Form("AbsoluteGain_%s",(const char*)fName),"AbsoluteGain", 200, 0.0, 700.0, nBinsPt, binLimPt);
311 fAbsoluteGain->SetYTitle("Momentum at TRD");
312 fAbsoluteGain->SetXTitle("charge deposit [a.u]");
313 fAbsoluteGain->SetZTitle("counts");
314 fAbsoluteGain->SetStats(0);
315 fAbsoluteGain->Sumw2();
316 fListHist->Add(fAbsoluteGain);
318 /////////////////////////////////////////
320 ///////////////////////////////////////
323 fLinearVdriftTest = new TH2S(Form("LFDV0testversion_%s",(const char*)fName),"LFDV0testversion",36,-0.9,0.9,48,-1.2,1.2);
324 fLinearVdriftTest->SetXTitle("tan(phi_{track})");
325 fLinearVdriftTest->SetYTitle("dy/dt");
326 fLinearVdriftTest->SetZTitle("Number of tracklets");
327 fLinearVdriftTest->SetStats(0);
328 fLinearVdriftTest->SetDirectory(0);
330 // Standart with AliESDfriend
331 fPH2dTest = new TProfile2D(Form("PH2dTest_%s",(const char*)fName),"Nz0Nrphi0"
332 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
334 fPH2dTest->SetYTitle("Det/pad groups");
335 fPH2dTest->SetXTitle("time [#mus]");
336 fPH2dTest->SetZTitle("<PH> [a.u.]");
337 fPH2dTest->SetStats(0);
339 fCH2dTest = new TH2I(Form("CH2dTest_%s",(const char*)fName),"Nz0Nrphi0",50,0,300,540,0,540);
340 fCH2dTest->SetYTitle("Det/pad groups");
341 fCH2dTest->SetXTitle("charge deposit [a.u]");
342 fCH2dTest->SetZTitle("counts");
343 fCH2dTest->SetStats(0);
347 fPH2dSM = new TProfile2D(Form("PH2dSM_%s",(const char*)fName),"Nz10Nrphi10"
348 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
350 fPH2dSM->SetYTitle("Det/pad groups");
351 fPH2dSM->SetXTitle("time [#mus]");
352 fPH2dSM->SetZTitle("<PH> [a.u.]");
353 fPH2dSM->SetStats(0);
355 fCH2dSM = new TH2I(Form("CH2dSM_%s",(const char*)fName),"Nz10Nrphi10",50,0,300,18,0,18);
356 fCH2dSM->SetYTitle("Det/pad groups");
357 fCH2dSM->SetXTitle("charge deposit [a.u]");
358 fCH2dSM->SetZTitle("counts");
359 fCH2dSM->SetStats(0);
362 fPH2dSum = new TProfile2D(Form("PH2dSum_%s",(const char*)fName),"Nz100Nrphi100"
363 ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
365 fPH2dSum->SetYTitle("Det/pad groups");
366 fPH2dSum->SetXTitle("time [#mus]");
367 fPH2dSum->SetZTitle("<PH> [a.u.]");
368 fPH2dSum->SetStats(0);
370 fCH2dSum = new TH2I(Form("CH2dSum_%s",(const char*)fName),"Nz100Nrphi100",50,0,300,1,0,1);
371 fCH2dSum->SetYTitle("Det/pad groups");
372 fCH2dSum->SetXTitle("charge deposit [a.u]");
373 fCH2dSum->SetZTitle("counts");
374 fCH2dSum->SetStats(0);
379 fListHist->Add(fLinearVdriftTest);
380 fListHist->Add(fPH2dTest);
381 fListHist->Add(fCH2dTest);
382 fListHist->Add(fPH2dSM);
383 fListHist->Add(fCH2dSM);
384 fListHist->Add(fPH2dSum);
385 fListHist->Add(fCH2dSum);
389 /////////////////////////////////////////
390 // Second debug level
391 ///////////////////////////////////////
394 fNbGoodTracks = new TH2F(Form("NbGoodTracks_%s",(const char*)fName),"NbGoodTracks",500,0.0,2500.0,200,0.0,100.0);
395 fNbGoodTracks->SetXTitle("Nb of good tracks");
396 fNbGoodTracks->SetYTitle("Centrality");
397 fNbGoodTracks->SetStats(0);
399 fNbTRDTrack = new TH1F(Form("TRDTrack_%s",(const char*)fName),"TRDTrack",50,0,50);
400 fNbTRDTrack->Sumw2();
401 fNbTRDTrackOffline = new TH1F(Form("TRDTrackOffline_%s",(const char*)fName),"TRDTrackOffline",50,0,50);
402 fNbTRDTrackOffline->Sumw2();
403 fNbTRDTrackStandalone = new TH1F(Form("TRDTrackStandalone_%s",(const char*)fName),"TRDTrackStandalone",50,0,50);
404 fNbTRDTrackStandalone->Sumw2();
405 fNbTPCTRDtrack = new TH2F(Form("NbTPCTRDtrack_%s",(const char*)fName),"NbTPCTRDtrack",100,0,100,100,0,100);
406 fNbTPCTRDtrack->Sumw2();
408 fNbTimeBin = new TH1F(Form("NbTimeBin_%s",(const char*)fName),"NbTimeBin",35,0,35);
410 fNbTimeBinOffline = new TH1F(Form("NbTimeBinOffline_%s",(const char*)fName),"NbTimeBinOffline",35,0,35);
411 fNbTimeBinOffline->Sumw2();
412 fNbTimeBinStandalone = new TH1F(Form("NbTimeBinStandalone_%s",(const char*)fName),"NbTimeBinStandalone",35,0,35);
413 fNbTimeBinStandalone->Sumw2();
415 fNbClusters = new TH1F(Form("NbClusters_%s",(const char*)fName),"",35,0,35);
416 fNbClusters->Sumw2();
417 fNbClustersOffline = new TH1F(Form("NbClustersOffline_%s",(const char*)fName),"",35,0,35);
418 fNbClustersOffline->Sumw2();
419 fNbClustersStandalone = new TH1F(Form("NbClustersStandalone_%s",(const char*)fName),"",35,0,35);
420 fNbClustersStandalone->Sumw2();
422 fNbTracklets = new TH1F(Form("NbTracklets_%s",(const char*)fName),"NbTracklets",540,0.,540.);
423 fNbTracklets->Sumw2();
424 fNbTrackletsOffline = new TH1F(Form("NbTrackletsOffline_%s",(const char*)fName),"NbTrackletsOffline",540,0.,540.);
425 fNbTrackletsOffline->Sumw2();
426 fNbTrackletsStandalone = new TH1F(Form("NbTrackletsStandalone_%s",(const char*)fName),"NbTrackletsStandalone",540,0.,540.);
427 fNbTrackletsStandalone->Sumw2();
429 fListHist->Add(fNbGoodTracks);
431 fListHist->Add(fNbTRDTrack);
432 fListHist->Add(fNbTRDTrackOffline);
433 fListHist->Add(fNbTRDTrackStandalone);
434 fListHist->Add(fNbTPCTRDtrack);
436 fListHist->Add(fNbTimeBin);
437 fListHist->Add(fNbTimeBinOffline);
438 fListHist->Add(fNbTimeBinStandalone);
439 fListHist->Add(fNbClusters);
440 fListHist->Add(fNbClustersOffline);
441 fListHist->Add(fNbClustersStandalone);
442 fListHist->Add(fNbTracklets);
443 fListHist->Add(fNbTrackletsOffline);
444 fListHist->Add(fNbTrackletsStandalone);
448 delete [] binLimLogPt;
451 PostData(1,fListHist);
453 //cout << "AliTRDCalibTask::UserCreateOutputObjects() OUT" << endl;
457 //________________________________________________________________________
458 void AliTRDCalibTask::UserExec(Option_t *)
461 // Filling of the histos
463 //cout << "AliTRDCalibTask::Exec() IN" << endl;
465 // Init Versions and subversions used
466 if((fFirstRunGain==-1) || (fVersionGainUsed==-1) || (fSubVersionGainUsed==-1) || (fFirstRunGainLocal==-1) || (fVersionGainLocalUsed==-1) || (fSubVersionGainLocalUsed==-1) || (fFirstRunVdrift==-1) || (fVersionVdriftUsed==-1) || (fSubVersionVdriftUsed==-1)) {
467 if(!SetVersionSubversion()) {
468 PostData(1, fListHist);
474 fTRDCalibraFillHisto->SetFirstRunGain(fFirstRunGain); // Gain Used
475 fTRDCalibraFillHisto->SetVersionGainUsed(fVersionGainUsed); // Gain Used
476 fTRDCalibraFillHisto->SetSubVersionGainUsed(fSubVersionGainUsed); // Gain Used
477 fTRDCalibraFillHisto->SetFirstRunGainLocal(fFirstRunGainLocal); // Gain Used
478 fTRDCalibraFillHisto->SetVersionGainLocalUsed(fVersionGainLocalUsed); // Gain Used
479 fTRDCalibraFillHisto->SetSubVersionGainLocalUsed(fSubVersionGainLocalUsed); // Gain Used
480 fTRDCalibraFillHisto->SetFirstRunVdrift(fFirstRunVdrift); // Vdrift Used
481 fTRDCalibraFillHisto->SetVersionVdriftUsed(fVersionVdriftUsed); // Vdrift Used
482 fTRDCalibraFillHisto->SetSubVersionVdriftUsed(fSubVersionVdriftUsed); // Vdrift Used
483 fTRDCalibraFillHisto->InitCalDet();
488 name += fVersionGainUsed;
490 name += fSubVersionGainUsed;
492 name += fFirstRunGain;
494 fCH2dTest->SetTitle(name);
496 TString namee("Ver");
497 namee += fVersionVdriftUsed;
499 namee += fSubVersionVdriftUsed;
501 namee += fFirstRunVdrift;
502 namee += "Nz0Nrphi0";
503 fPH2dTest->SetTitle(namee);
507 // AliLog::SetGlobalLogLevel(AliLog::kError);
508 // cout << "AliTRDCalibTask::Exec() 1" << endl;
509 fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
511 AliError("ESD Event missing");
512 PostData(1, fListHist);
516 const char* type = fESD->GetBeamType();
519 //printf("Counter %d\n",fCounter);
522 fNEventsInput->Fill(1);
524 //cout << "maxEvent = " << fMaxEvent << endl;
525 //if(fCounter%100==0) cout << "fCounter = " << fCounter << endl;
526 if((fMaxEvent != 0) && (fMaxEvent < fCounter)) {
527 PostData(1, fListHist);
530 //if(fCounter%100==0) cout << "fCounter1 = " << fCounter << endl;
531 //cout << "event = " << fCounter << endl;
533 //printf("Counter %d\n",fCounter);
540 if (strstr(type,"p-p")) {
542 //printf("Will check the triggers\n");
544 Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
545 //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
548 for(Int_t k = 0; k < numberOfTriggerSelected; k++){
549 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
550 const TString tString=obString->GetString();
551 if(fESD->IsTriggerClassFired((const char*)tString)) {
558 for(Int_t k = 0; k < numberOfTriggerSelected; k++){
559 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
560 const TString tString=obString->GetString();
561 if(fESD->IsTriggerClassFired((const char*)tString)) {
567 PostData(1, fListHist);
573 //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
574 //printf("Trigger passed\n");
576 ///////////////////////////////
577 // Require a primary vertex
578 //////////////////////////////
579 if(fRequirePrimaryVertex) {
580 const AliESDVertex* vtxESD = 0x0;
581 if (fVtxTPC) vtxESD = fESD->GetPrimaryVertexTPC() ;
582 else if (fVtxSPD) vtxESD = fESD->GetPrimaryVertexSPD() ;
583 else vtxESD = fESD->GetPrimaryVertexTracks() ;
585 PostData(1, fListHist);
588 Int_t nCtrb = vtxESD->GetNContributors();
589 if(nCtrb < fMinNbContributors) {
590 PostData(1, fListHist);
593 Double_t zPosition = vtxESD->GetZ();
594 if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) {
595 PostData(1, fListHist);
601 //printf("Primary vertex passed\n");
603 //////////////////////////////////////
604 // Requirement on number of good tracks
605 //////////////////////////////////////
606 Int_t nGoodParticles = 0;
607 Double_t nbTracks = fESD->GetNumberOfTracks();
608 for(Int_t itrack = 0; itrack < nbTracks; itrack++) {
609 if(ParticleGood(itrack)) nGoodParticles++;
613 AliCentrality *esdCentrality = fESD->GetCentrality();
614 Float_t centrality = esdCentrality->GetCentralityPercentile("V0M");
615 //Float_t centralityb = esdCentrality->GetCentralityPercentile("CL1");
616 fNbGoodTracks->Fill(nGoodParticles,centrality);
617 //printf("centrality %f, centralityb %f\n",centrality,centralityb);
620 if (strstr(type,"Pb-Pb")) {
621 //printf("Will check the number of good tracks\n");
622 if((nGoodParticles < fMinNbTracks) || (nGoodParticles > fMaxNbTracks)) {
623 PostData(1, fListHist);
631 Int_t nbTrdTracks = 0;
633 Int_t nbTrdTracksStandalone = 0;
635 Int_t nbTrdTracksOffline = 0;
637 Int_t nbtrackTPC = 0;
641 if (nbTracks <= 0.0) {
644 fNbTRDTrack->Fill(nbTrdTracks);
645 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
646 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
648 PostData(1, fListHist);
653 fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
655 AliError("fESDfriend not available");
656 PostData(1, fListHist);
660 if(fESDfriend->TestSkipBit()) {
661 PostData(1, fListHist);
665 //printf("has friends\n");
667 /////////////////////////////////////
668 // Loop on AliESDtrack
669 ////////////////////////////////////
670 //printf("Nb of tracks %f\n",nbTracks);
671 for(int itrk=0; itrk < nbTracks; ++itrk){
674 fkEsdTrack = fESD->GetTrack(itrk);
675 if(!fkEsdTrack) continue;
676 ULong_t status = fkEsdTrack->GetStatus();
677 if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
679 // Quality cuts on the AliESDtrack
680 if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
681 //printf("Not a good track\n");
685 // First Absolute gain calibration
686 Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
687 Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID();
688 //printf("Number of trd tracklets %d and PID trd tracklets %d\n",trdNTracklets,trdNTrackletsPID);
689 if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
690 for(Int_t iPlane = 0; iPlane < 6; ++iPlane){
691 //Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
692 //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
693 //Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
694 //printf("momentum %f, slide %f\n",momentum,slide);
695 if(fkEsdTrack->GetTRDslice(iPlane) > 0.0)
696 fAbsoluteGain->Fill(fkEsdTrack->GetTRDslice(iPlane)*8.0/100.0,
697 fkEsdTrack->GetTRDmomentum(iPlane));
703 Bool_t standalonetrack = kFALSE;
704 Bool_t offlinetrack = kFALSE;
705 //ULong_t status = fkEsdTrack->GetStatus();
707 fFriendTrack = fESDfriend->GetTrack(itrk);
709 //printf("No friend track %d\n",itrk);
712 //////////////////////////////////////
713 // Loop on calibration objects
714 //////////////////////////////////////
717 while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
718 //printf("Name %s\n",fCalibObject->IsA()->GetName());
719 if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
720 //printf("Find the calibration object\n");
723 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
724 standalonetrack = kTRUE;
726 if((status&(AliESDtrack::kTRDin))) {
727 offlinetrack = kTRUE;
734 else if(fStandaloneTracks){
735 if(!standalonetrack){
740 fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
741 if(good && fOnInstance) {
742 //cout << "good" << endl;
743 fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);
744 //printf("Fill fTRDCalibraFillHisto\n");
747 //////////////////////////////////
749 ////////////////////////////////
753 //printf("Enter debug\n");
755 Int_t nbtracklets = 0;
758 Bool_t standalonetracklet = kFALSE;
759 const AliTRDseedV1 *tracklet = 0x0;
760 //////////////////////////////////////
762 /////////////////////////////////////
764 Double_t phtb[AliTRDseedV1::kNtb];
765 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
767 Float_t normalisation = 6.67;
770 for(Int_t itr = 0; itr < 6; ++itr){
772 if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
773 if(!tracklet->IsOK()) continue;
775 standalonetracklet = kFALSE;
776 if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
779 memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
783 //Int_t crossrow = 0;
785 // Check no shared clusters
786 //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
787 // if((fcl = tracklet->GetClusters(icc))) crossrow = 1;
794 for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
796 if(!(fCl = tracklet->GetClusters(ic))) continue;
798 time = fCl->GetPadTime();
799 ch = tracklet->GetdQdl(ic);
800 qcl = TMath::Abs(fCl->GetQ());
801 detector = fCl->GetDetector();
802 // Add the charge if shared cluster
803 if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
804 if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
805 qcl += TMath::Abs(fCl->GetQ());
806 //printf("Add the cluster charge\n");
809 if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
810 if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;
811 else sum += ch/normalisation;
814 fNbTimeBin->Fill(time);
815 if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
816 else fNbTimeBinOffline->Fill(time);
819 sector = AliTRDgeometry::GetSector(detector);
822 fNbTracklets->Fill(detector);
823 if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
824 else fNbTrackletsOffline->Fill(detector);
826 fNbClusters->Fill(nbclusters);
827 if(tracklet->IsStandAlone()) fNbClustersStandalone->Fill(nbclusters);
828 else fNbClustersOffline->Fill(nbclusters);
831 if((nbclusters > fLow) && (nbclusters < fHigh)){
832 if(fRelativeScale > 0.0) sum = sum/fRelativeScale;
833 fCH2dTest->Fill(sum,detector+0.5);
834 fCH2dSM->Fill(sum,sector+0.5);
835 fCH2dSum->Fill(sum,0.5);
836 Bool_t checknoise = kTRUE;
837 if(fMaxCluster > 0) {
838 if(phtb[0] > fMaxCluster) checknoise = kFALSE;
839 if(fNbTimeBins > fNbMaxCluster) {
840 for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
841 if(phtb[k] > fMaxCluster) checknoise = kFALSE;
846 for(int ic=0; ic<fNbTimeBins; ic++){
848 fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
849 fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
850 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
854 fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
855 fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
856 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
862 if(detector == 0) FindP1TrackPHtrackletV1Test(tracklet,nbclusters);
864 } // loop on tracklets
869 }// while calibration objects
870 if(nTRDtrackV1 > 0) {
872 if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
873 ++nbTrdTracksStandalone;
875 if((status&(AliESDtrack::kTRDin))) {
876 ++nbTrdTracksOffline;
879 //delete fFriendTrack;
883 fNbTRDTrack->Fill(nbTrdTracks);
884 fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
885 fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
886 fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
890 PostData(1, fListHist);
891 //cout << "AliTRDCalibTask::Exec() OUT" << endl;
894 //________________________________________________________________________
895 void AliTRDCalibTask::Terminate(Option_t *)
901 if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
905 //_______________________________________________________
906 Bool_t AliTRDCalibTask::Load(const Char_t *filename)
909 // Generic container loader
912 if(!TFile::Open(filename)){
913 //AliWarning(Form("Couldn't open file %s.", filename));
917 if(!(o = (TList*)gFile->Get(GetName()))){
918 //AliWarning("Missing histogram container.");
921 fListHist = (TList*)o->Clone(GetName());
925 //_______________________________________________________
926 Bool_t AliTRDCalibTask::Load(TList *lister)
929 // Generic container loader
932 fListHist = (TList*)lister->Clone(GetName());
935 //________________________________________________________________________
936 void AliTRDCalibTask::Plot()
939 // Plot the histos stored in the TList
942 if(!fListHist) return;
944 /////////////////////////////////////
945 // Take the debug stuff
946 /////////////////////////////////////
948 TH1I *nEvents = (TH1I *) fListHist->FindObject("NEvents");
950 TH2F *absoluteGain = (TH2F *) fListHist->FindObject("AbsoluteGain");
952 TH1F *trdTrack = (TH1F *) fListHist->FindObject("TRDTrack");
953 TH1F *trdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");
954 TH1F *trdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");
956 TH2F *tpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");
958 TH1F *nbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");
959 TH1F *nbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");
960 TH1F *nbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");
962 TH1F *nbClusters = (TH1F *) fListHist->FindObject("NbClusters");
963 TH1F *nbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");
964 TH1F *nbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");
966 TH1F *nbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");
967 TH1F *nbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");
968 TH1F *nbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");
970 /////////////////////////////////////
971 // Take the calibration objects
972 /////////////////////////////////////
974 TH2I *ch2d = (TH2I *) fListHist->FindObject("CH2d");
975 TProfile2D *ph2d = (TProfile2D *) fListHist->FindObject("PH2d");
977 TH2I *ch2dSum = (TH2I *) fListHist->FindObject("CH2dSum");
978 TProfile2D *ph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");
980 TH2I *ch2dSM = (TH2I *) fListHist->FindObject("CH2dSM");
981 TProfile2D *ph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");
983 AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
985 ////////////////////////////////////////////////
986 // Add the AliTRDCalibraVdriftLinearFit
987 ///////////////////////////////////////////////
990 TH2S *histolinearfitsum = 0x0;
993 for(Int_t det = 0; det < 540; det++) {
994 if(linearfit->GetLinearFitterHisto(det)){
995 if(TMath::Abs(first)<0.0001){
996 histolinearfitsum = linearfit->GetLinearFitterHisto(det);
1000 if (histolinearfitsum) {
1001 histolinearfitsum->Add(linearfit->GetLinearFitterHisto(det));
1008 ///////////////////////////
1010 //////////////////////////
1012 gStyle->SetPalette(1);
1013 gStyle->SetOptStat(1111);
1014 gStyle->SetOptFit(1111);
1015 gStyle->SetPadBorderMode(0);
1016 gStyle->SetCanvasColor(10);
1017 gStyle->SetPadLeftMargin(0.13);
1018 gStyle->SetPadRightMargin(0.13);
1020 /////////////////////////
1022 ////////////////////////
1026 TCanvas *debugEvents = new TCanvas("cNEvents","cNEvents",10,10,510,510);
1028 if(nEvents) nEvents->Draw();
1034 TCanvas *debugAbsoluteGain = new TCanvas("cAbsoluteGain","cAbsoluteGain",10,10,510,510);
1035 debugAbsoluteGain->cd(1);
1036 if(absoluteGain) absoluteGain->Draw();
1040 if(trdTrack || tpctrdTrack) {
1042 TCanvas *debugtrdtpcTrack = new TCanvas("TRDtracktpctrdtrack","TRDtracktpctrdtrack",10,10,510,510);
1043 debugtrdtpcTrack->Divide(2,1);
1044 debugtrdtpcTrack->cd(1);
1045 if(trdTrack) trdTrack->Draw();
1046 if(trdTrackOffline) trdTrackOffline->Draw("same");
1047 if(trdTrackStandalone) trdTrackStandalone->Draw("same");
1048 TLegend *leg = new TLegend(0.4,0.6,0.89,0.89);
1049 if(trdTrack) leg->AddEntry(trdTrack,"All","p");
1050 if(trdTrackOffline) leg->AddEntry(trdTrackOffline,"Offline","p");
1051 if(trdTrackStandalone) leg->AddEntry(trdTrackStandalone,"Standalone","p");
1053 debugtrdtpcTrack->cd(2);
1054 if(tpctrdTrack) tpctrdTrack->Draw();
1055 TLine *line = new TLine(0.0,0.0,100.0,100.0);
1060 if(nbTimeBin || nbTracklets || nbClusters) {
1062 TCanvas *debugTracklets = new TCanvas("TRDtimebintrackletcluster","TRDtimebintrackletcluster",10,10,510,510);
1063 debugTracklets->Divide(3,1);
1064 debugTracklets->cd(1);
1065 if(nbTimeBin) nbTimeBin->Draw();
1066 if(nbTimeBinOffline) nbTimeBinOffline->Draw("same");
1067 if(nbTimeBinStandalone) nbTimeBinStandalone->Draw("same");
1068 TLegend *lega = new TLegend(0.4,0.6,0.89,0.89);
1069 if(nbTimeBin) lega->AddEntry(nbTimeBin,"All","p");
1070 if(nbTimeBinOffline) lega->AddEntry(nbTimeBinOffline,"Offline","p");
1071 if(nbTimeBinStandalone) lega->AddEntry(nbTimeBinStandalone,"Standalone","p");
1073 debugTracklets->cd(2);
1074 if(nbTracklets) nbTracklets->Draw();
1075 if(nbTrackletsOffline) nbTrackletsOffline->Draw("same");
1076 if(nbTrackletsStandalone) nbTrackletsStandalone->Draw("same");
1077 TLegend *legb = new TLegend(0.4,0.6,0.89,0.89);
1078 if(nbTracklets) legb->AddEntry(nbTracklets,"All","p");
1079 if(nbTrackletsOffline) legb->AddEntry(nbTrackletsOffline,"Offline","p");
1080 if(nbTrackletsStandalone) legb->AddEntry(nbTrackletsStandalone,"Standalone","p");
1082 debugTracklets->cd(3);
1083 if(nbClusters) nbClusters->Draw();
1084 if(nbClustersOffline) nbClustersOffline->Draw("same");
1085 if(nbClustersStandalone) nbClustersStandalone->Draw("same");
1086 TLegend *legc = new TLegend(0.4,0.6,0.89,0.89);
1087 if(nbClusters) legc->AddEntry(nbClusters,"All","p");
1088 if(nbClustersOffline) legc->AddEntry(nbClustersOffline,"Offline","p");
1089 if(nbClustersStandalone) legc->AddEntry(nbClustersStandalone,"Standalone","p");
1094 if(ch2dSum || ph2dSum || histolinearfitsum) {
1096 TCanvas *debugSum = new TCanvas("SumCalibrationObjects","SumCalibrationObjects",10,10,510,510);
1097 debugSum->Divide(3,1);
1099 if(ch2dSum) ch2dSum->Draw("lego");
1101 if(ph2dSum) ph2dSum->Draw("lego");
1103 if(histolinearfitsum) histolinearfitsum->Draw();
1107 if(ch2dSM || ph2dSM) {
1109 TCanvas *debugSM = new TCanvas("SMCalibrationObjects","SMCalibrationObjects",10,10,510,510);
1110 debugSM->Divide(2,1);
1112 if(ch2dSM) ch2dSM->Draw("lego");
1114 if(ph2dSM) ph2dSM->Draw("lego");
1120 TCanvas *debug = new TCanvas("CalibrationObjects","CalibrationObjects",10,10,510,510);
1123 if(ch2d) ch2d->Draw("lego");
1125 if(ph2d) ph2d->Draw("lego");
1130 //_______________________________________________________________________________________
1131 void AliTRDCalibTask::AddTask(const AliTRDCalibTask * calibTask) {
1137 TList *listcalibTask = calibTask->GetList();
1138 if(!listcalibTask) return;
1140 TH1I *nEvents = (TH1I *) listcalibTask->FindObject("NEvents");
1141 TH1I *nEventsInput = (TH1I *) listcalibTask->FindObject("NEventsInput");
1142 TH2F *absoluteGain = (TH2F *) listcalibTask->FindObject("AbsoluteGain");
1144 TH1F *trdTrack = (TH1F *) listcalibTask->FindObject("TRDTrack");
1145 TH1F *trdTrackOffline = (TH1F *) listcalibTask->FindObject("TRDTrackOffline");
1146 TH1F *trdTrackStandalone = (TH1F *) listcalibTask->FindObject("TRDTrackStandalone");
1148 TH2F *tpctrdTrack = (TH2F *) listcalibTask->FindObject("NbTPCTRDtrack");
1150 TH1F *nbTimeBin = (TH1F *) listcalibTask->FindObject("NbTimeBin");
1151 TH1F *nbTimeBinOffline = (TH1F *) listcalibTask->FindObject("NbTimeBinOffline");
1152 TH1F *nbTimeBinStandalone = (TH1F *) listcalibTask->FindObject("NbTimeBinStandalone");
1154 TH1F *nbClusters = (TH1F *) listcalibTask->FindObject("NbClusters");
1155 TH1F *nbClustersOffline = (TH1F *) listcalibTask->FindObject("NbClustersOffline");
1156 TH1F *nbClustersStandalone = (TH1F *) listcalibTask->FindObject("NbClustersStandalone");
1158 TH1F *nbTracklets = (TH1F *) listcalibTask->FindObject("NbTracklets");
1159 TH1F *nbTrackletsOffline = (TH1F *) listcalibTask->FindObject("NbTrackletsOffline");
1160 TH1F *nbTrackletsStandalone = (TH1F *) listcalibTask->FindObject("NbTrackletsStandalone");
1162 TH2I *ch2d = (TH2I *) listcalibTask->FindObject("CH2d");
1163 TProfile2D *ph2d = (TProfile2D *) listcalibTask->FindObject("PH2d");
1164 TProfile2D *prf2d = (TProfile2D *) listcalibTask->FindObject("PRF2d");
1166 TH2I *ch2dSum = (TH2I *) listcalibTask->FindObject("CH2dSum");
1167 TProfile2D *ph2dSum = (TProfile2D *) listcalibTask->FindObject("PH2dSum");
1169 TH2I *ch2dSM = (TH2I *) listcalibTask->FindObject("CH2dSM");
1170 TProfile2D *ph2dSM = (TProfile2D *) listcalibTask->FindObject("PH2dSM");
1172 AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) listcalibTask->FindObject("AliTRDCalibraVdriftLinearFit");
1173 AliTRDCalibraVector *calibraVector = (AliTRDCalibraVector *) listcalibTask->FindObject("AliTRDCalibraVector");
1177 TH1I *inEventsInput = (TH1I *) fListHist->FindObject("NEventsInput");
1178 TH1I *inEvents = (TH1I *) fListHist->FindObject("NEvents");
1179 TH2F *iabsoluteGain = (TH2F *) fListHist->FindObject("AbsoluteGain");
1181 TH1F *itrdTrack = (TH1F *) fListHist->FindObject("TRDTrack");
1182 TH1F *itrdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");
1183 TH1F *itrdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");
1185 TH2F *itpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");
1187 TH1F *inbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");
1188 TH1F *inbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");
1189 TH1F *inbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");
1191 TH1F *inbClusters = (TH1F *) fListHist->FindObject("NbClusters");
1192 TH1F *inbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");
1193 TH1F *inbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");
1195 TH1F *inbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");
1196 TH1F *inbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");
1197 TH1F *inbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");
1199 TH2I *ich2d = (TH2I *) fListHist->FindObject("CH2d");
1200 TProfile2D *iph2d = (TProfile2D *) fListHist->FindObject("PH2d");
1201 TProfile2D *iprf2d = (TProfile2D *) fListHist->FindObject("PRF2d");
1203 TH2I *ich2dSum = (TH2I *) fListHist->FindObject("CH2dSum");
1204 TProfile2D *iph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");
1206 TH2I *ich2dSM = (TH2I *) fListHist->FindObject("CH2dSM");
1207 TProfile2D *iph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");
1209 AliTRDCalibraVdriftLinearFit *ilinearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
1210 AliTRDCalibraVector *icalibraVector = (AliTRDCalibraVector *) fListHist->FindObject("AliTRDCalibraVector");
1217 inEventsInput->Add(nEventsInput);
1218 //printf("Add Events\n");
1221 //printf("Create new Events\n");
1222 inEventsInput = new TH1I(*nEventsInput);
1223 fListHist->Add(inEventsInput);
1229 inEvents->Add(nEvents);
1230 //printf("Add Events\n");
1233 //printf("Create new Events\n");
1234 inEvents = new TH1I(*nEvents);
1235 fListHist->Add(inEvents);
1240 if(iabsoluteGain) iabsoluteGain->Add(absoluteGain);
1242 iabsoluteGain = new TH2F(*absoluteGain);
1243 fListHist->Add(iabsoluteGain);
1248 if(itrdTrack) itrdTrack->Add(trdTrack);
1250 itrdTrack = new TH1F(*trdTrack);
1251 fListHist->Add(itrdTrack);
1255 if(trdTrackOffline) {
1256 if(itrdTrackOffline) itrdTrackOffline->Add(trdTrackOffline);
1258 itrdTrackOffline = new TH1F(*trdTrackOffline);
1259 fListHist->Add(itrdTrackOffline);
1263 if(trdTrackStandalone) {
1264 if(itrdTrackStandalone) itrdTrackStandalone->Add(trdTrackStandalone);
1266 itrdTrackStandalone = new TH1F(*trdTrackStandalone);
1267 fListHist->Add(itrdTrackStandalone);
1272 if(itpctrdTrack) itpctrdTrack->Add(tpctrdTrack);
1274 itpctrdTrack = new TH2F(*tpctrdTrack);
1275 fListHist->Add(itpctrdTrack);
1280 if(inbTimeBin) inbTimeBin->Add(nbTimeBin);
1282 inbTimeBin = new TH1F(*inbTimeBin);
1283 fListHist->Add(inbTimeBin);
1287 if(nbTimeBinOffline) {
1288 if(inbTimeBinOffline) inbTimeBinOffline->Add(nbTimeBinOffline);
1290 inbTimeBinOffline = new TH1F(*nbTimeBinOffline);
1291 fListHist->Add(inbTimeBinOffline);
1295 if(nbTimeBinStandalone) {
1296 if(inbTimeBinStandalone) inbTimeBinStandalone->Add(nbTimeBinStandalone);
1298 inbTimeBinStandalone = new TH1F(*nbTimeBinStandalone);
1299 fListHist->Add(inbTimeBinStandalone);
1304 if(inbClusters) inbClusters->Add(nbClusters);
1306 inbClusters = new TH1F(*nbClusters);
1307 fListHist->Add(inbClusters);
1311 if(nbClustersOffline) {
1312 if(inbClustersOffline) inbClustersOffline->Add(nbClustersOffline);
1314 inbClustersOffline = new TH1F(*nbClustersOffline);
1315 fListHist->Add(inbClustersOffline);
1319 if(nbClustersStandalone) {
1320 if(inbClustersStandalone) inbClustersStandalone->Add(nbClustersStandalone);
1322 inbClustersStandalone = new TH1F(*nbClustersStandalone);
1323 fListHist->Add(inbClustersStandalone);
1328 if(inbTracklets) inbTracklets->Add(nbTracklets);
1330 inbTracklets = new TH1F(*nbTracklets);
1331 fListHist->Add(inbTracklets);
1335 if(nbTrackletsOffline) {
1336 if(inbTrackletsOffline) inbTrackletsOffline->Add(nbTrackletsOffline);
1338 inbTrackletsOffline = new TH1F(*nbTrackletsOffline);
1339 fListHist->Add(inbTrackletsOffline);
1343 if(nbTrackletsStandalone) {
1344 if(inbTrackletsStandalone) inbTrackletsStandalone->Add(nbTrackletsStandalone);
1346 inbTrackletsStandalone = new TH1F(*nbTrackletsStandalone);
1347 fListHist->Add(inbTrackletsStandalone);
1352 if(ich2d) ich2d->Add(ch2d);
1354 ich2d = new TH2I(*ch2d);
1355 fListHist->Add(ich2d);
1360 if(iph2d) iph2d->Add(ph2d);
1362 iph2d = new TProfile2D(*ph2d);
1363 fListHist->Add(iph2d);
1368 if(iprf2d) iprf2d->Add(prf2d);
1370 iprf2d = new TProfile2D(*prf2d);
1371 fListHist->Add(iprf2d);
1376 if(ich2dSum) ich2dSum->Add(ch2dSum);
1378 ich2dSum = new TH2I(*ch2dSum);
1379 fListHist->Add(ich2dSum);
1384 if(iph2dSum) iph2dSum->Add(ph2dSum);
1386 iph2dSum = new TProfile2D(*ph2dSum);
1387 fListHist->Add(iph2dSum);
1392 if(ich2dSM) ich2dSM->Add(ch2dSM);
1394 ich2dSM = new TH2I(*ch2dSM);
1395 fListHist->Add(ich2dSM);
1400 if(iph2dSM) iph2dSM->Add(ph2dSM);
1402 iph2dSM = new TProfile2D(*ph2dSM);
1403 fListHist->Add(iph2dSM);
1408 if(ilinearfit) ilinearfit->Add(linearfit);
1410 ilinearfit = new AliTRDCalibraVdriftLinearFit(*linearfit);
1411 fListHist->Add(ilinearfit);
1416 if(icalibraVector) icalibraVector->Add(calibraVector);
1418 icalibraVector = new AliTRDCalibraVector(*calibraVector);
1419 fListHist->Add(icalibraVector);
1424 //________________________________________________________________________________
1425 Long64_t AliTRDCalibTask::Merge(TCollection *li) {
1431 TIterator* iter = li->MakeIterator();
1432 AliTRDCalibTask* cal = 0;
1434 while ((cal = (AliTRDCalibTask*)iter->Next())) {
1435 if (!cal->InheritsFrom(AliTRDCalibTask::Class())) {
1436 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
1440 // add histograms here...
1448 //_____________________________________________________
1449 Bool_t AliTRDCalibTask::SetVersionSubversion(){
1451 // Load Chamber Gain factors into the Tender supply
1454 printf("SetVersionSubversion\n");
1456 //find previous entry from the UserInfo
1457 TTree *tree=((TChain*)GetInputData(0))->GetTree();
1459 AliError("Tree not found in ESDhandler");
1463 TList *userInfo=(TList*)tree->GetUserInfo();
1465 AliError("No UserInfo found in tree");
1469 TList *cdbList=(TList*)userInfo->FindObject("cdbList");
1471 AliError("No cdbList found in UserInfo");
1472 if (AliLog::GetGlobalLogLevel()>=AliLog::kError) userInfo->Print();
1476 TIter nextCDB(cdbList);
1478 while ( (os=(TObjString*)nextCDB()) ){
1479 if(os->GetString().Contains("TRD/Calib/ChamberGainFactor")){
1480 // Get Old gain calibration
1481 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1482 fFirstRunGain = id->GetFirstRun();
1483 fVersionGainUsed = id->GetVersion();
1484 fSubVersionGainUsed = id->GetSubVersion();
1485 } else if(os->GetString().Contains("TRD/Calib/ChamberVdrift")){
1486 // Get Old drift velocity calibration
1487 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1488 fFirstRunVdrift = id->GetFirstRun();
1489 fVersionVdriftUsed = id->GetVersion();
1490 fSubVersionVdriftUsed = id->GetSubVersion();
1491 } else if(os->GetString().Contains("TRD/Calib/LocalGainFactor")){
1492 // Get Old drift velocity calibration
1493 AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1494 fFirstRunGainLocal = id->GetFirstRun();
1495 fVersionGainLocalUsed = id->GetVersion();
1496 fSubVersionGainLocalUsed = id->GetSubVersion();
1500 //printf("VersionGain %d, SubversionGain %d, VersionLocalGain %d, Subversionlocalgain %d, Versionvdrift %d, Subversionvdrift %d\n",fVersionGainUsed,fSubVersionGainUsed,fVersionGainLocalUsed,fSubVersionGainLocalUsed,fVersionVdriftUsed,fSubVersionVdriftUsed);
1503 if((fFirstRunGain < 0) ||
1504 (fFirstRunGainLocal < 0) ||
1505 (fFirstRunVdrift < 0) ||
1506 (fVersionGainUsed < 0) ||
1507 (fVersionGainLocalUsed < 0) ||
1508 (fSubVersionGainUsed < 0) ||
1509 (fSubVersionGainLocalUsed < 0) ||
1510 (fVersionVdriftUsed < 0) ||
1511 (fSubVersionVdriftUsed < 0)) {
1512 AliError("No recent calibration found");
1518 //_________________________________________________________________________________________________________________________
1519 Bool_t AliTRDCalibTask::ParticleGood(int i) const {
1522 // Definition of good tracks
1526 AliESDtrack *track = fESD->GetTrack(i);
1527 if (!track->IsOn(AliESDtrack::kTPCrefit)) return 0; // TPC refit
1528 if (track->GetTPCNcls() < 90) return 0; // number of TPC clusters
1529 if (fabs(track->Eta())>0.8) return 0; // fiducial pseudorapidity
1531 track->GetImpactParametersTPC(r,z);
1532 if (fabs(z)>2.0) return 0; // impact parameter in z
1533 if (fabs(r)>2.0) return 0; // impact parameter in xy
1539 //______________________________________________________________________________________________________________________
1540 Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1543 // Drift velocity calibration:
1544 // Fit the clusters with a straight line
1545 // From the slope find the drift velocity
1548 ////////////////////////////////////////////////
1549 //Number of points: if less than 3 return kFALSE
1550 /////////////////////////////////////////////////
1551 if(nbclusters <= 2) return kFALSE;
1556 // results of the linear fit
1557 Double_t dydt = 0.0; // dydt tracklet after straight line fit
1558 Double_t errorpar = 0.0; // error after straight line fit on dy/dt
1559 Double_t pointError = 0.0; // error after straight line fit
1560 // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
1561 Int_t crossrow = 0; // if it crosses a pad row
1562 Int_t rowp = -1; // if it crosses a pad row
1563 Float_t tnt = tracklet->GetTilt(); // tan tiltingangle
1564 TLinearFitter linearFitterTracklet(2,"pol1");
1565 linearFitterTracklet.StoreData(kTRUE);
1568 ///////////////////////////////////////////
1569 // Take the parameters of the track
1570 //////////////////////////////////////////
1571 // take now the snp, tnp and tgl from the track
1572 Double_t snp = tracklet->GetSnp(); // sin dy/dx at the end of the chamber
1573 Double_t tnp = 0.0; // dy/dx at the end of the chamber
1574 if( TMath::Abs(snp) < 1.){
1575 tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1577 Double_t tgl = tracklet->GetTgl(); // dz/dl
1578 Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp); // dz/dx calculated from dz/dl
1580 //Double_t tnp = tracklet->GetYref(1); // dy/dx at the entrance of the chamber
1581 //Double_t tgl = tracklet->GetZref(1); // dz/dl at the entrance of the chamber
1582 //Double_t dzdx = tgl; //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1583 // at the end with correction due to linear fit
1584 //Double_t tnp = tracklet->GetYfit(1); // dy/dx at the end of the chamber after fit correction
1585 //Double_t tgl = tracklet->GetZfit(1); // dz/dl at the end of the chamber after fit correction
1588 ////////////////////////////
1589 // loop over the clusters
1590 ////////////////////////////
1592 AliTRDcluster *cl = 0x0;
1593 //////////////////////////////
1594 // Check no shared clusters
1595 //////////////////////////////
1596 for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
1597 cl = tracklet->GetClusters(icc);
1598 if(cl) crossrow = 1;
1600 //////////////////////////////////
1602 //////////////////////////////////
1603 for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1604 if(!(cl = tracklet->GetClusters(ic))) continue;
1605 //if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
1607 Double_t ycluster = cl->GetY();
1608 Int_t time = cl->GetPadTime();
1609 Double_t timeis = time/10.0;
1610 //See if cross two pad rows
1611 Int_t row = cl->GetPadRow();
1612 if(rowp==-1) rowp = row;
1613 if(row != rowp) crossrow = 1;
1615 linearFitterTracklet.AddPoint(&timeis,ycluster,1);
1621 ////////////////////////////////////
1622 // Do the straight line fit now
1623 ///////////////////////////////////
1625 linearFitterTracklet.ClearPoints();
1629 linearFitterTracklet.Eval();
1630 linearFitterTracklet.GetParameters(pars);
1631 pointError = TMath::Sqrt(linearFitterTracklet.GetChisquare()/(nbli-2));
1632 errorpar = linearFitterTracklet.GetParError(1)*pointError;
1634 //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",linearFitterTracklet->GetChisquare(),nbli,pointError,linearFitterTracklet->GetParError(1),errorpar);
1635 linearFitterTracklet.ClearPoints();
1637 /////////////////////////
1639 ////////////////////////
1641 if(nbclusters < fLow) return kFALSE;
1642 if(nbclusters > fHigh) return kFALSE;
1643 if(pointError >= 0.3) return kFALSE;
1644 if(crossrow == 1) return kTRUE;
1646 ///////////////////////
1648 //////////////////////
1651 //Add to the linear fitter of the detector
1652 if( TMath::Abs(snp) < 1.){
1653 Double_t x = tnp-dzdx*tnt;
1654 //if(!fLinearVdriftTest) printf("Not there\n");
1655 Double_t nbentries = fLinearVdriftTest->GetEntries();
1656 if(nbentries < (5.0*32767)) fLinearVdriftTest->Fill(x,dydt);