Bug fix in the constructor (thanks to A.Marin)
[u/mrichter/AliRoot.git] / TRD / AliTRDCalibTask.cxx
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>
29 using 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
75 #include "AliTRDcalibDB.h"
76 #include "AliCDBId.h"
77 #include "AliLog.h"
78
79
80 #include "AliTRDCalibTask.h"
81
82
83 ClassImp(AliTRDCalibTask)
84
85 //________________________________________________________________________
86   AliTRDCalibTask::AliTRDCalibTask(const char *name) 
87     : AliAnalysisTaskSE(name), fESD(0),
88       fESDfriend(0),
89       fkEsdTrack(0),
90       fFriendTrack(0),
91       fCalibObject(0),
92       fTrdTrack(0),
93       fCl(0),
94       fListHist(0),
95       fTRDCalibraFillHisto(0),
96       fNEvents(0),
97       fNEventsInput(0),
98       fNbTRDTrack(0),
99       fNbTRDTrackOffline(0),
100       fNbTRDTrackStandalone(0),
101       fNbTPCTRDtrack(0),
102       fNbGoodTracks(0),
103       fNbTimeBin(0),
104       fNbTimeBinOffline(0),
105       fNbTimeBinStandalone(0),
106       fNbClusters(0),
107       fNbClustersOffline(0),
108       fNbClustersStandalone(0),
109       fNbTracklets(0),
110       fNbTrackletsOffline(0),
111       fNbTrackletsStandalone(0),
112       fAbsoluteGain(0),
113       fCH2dSum(0),
114       fPH2dSum(0),
115       fCH2dSM(0),
116       fPH2dSM(0),
117       fCH2dTest(0),
118       fPH2dTest(0),
119       fLinearVdriftTest(0),
120       fOnInstance(kTRUE),
121       fHisto2d(kTRUE),
122       fVector2d(kFALSE),
123       fVdriftLinear(kTRUE),
124       fNbTimeBins(0),
125       fSelectedTrigger(new TObjArray()),
126       fRejected(kTRUE),
127       fEsdTrackCuts(0),
128       fRequirePrimaryVertex(kFALSE),
129       fVtxTPC(kFALSE),
130       fVtxSPD(kFALSE),
131       fMinNbContributors(0),
132       fRangePrimaryVertexZ(9999999.0),
133       fMinNbTracks(9),
134       fMaxNbTracks(500),
135       fLow(0),
136       fHigh(30),
137       fFillZero(kFALSE),
138       fNormalizeNbOfCluster(kFALSE),
139       fRelativeScale(0.0),
140       fMaxCluster(100.0),
141       fNbMaxCluster(2),
142       fOfflineTracks(kFALSE),
143       fStandaloneTracks(kFALSE),
144       fFirstRunGain(-1),
145       fVersionGainUsed(-1),
146       fSubVersionGainUsed(-1),
147       fFirstRunGainLocal(-1),
148       fVersionGainLocalUsed(-1),
149       fSubVersionGainLocalUsed(-1),
150       fFirstRunVdrift(-1),
151       fVersionVdriftUsed(-1), 
152       fSubVersionVdriftUsed(-1),
153       fCalDetGain(0x0),
154       fMaxEvent(0),
155       fCounter(0),
156       fDebug(0)
157 {
158   //
159   // Default constructor
160   //
161
162   fNz[0] = 0;
163   fNz[1] = 0;
164   fNz[2] = 0;
165   
166   fNrphi[0] = 0;
167   fNrphi[1] = 0;
168   fNrphi[2] = 0;
169
170   // Define input and output slots here
171   // Input slot #0 works with a TChain
172   DefineInput(0, TChain::Class());
173         
174   // Output slot #0 writes into a TList container
175   DefineOutput(1, TList::Class());
176   
177    
178 }
179 //____________________________________________________________________________________
180 AliTRDCalibTask::~AliTRDCalibTask()
181 {
182   //
183   // AliTRDCalibTask destructor
184   //
185
186   // Pointeur
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;
212   
213   if(fSelectedTrigger) {
214     fSelectedTrigger->Delete();
215     delete fSelectedTrigger;
216   }
217   if(fEsdTrackCuts) {
218     delete fEsdTrackCuts;
219   }
220   
221 }
222
223 //________________________________________________________________________
224 void AliTRDCalibTask::UserCreateOutputObjects() 
225 {
226   //
227   // Create the histos
228   //
229   //cout << "AliTRDCalibTask::CreateOutputObjects() IN" << endl;
230
231   // Number of time bins
232   if(fNbTimeBins==0) {
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));
237       fNbTimeBins = 30;
238     }
239   }
240
241   // output list
242   fListHist = new TList();
243   fListHist->SetOwner();
244
245   // instance calibration 
246   fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
247   if(fOnInstance) {
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
259       }
260       else {
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);
264         }
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);
268         }
269       }
270     }
271     // Variables for how to fill
272     fTRDCalibraFillHisto->SetFillWithZero(fFillZero);
273     fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fNormalizeNbOfCluster); 
274     fTRDCalibraFillHisto->SetMaxCluster(fMaxCluster);
275     fTRDCalibraFillHisto->SetNbMaxCluster(fNbMaxCluster);
276   
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
281     
282     // For testing only
283     if(fDebug > 2) fTRDCalibraFillHisto->SetDebugLevel(1); //debug stuff
284     
285     if(fHisto2d) {  
286       fListHist->Add(fTRDCalibraFillHisto->GetCH2d());
287       fListHist->Add(fTRDCalibraFillHisto->GetPH2d()); 
288       fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());
289     } 
290     if(fVdriftLinear) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetVdriftLinearFit());
291     if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector  
292   }
293   fRelativeScale = fTRDCalibraFillHisto->GetRelativeScale(); // Get the relative scale for the gain
294   
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);
299   
300   // absolute gain calibration even without AliESDfriend
301   Int_t nBinsPt = 25;
302   Double_t minPt = 0.001;
303   Double_t maxPt = 10.0;
304   
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]);
309   
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);
317   
318   /////////////////////////////////////////
319   // First debug level
320   ///////////////////////////////////////
321   if(fDebug > 0) {
322
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);
329     
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)
333                            ,540,0,540);
334     fPH2dTest->SetYTitle("Det/pad groups");
335     fPH2dTest->SetXTitle("time [#mus]");
336     fPH2dTest->SetZTitle("<PH> [a.u.]");
337     fPH2dTest->SetStats(0);
338     //
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);
344     fCH2dTest->Sumw2();
345
346     //
347     fPH2dSM = new TProfile2D(Form("PH2dSM_%s",(const char*)fName),"Nz10Nrphi10"
348                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
349                            ,18,0,18);
350     fPH2dSM->SetYTitle("Det/pad groups");
351     fPH2dSM->SetXTitle("time [#mus]");
352     fPH2dSM->SetZTitle("<PH> [a.u.]");
353     fPH2dSM->SetStats(0);
354     //
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);
360     fCH2dSM->Sumw2();
361     //
362     fPH2dSum = new TProfile2D(Form("PH2dSum_%s",(const char*)fName),"Nz100Nrphi100"
363                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
364                             ,1,0,1);
365     fPH2dSum->SetYTitle("Det/pad groups");
366     fPH2dSum->SetXTitle("time [#mus]");
367     fPH2dSum->SetZTitle("<PH> [a.u.]");
368     fPH2dSum->SetStats(0);
369     //
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);
375     fCH2dSum->Sumw2();
376     
377     
378     // Add them
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);
386
387   }
388
389   /////////////////////////////////////////
390   // Second debug level
391   ///////////////////////////////////////
392   if(fDebug > 1) {
393
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);
398
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();
407     //
408     fNbTimeBin = new TH1F(Form("NbTimeBin_%s",(const char*)fName),"NbTimeBin",35,0,35);
409     fNbTimeBin->Sumw2();
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();
414     //
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();
421     //
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();
428    
429     fListHist->Add(fNbGoodTracks);
430    
431     fListHist->Add(fNbTRDTrack);
432     fListHist->Add(fNbTRDTrackOffline);
433     fListHist->Add(fNbTRDTrackStandalone);
434     fListHist->Add(fNbTPCTRDtrack);
435     
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);
445     
446   }
447
448   delete [] binLimLogPt;
449   delete [] binLimPt;
450
451   PostData(1,fListHist);
452
453   //cout << "AliTRDCalibTask::UserCreateOutputObjects() OUT" << endl;
454
455 }
456
457 //________________________________________________________________________
458 void AliTRDCalibTask::UserExec(Option_t *) 
459 {
460   //
461   // Filling of the histos
462   //
463   //cout << "AliTRDCalibTask::Exec() IN" << endl;
464   
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);
469       return;
470     }
471   }
472   if(fCounter==0) {
473     if(fOnInstance) {
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();
484     }
485     if(fDebug > 1){
486       // title CH2dTest
487       TString name("Ver");
488       name += fVersionGainUsed;
489       name += "Subver";
490       name += fSubVersionGainUsed;
491       name += "FirstRun";
492       name += fFirstRunGain;
493       name += "Nz0Nrphi0";
494       fCH2dTest->SetTitle(name);  
495       // title PH2dTest
496       TString namee("Ver");
497       namee += fVersionVdriftUsed;
498       namee += "Subver";
499       namee += fSubVersionVdriftUsed;
500       namee += "FirstRun";
501       namee += fFirstRunVdrift;
502       namee += "Nz0Nrphi0";
503       fPH2dTest->SetTitle(namee); 
504     }
505   }
506   
507   //  AliLog::SetGlobalLogLevel(AliLog::kError);
508   //  cout << "AliTRDCalibTask::Exec() 1" << endl;
509   fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
510   if(!fESD){
511     AliError("ESD Event missing");
512     PostData(1, fListHist);
513     return;
514   }
515
516   const char* type = fESD->GetBeamType();
517   
518   
519   //printf("Counter %d\n",fCounter);
520   
521   fCounter++;
522   fNEventsInput->Fill(1);
523
524   //cout << "maxEvent = " << fMaxEvent << endl;
525   //if(fCounter%100==0) cout << "fCounter = " << fCounter << endl;
526   if((fMaxEvent != 0) && (fMaxEvent < fCounter)) {
527     PostData(1, fListHist);
528     return;
529   }
530   //if(fCounter%100==0) cout << "fCounter1 = " << fCounter << endl;
531   //cout << "event = " << fCounter << endl;
532   
533   //printf("Counter %d\n",fCounter);
534   
535   ///////////////////
536   // Check trigger
537   ///////////////////
538   Bool_t pass = kTRUE;
539
540   if (strstr(type,"p-p")) {
541    
542     //printf("Will check the triggers\n");
543
544     Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
545     //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
546     if(fRejected) {
547       pass = kTRUE;
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)) {
552           pass = kFALSE;
553         }
554       }
555     }
556     else {
557       pass = kFALSE;
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)) {
562           pass = kTRUE;
563         }
564       }
565     }
566     if(!pass) {
567       PostData(1, fListHist);
568       return;
569     }   
570
571   }
572     
573   //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
574   //printf("Trigger passed\n");
575   
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() ;
584     if(!vtxESD){
585       PostData(1, fListHist);
586       return;
587     }
588     Int_t nCtrb = vtxESD->GetNContributors();
589     if(nCtrb < fMinNbContributors) {
590       PostData(1, fListHist);     
591       return;
592     }
593     Double_t zPosition = vtxESD->GetZ();
594     if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) {
595       PostData(1, fListHist);
596       return;
597     }     
598     
599   }
600   
601   //printf("Primary vertex passed\n");
602   
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++;  
610   }
611   if(fDebug > 1)  {
612     // Centrality
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);
618   }
619   
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);
624       return;
625     }
626   }
627   
628   fNEvents->Fill(1);
629   
630   // In total
631   Int_t nbTrdTracks = 0;
632   // standalone
633   Int_t nbTrdTracksStandalone = 0;
634   // offline
635   Int_t nbTrdTracksOffline = 0;
636   // TPC
637   Int_t nbtrackTPC = 0;
638   
639
640   
641   if (nbTracks <= 0.0) {
642     
643     if(fDebug > 1) {
644       fNbTRDTrack->Fill(nbTrdTracks);
645       fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
646       fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
647     }
648     PostData(1, fListHist);
649     return;
650   }
651   
652   
653   fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
654   if(!fESDfriend){
655     AliError("fESDfriend not available");
656     PostData(1, fListHist);
657     return;
658   }
659
660   if(fESDfriend->TestSkipBit()) {
661     PostData(1, fListHist);
662     return;
663   }
664   
665   //printf("has friends\n");
666
667   /////////////////////////////////////
668   // Loop on AliESDtrack
669   ////////////////////////////////////
670   //printf("Nb of tracks %f\n",nbTracks);      
671   for(int itrk=0; itrk < nbTracks; ++itrk){
672     
673     // Get ESD track
674     fkEsdTrack = fESD->GetTrack(itrk);
675     if(!fkEsdTrack) continue;
676     ULong_t status = fkEsdTrack->GetStatus(); 
677     if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
678     
679     // Quality cuts on the AliESDtrack
680     if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
681       //printf("Not a good track\n");
682       continue;
683     }
684     
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)); 
698       }
699     }     
700     
701     // Other cuts
702     Bool_t good = kTRUE;
703     Bool_t standalonetrack = kFALSE;
704     Bool_t offlinetrack = kFALSE;
705     //ULong_t status = fkEsdTrack->GetStatus();
706     
707     fFriendTrack = fESDfriend->GetTrack(itrk);
708     if(!fFriendTrack)  {
709       //printf("No friend track %d\n",itrk);
710       continue;
711     }
712     //////////////////////////////////////
713     // Loop on calibration objects
714     //////////////////////////////////////
715     Int_t icalib=0;
716     Int_t nTRDtrackV1=0;
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");
721       ++nTRDtrackV1;
722       
723       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
724         standalonetrack = kTRUE;
725       }
726       if((status&(AliESDtrack::kTRDin))) {
727         offlinetrack = kTRUE;
728       }
729       if(fOfflineTracks){
730         if(!offlinetrack){
731           good = kFALSE;
732         }
733       }
734       else if(fStandaloneTracks){
735         if(!standalonetrack){
736           good = kFALSE;
737         }
738       }
739       
740       fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
741       if(good && fOnInstance) {
742         //cout << "good" << endl;
743         fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);
744         //printf("Fill fTRDCalibraFillHisto\n");
745       }
746       
747       //////////////////////////////////
748       // Debug 
749       ////////////////////////////////
750       
751       if(fDebug > 0) {
752         
753         //printf("Enter debug\n");
754         
755         Int_t nbtracklets = 0;
756         
757         // Check some stuff
758         Bool_t standalonetracklet = kFALSE;  
759         const AliTRDseedV1 *tracklet = 0x0;
760         //////////////////////////////////////
761         // Loop tracklets
762         ///////////////////////////////////// 
763         Int_t nbclusters=0;
764         Double_t phtb[AliTRDseedV1::kNtb];
765         memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
766         Double_t sum = 0.0;
767         Float_t normalisation = 6.67;
768         Int_t detector = 0;
769         Int_t sector = 0;
770         for(Int_t itr = 0; itr < 6; ++itr){
771           
772           if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
773           if(!tracklet->IsOK()) continue;
774           ++nbtracklets;
775           standalonetracklet = kFALSE; 
776           if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
777           
778           nbclusters = 0;
779           memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
780           sum = 0.0;
781           detector = 0;
782           sector = 0;
783           //Int_t crossrow = 0;
784           
785           // Check no shared clusters
786           //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
787           //  if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;
788           // }
789           
790           // Loop on clusters
791           Int_t time = 0;
792           Float_t ch = 0;
793           Float_t qcl = 0;
794           for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
795             
796             if(!(fCl = tracklet->GetClusters(ic))) continue;
797             ++nbclusters;
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");
807               }
808             }
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;
812             
813             if(fDebug > 1) {
814               fNbTimeBin->Fill(time);
815               if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
816               else fNbTimeBinOffline->Fill(time);
817             }
818           }
819           sector = AliTRDgeometry::GetSector(detector);
820           
821           if(fDebug > 1) {
822             fNbTracklets->Fill(detector);
823             if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
824             else fNbTrackletsOffline->Fill(detector);
825             
826             fNbClusters->Fill(nbclusters);
827             if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);
828             else  fNbClustersOffline->Fill(nbclusters);
829           }        
830           
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;
842                 }
843               }
844             }
845             if(checknoise) {           
846               for(int ic=0; ic<fNbTimeBins; ic++){
847                 if(fFillZero) {
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]);
851                 }
852                 else {
853                   if(phtb[ic] > 0.0) {
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]);
857                   }
858                 }
859               }
860             }
861           }
862           if(detector == 0) FindP1TrackPHtrackletV1Test(tracklet,nbclusters);
863
864         } // loop on tracklets
865
866         
867       } // debug
868       
869     }// while calibration objects
870     if(nTRDtrackV1 > 0) {
871       ++nbTrdTracks;      
872       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
873         ++nbTrdTracksStandalone;
874       }
875       if((status&(AliESDtrack::kTRDin))) {
876         ++nbTrdTracksOffline;
877       }
878     }
879     //delete fFriendTrack;
880   } // loop ESD track
881   
882   if(fDebug > 1) {
883     fNbTRDTrack->Fill(nbTrdTracks);
884     fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
885     fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
886     fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
887   }
888   
889   // Post output data
890   PostData(1, fListHist);
891   //cout << "AliTRDCalibTask::Exec() OUT" << endl;
892 }
893      
894 //________________________________________________________________________
895 void AliTRDCalibTask::Terminate(Option_t *) 
896 {
897   //
898   // Terminate
899   //
900   
901   if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
902
903  
904 }
905 //_______________________________________________________
906 Bool_t AliTRDCalibTask::Load(const Char_t *filename)
907 {
908   //
909   // Generic container loader
910   //
911
912   if(!TFile::Open(filename)){
913     //AliWarning(Form("Couldn't open file %s.", filename));
914     return kFALSE;
915   }
916   TList *o = 0x0;
917   if(!(o = (TList*)gFile->Get(GetName()))){
918     //AliWarning("Missing histogram container.");
919     return kFALSE;
920   }
921   fListHist = (TList*)o->Clone(GetName());
922   gFile->Close();
923   return kTRUE;
924 }
925 //_______________________________________________________
926 Bool_t AliTRDCalibTask::Load(TList *lister)
927 {
928   //
929   // Generic container loader
930   //
931
932   fListHist = (TList*)lister->Clone(GetName());
933   return kTRUE;
934 }
935 //________________________________________________________________________
936 void AliTRDCalibTask::Plot() 
937 {
938   //
939   // Plot the histos stored in the TList
940   //
941  
942   if(!fListHist) return;
943
944   /////////////////////////////////////
945   // Take the debug stuff
946   /////////////////////////////////////
947
948   TH1I *nEvents  = (TH1I *) fListHist->FindObject("NEvents");
949
950   TH2F *absoluteGain  = (TH2F *) fListHist->FindObject("AbsoluteGain");
951
952   TH1F *trdTrack = (TH1F *) fListHist->FindObject("TRDTrack");
953   TH1F *trdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");
954   TH1F *trdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");
955
956   TH2F *tpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");
957
958   TH1F *nbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");
959   TH1F *nbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");
960   TH1F *nbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");
961
962   TH1F *nbClusters = (TH1F *) fListHist->FindObject("NbClusters");
963   TH1F *nbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");
964   TH1F *nbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");
965
966   TH1F *nbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");
967   TH1F *nbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");
968   TH1F *nbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");
969   
970   /////////////////////////////////////
971   // Take the calibration objects
972   /////////////////////////////////////
973
974   TH2I *ch2d = (TH2I *) fListHist->FindObject("CH2d");
975   TProfile2D *ph2d = (TProfile2D *) fListHist->FindObject("PH2d");
976
977   TH2I *ch2dSum = (TH2I *) fListHist->FindObject("CH2dSum");
978   TProfile2D *ph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");
979
980   TH2I *ch2dSM = (TH2I *) fListHist->FindObject("CH2dSM");
981   TProfile2D *ph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");
982   
983   AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
984   
985   ////////////////////////////////////////////////
986   // Add the AliTRDCalibraVdriftLinearFit
987   ///////////////////////////////////////////////
988   
989   Int_t first = 0;
990   TH2S *histolinearfitsum = 0x0;
991   
992   if(linearfit) {
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);
997           first += 1;
998         }
999         else {
1000           if (histolinearfitsum) {
1001             histolinearfitsum->Add(linearfit->GetLinearFitterHisto(det));
1002           }
1003         }
1004       }
1005     }
1006   }
1007
1008   ///////////////////////////
1009   // Style
1010   //////////////////////////
1011
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);
1019
1020   /////////////////////////
1021   // Plot
1022   ////////////////////////
1023
1024  if(nEvents) {
1025    
1026     TCanvas *debugEvents = new TCanvas("cNEvents","cNEvents",10,10,510,510);
1027     debugEvents->cd(1);
1028     if(nEvents) nEvents->Draw();
1029       
1030   }
1031
1032  if(absoluteGain) {
1033    
1034     TCanvas *debugAbsoluteGain = new TCanvas("cAbsoluteGain","cAbsoluteGain",10,10,510,510);
1035     debugAbsoluteGain->cd(1);
1036     if(absoluteGain) absoluteGain->Draw();
1037       
1038   }
1039
1040   if(trdTrack || tpctrdTrack) {
1041     
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");
1052     leg->Draw("same");
1053     debugtrdtpcTrack->cd(2);
1054     if(tpctrdTrack) tpctrdTrack->Draw();
1055     TLine *line = new TLine(0.0,0.0,100.0,100.0);
1056     line->Draw("same");
1057     
1058   }
1059  
1060   if(nbTimeBin || nbTracklets || nbClusters) {
1061     
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");
1072     lega->Draw("same");
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");
1081     legb->Draw("same");
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");
1090     legc->Draw("same");
1091   
1092   }
1093
1094   if(ch2dSum || ph2dSum || histolinearfitsum) {
1095     
1096     TCanvas *debugSum = new TCanvas("SumCalibrationObjects","SumCalibrationObjects",10,10,510,510);
1097     debugSum->Divide(3,1);
1098     debugSum->cd(1);
1099     if(ch2dSum) ch2dSum->Draw("lego");
1100     debugSum->cd(2);
1101     if(ph2dSum) ph2dSum->Draw("lego");
1102     debugSum->cd(3);
1103     if(histolinearfitsum) histolinearfitsum->Draw();
1104   
1105   }
1106
1107   if(ch2dSM || ph2dSM) {
1108     
1109     TCanvas *debugSM = new TCanvas("SMCalibrationObjects","SMCalibrationObjects",10,10,510,510);
1110     debugSM->Divide(2,1);
1111     debugSM->cd(1);
1112     if(ch2dSM) ch2dSM->Draw("lego");
1113     debugSM->cd(2);
1114     if(ph2dSM) ph2dSM->Draw("lego");
1115      
1116   }
1117
1118   if(ch2d || ph2d) {
1119     
1120     TCanvas *debug = new TCanvas("CalibrationObjects","CalibrationObjects",10,10,510,510);
1121     debug->Divide(2,1);
1122     debug->cd(1);
1123     if(ch2d) ch2d->Draw("lego");
1124     debug->cd(2);
1125     if(ph2d) ph2d->Draw("lego");
1126      
1127   }
1128  
1129 }
1130 //_______________________________________________________________________________________
1131 void  AliTRDCalibTask::AddTask(const AliTRDCalibTask * calibTask) {
1132
1133   //
1134   // Add stats
1135   //
1136
1137   TList *listcalibTask = calibTask->GetList();
1138   if(!listcalibTask) return;
1139
1140   TH1I *nEvents  = (TH1I *) listcalibTask->FindObject("NEvents");
1141   TH1I *nEventsInput  = (TH1I *) listcalibTask->FindObject("NEventsInput");
1142   TH2F *absoluteGain  = (TH2F *) listcalibTask->FindObject("AbsoluteGain");
1143
1144   TH1F *trdTrack = (TH1F *) listcalibTask->FindObject("TRDTrack");
1145   TH1F *trdTrackOffline = (TH1F *) listcalibTask->FindObject("TRDTrackOffline");
1146   TH1F *trdTrackStandalone = (TH1F *) listcalibTask->FindObject("TRDTrackStandalone");
1147
1148   TH2F *tpctrdTrack = (TH2F *) listcalibTask->FindObject("NbTPCTRDtrack");
1149
1150   TH1F *nbTimeBin = (TH1F *) listcalibTask->FindObject("NbTimeBin");
1151   TH1F *nbTimeBinOffline = (TH1F *) listcalibTask->FindObject("NbTimeBinOffline");
1152   TH1F *nbTimeBinStandalone = (TH1F *) listcalibTask->FindObject("NbTimeBinStandalone");
1153
1154   TH1F *nbClusters = (TH1F *) listcalibTask->FindObject("NbClusters");
1155   TH1F *nbClustersOffline = (TH1F *) listcalibTask->FindObject("NbClustersOffline");
1156   TH1F *nbClustersStandalone = (TH1F *) listcalibTask->FindObject("NbClustersStandalone");
1157
1158   TH1F *nbTracklets = (TH1F *) listcalibTask->FindObject("NbTracklets");
1159   TH1F *nbTrackletsOffline = (TH1F *) listcalibTask->FindObject("NbTrackletsOffline");
1160   TH1F *nbTrackletsStandalone = (TH1F *) listcalibTask->FindObject("NbTrackletsStandalone");
1161   
1162   TH2I *ch2d = (TH2I *) listcalibTask->FindObject("CH2d");
1163   TProfile2D *ph2d = (TProfile2D *) listcalibTask->FindObject("PH2d");
1164   TProfile2D *prf2d = (TProfile2D *) listcalibTask->FindObject("PRF2d");
1165
1166   TH2I *ch2dSum = (TH2I *) listcalibTask->FindObject("CH2dSum");
1167   TProfile2D *ph2dSum = (TProfile2D *) listcalibTask->FindObject("PH2dSum");
1168
1169   TH2I *ch2dSM = (TH2I *) listcalibTask->FindObject("CH2dSM");
1170   TProfile2D *ph2dSM = (TProfile2D *) listcalibTask->FindObject("PH2dSM");
1171   
1172   AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) listcalibTask->FindObject("AliTRDCalibraVdriftLinearFit");  
1173   AliTRDCalibraVector *calibraVector = (AliTRDCalibraVector *) listcalibTask->FindObject("AliTRDCalibraVector");  
1174
1175   //
1176
1177   TH1I *inEventsInput  = (TH1I *) fListHist->FindObject("NEventsInput");
1178   TH1I *inEvents  = (TH1I *) fListHist->FindObject("NEvents");
1179   TH2F *iabsoluteGain  = (TH2F *) fListHist->FindObject("AbsoluteGain");
1180
1181   TH1F *itrdTrack = (TH1F *) fListHist->FindObject("TRDTrack");
1182   TH1F *itrdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");
1183   TH1F *itrdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");
1184
1185   TH2F *itpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");
1186
1187   TH1F *inbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");
1188   TH1F *inbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");
1189   TH1F *inbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");
1190
1191   TH1F *inbClusters = (TH1F *) fListHist->FindObject("NbClusters");
1192   TH1F *inbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");
1193   TH1F *inbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");
1194
1195   TH1F *inbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");
1196   TH1F *inbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");
1197   TH1F *inbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");
1198   
1199   TH2I *ich2d = (TH2I *) fListHist->FindObject("CH2d");
1200   TProfile2D *iph2d = (TProfile2D *) fListHist->FindObject("PH2d");
1201   TProfile2D *iprf2d = (TProfile2D *) fListHist->FindObject("PRF2d");
1202
1203   TH2I *ich2dSum = (TH2I *) fListHist->FindObject("CH2dSum");
1204   TProfile2D *iph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");
1205
1206   TH2I *ich2dSM = (TH2I *) fListHist->FindObject("CH2dSM");
1207   TProfile2D *iph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");
1208   
1209   AliTRDCalibraVdriftLinearFit *ilinearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");  
1210   AliTRDCalibraVector *icalibraVector = (AliTRDCalibraVector *) fListHist->FindObject("AliTRDCalibraVector");  
1211
1212
1213   // Add
1214
1215   if(nEventsInput) {
1216     if(inEventsInput) {
1217       inEventsInput->Add(nEventsInput);
1218       //printf("Add Events\n");
1219     }
1220     else {
1221       //printf("Create new Events\n");
1222       inEventsInput = new TH1I(*nEventsInput);
1223       fListHist->Add(inEventsInput);
1224     }
1225   }
1226   
1227   if(nEvents) {
1228     if(inEvents) {
1229       inEvents->Add(nEvents);
1230       //printf("Add Events\n");
1231     }
1232     else {
1233       //printf("Create new Events\n");
1234       inEvents = new TH1I(*nEvents);
1235       fListHist->Add(inEvents);
1236     }
1237   }
1238   
1239   if(absoluteGain) {
1240     if(iabsoluteGain) iabsoluteGain->Add(absoluteGain);
1241     else {
1242       iabsoluteGain = new TH2F(*absoluteGain);
1243       fListHist->Add(iabsoluteGain);
1244     }
1245   }
1246   
1247   if(trdTrack) {
1248     if(itrdTrack) itrdTrack->Add(trdTrack);
1249     else {
1250      itrdTrack = new TH1F(*trdTrack);
1251      fListHist->Add(itrdTrack);
1252     }
1253   }
1254
1255   if(trdTrackOffline) {
1256     if(itrdTrackOffline) itrdTrackOffline->Add(trdTrackOffline);
1257     else {
1258       itrdTrackOffline = new TH1F(*trdTrackOffline);
1259       fListHist->Add(itrdTrackOffline);
1260     }
1261   }
1262
1263   if(trdTrackStandalone) {
1264     if(itrdTrackStandalone) itrdTrackStandalone->Add(trdTrackStandalone);
1265     else {
1266       itrdTrackStandalone = new TH1F(*trdTrackStandalone);
1267       fListHist->Add(itrdTrackStandalone);
1268     }
1269   }
1270
1271   if(tpctrdTrack) {
1272     if(itpctrdTrack) itpctrdTrack->Add(tpctrdTrack);
1273     else {
1274       itpctrdTrack = new TH2F(*tpctrdTrack);
1275       fListHist->Add(itpctrdTrack);
1276     }
1277   }
1278
1279   if(nbTimeBin) {
1280     if(inbTimeBin) inbTimeBin->Add(nbTimeBin);
1281     else {
1282       inbTimeBin = new TH1F(*inbTimeBin);
1283       fListHist->Add(inbTimeBin);
1284     }
1285   }
1286
1287   if(nbTimeBinOffline) {
1288     if(inbTimeBinOffline) inbTimeBinOffline->Add(nbTimeBinOffline);
1289     else {
1290       inbTimeBinOffline = new TH1F(*nbTimeBinOffline);
1291       fListHist->Add(inbTimeBinOffline);
1292     }
1293   }
1294   
1295   if(nbTimeBinStandalone) {
1296     if(inbTimeBinStandalone) inbTimeBinStandalone->Add(nbTimeBinStandalone);
1297     else {
1298       inbTimeBinStandalone = new TH1F(*nbTimeBinStandalone);
1299       fListHist->Add(inbTimeBinStandalone);
1300     }
1301   }
1302
1303   if(nbClusters) {
1304     if(inbClusters) inbClusters->Add(nbClusters);
1305     else {
1306       inbClusters = new TH1F(*nbClusters);
1307       fListHist->Add(inbClusters);
1308     }
1309   }
1310   
1311   if(nbClustersOffline) {
1312     if(inbClustersOffline) inbClustersOffline->Add(nbClustersOffline);
1313     else {
1314       inbClustersOffline = new TH1F(*nbClustersOffline);
1315       fListHist->Add(inbClustersOffline);
1316     }
1317   }
1318   
1319   if(nbClustersStandalone) {
1320     if(inbClustersStandalone) inbClustersStandalone->Add(nbClustersStandalone);
1321     else {
1322       inbClustersStandalone = new TH1F(*nbClustersStandalone);
1323       fListHist->Add(inbClustersStandalone);
1324     }
1325   }
1326
1327   if(nbTracklets) {
1328     if(inbTracklets) inbTracklets->Add(nbTracklets);
1329     else {
1330       inbTracklets = new TH1F(*nbTracklets);
1331       fListHist->Add(inbTracklets);
1332     }
1333   }
1334
1335   if(nbTrackletsOffline) {
1336     if(inbTrackletsOffline) inbTrackletsOffline->Add(nbTrackletsOffline);
1337     else {
1338       inbTrackletsOffline = new TH1F(*nbTrackletsOffline);
1339       fListHist->Add(inbTrackletsOffline);
1340     }
1341   }
1342   
1343   if(nbTrackletsStandalone) {
1344     if(inbTrackletsStandalone) inbTrackletsStandalone->Add(nbTrackletsStandalone);
1345     else {
1346       inbTrackletsStandalone = new TH1F(*nbTrackletsStandalone);
1347       fListHist->Add(inbTrackletsStandalone);
1348     }
1349   }
1350   
1351   if(ch2d) {
1352     if(ich2d) ich2d->Add(ch2d);
1353     else {
1354       ich2d = new TH2I(*ch2d);
1355       fListHist->Add(ich2d);
1356     }
1357   }
1358
1359   if(ph2d) {
1360     if(iph2d) iph2d->Add(ph2d);
1361     else {
1362       iph2d = new TProfile2D(*ph2d);
1363       fListHist->Add(iph2d);
1364     }
1365   }
1366
1367   if(prf2d) {
1368     if(iprf2d) iprf2d->Add(prf2d);
1369     else {
1370       iprf2d = new TProfile2D(*prf2d);
1371       fListHist->Add(iprf2d);
1372     }
1373   }
1374
1375   if(ch2dSum) {
1376     if(ich2dSum) ich2dSum->Add(ch2dSum);
1377     else {
1378       ich2dSum = new TH2I(*ch2dSum);
1379       fListHist->Add(ich2dSum);
1380     }
1381   }
1382
1383   if(ph2dSum) {
1384     if(iph2dSum) iph2dSum->Add(ph2dSum);
1385     else {
1386       iph2dSum = new TProfile2D(*ph2dSum);
1387       fListHist->Add(iph2dSum);
1388     }
1389   }
1390
1391   if(ch2dSM) {
1392     if(ich2dSM) ich2dSM->Add(ch2dSM);
1393     else {
1394       ich2dSM = new TH2I(*ch2dSM);
1395       fListHist->Add(ich2dSM);
1396     }
1397   }
1398
1399   if(ph2dSM) {
1400     if(iph2dSM) iph2dSM->Add(ph2dSM);
1401     else {
1402       iph2dSM = new TProfile2D(*ph2dSM);
1403       fListHist->Add(iph2dSM);
1404     }
1405   }
1406   
1407   if(linearfit) {
1408     if(ilinearfit) ilinearfit->Add(linearfit);
1409     else {
1410       ilinearfit = new AliTRDCalibraVdriftLinearFit(*linearfit);
1411       fListHist->Add(ilinearfit);
1412     }
1413   }
1414
1415   if(calibraVector) {
1416     if(icalibraVector) icalibraVector->Add(calibraVector);
1417     else {
1418       icalibraVector = new AliTRDCalibraVector(*calibraVector);
1419       fListHist->Add(icalibraVector);
1420     }
1421   }
1422   
1423 }
1424 //________________________________________________________________________________
1425 Long64_t AliTRDCalibTask::Merge(TCollection *li) {
1426   
1427   //
1428   // merge component
1429   //
1430   
1431   TIterator* iter = li->MakeIterator();
1432   AliTRDCalibTask* cal = 0;
1433
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());
1437       return -1;
1438     }
1439
1440     // add histograms here...
1441     this->AddTask(cal);
1442     
1443   }
1444   
1445   return 0;
1446   
1447 }
1448 //_____________________________________________________
1449 Bool_t AliTRDCalibTask::SetVersionSubversion(){
1450   //
1451   // Load Chamber Gain factors into the Tender supply
1452   //
1453   
1454   printf("SetVersionSubversion\n");
1455
1456   //find previous entry from the UserInfo
1457   TTree *tree=((TChain*)GetInputData(0))->GetTree();
1458   if (!tree) {
1459     AliError("Tree not found in ESDhandler");
1460     return kFALSE;
1461   }
1462          
1463   TList *userInfo=(TList*)tree->GetUserInfo();
1464   if (!userInfo) {
1465     AliError("No UserInfo found in tree");
1466     return kFALSE;
1467   }
1468
1469   TList *cdbList=(TList*)userInfo->FindObject("cdbList");
1470   if (!cdbList) {
1471     AliError("No cdbList found in UserInfo");
1472     if (AliLog::GetGlobalLogLevel()>=AliLog::kError) userInfo->Print();
1473     return kFALSE;
1474   }
1475         
1476   TIter nextCDB(cdbList);
1477   TObjString *os=0x0;
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();
1497     }
1498   }
1499
1500   //printf("VersionGain %d, SubversionGain %d, VersionLocalGain %d, Subversionlocalgain %d, Versionvdrift %d, Subversionvdrift %d\n",fVersionGainUsed,fSubVersionGainUsed,fVersionGainLocalUsed,fSubVersionGainLocalUsed,fVersionVdriftUsed,fSubVersionVdriftUsed);
1501
1502   // Check
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");
1513     return kFALSE;
1514   }
1515   else return kTRUE;
1516
1517 }
1518 //_________________________________________________________________________________________________________________________
1519 Bool_t AliTRDCalibTask::ParticleGood(int i) const {
1520
1521   //
1522   // Definition of good tracks
1523   //
1524
1525   
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
1530   Float_t r,z;
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
1534   if (r==0) return 0;
1535   return 1;   
1536
1537
1538 }
1539 //______________________________________________________________________________________________________________________
1540 Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1541 {
1542   //
1543   // Drift velocity calibration:
1544   // Fit the clusters with a straight line
1545   // From the slope find the drift velocity
1546   //
1547
1548   ////////////////////////////////////////////////
1549   //Number of points: if less than 3 return kFALSE
1550   /////////////////////////////////////////////////
1551   if(nbclusters <= 2) return kFALSE;
1552
1553   ////////////
1554   //Variables
1555   ////////////
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);  
1566  
1567   
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));
1576   } 
1577   Double_t tgl  = tracklet->GetTgl();           // dz/dl
1578   Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp);   // dz/dx calculated from dz/dl
1579   // at the entrance
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 
1586
1587
1588   ////////////////////////////
1589   // loop over the clusters
1590   ////////////////////////////
1591   Int_t  nbli = 0;
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;
1599   }
1600   //////////////////////////////////
1601   // Loop clusters
1602   //////////////////////////////////
1603   for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1604     if(!(cl = tracklet->GetClusters(ic))) continue;
1605     //if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
1606     
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;
1614
1615     linearFitterTracklet.AddPoint(&timeis,ycluster,1);
1616     nbli++;  
1617
1618     
1619   }
1620   
1621   ////////////////////////////////////
1622   // Do the straight line fit now
1623   ///////////////////////////////////
1624   if(nbli <= 2){ 
1625     linearFitterTracklet.ClearPoints();  
1626     return kFALSE; 
1627   }
1628   TVectorD pars;
1629   linearFitterTracklet.Eval();
1630   linearFitterTracklet.GetParameters(pars);
1631   pointError  =  TMath::Sqrt(linearFitterTracklet.GetChisquare()/(nbli-2));
1632   errorpar    =  linearFitterTracklet.GetParError(1)*pointError;
1633   dydt        = pars[1]; 
1634   //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",linearFitterTracklet->GetChisquare(),nbli,pointError,linearFitterTracklet->GetParError(1),errorpar);
1635   linearFitterTracklet.ClearPoints();  
1636  
1637   /////////////////////////
1638   // Cuts quality
1639   ////////////////////////
1640   
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;
1645   
1646   ///////////////////////
1647   // Fill
1648   //////////////////////
1649
1650   if(fDebug > 0){
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);
1657     }
1658   }
1659   
1660   return kTRUE;
1661 }
1662