Fix compilation warnings
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibTime.cxx
CommitLineData
c74c5f6c 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/*
17 Comments to be written here:
18 1. What do we calibrate.
19
20Time dependence of gain and drift velocity in order to account for changes in: temperature, pressure, gas composition.
21
22 AliTPCcalibTime *calibTime = new AliTPCcalibTime("cosmicTime","cosmicTime",0, 1213.9e+06, 1213.96e+06, 0.04e+04, 0.04e+04);
23
24
25 2. How to interpret results
26 3. Simple example
27
28a.) determine the required time range:
29
30AliXRDPROOFtoolkit tool;
31TChain * chain = tool.MakeChain("pass2.txt","esdTree",0,6000);
32chain->Draw("GetTimeStamp()")
33
34b) analyse calibration object on Proof in calibration train
35
36AliTPCcalibTime *calibTime = new AliTPCcalibTime("cosmicTime","cosmicTime", StartTimeStamp, EndTimeStamp, IntegrationTimeVdrift, IntegrationTimeDeDx);
37
38c) plot results
39
40TFile f("CalibObjects.root");
828383be 41AliTPCcalibTime *cal = (AliTPCcalibTime *)f->Get("TPCCalib")->FindObject("cosmicTime");
c74c5f6c 42cal->GetHistVdrift()->Projection(1,0)->Draw()
43
44 4. Analysis using debug streamers.
45
46*/
47
48
49#include "Riostream.h"
50#include "TChain.h"
51#include "TTree.h"
52#include "TH1F.h"
53#include "TH2F.h"
54#include "TH3F.h"
55#include "THnSparse.h"
56#include "TList.h"
57#include "TMath.h"
58#include "TCanvas.h"
59#include "TFile.h"
60#include "TF1.h"
61#include "TVectorD.h"
62#include "TProfile.h"
63#include "TGraphErrors.h"
64#include "TCanvas.h"
65
66#include "AliTPCclusterMI.h"
67#include "AliTPCseed.h"
68#include "AliESDVertex.h"
69#include "AliESDEvent.h"
70#include "AliESDfriend.h"
71#include "AliESDInputHandler.h"
72#include "AliAnalysisManager.h"
73
74#include "AliTracker.h"
75#include "AliMagFMaps.h"
76#include "AliTPCCalROC.h"
77
78#include "AliLog.h"
79
80#include "AliTPCcalibTime.h"
81
82#include "TTreeStream.h"
83#include "AliTPCTracklet.h"
84
85ClassImp(AliTPCcalibTime)
86
87
88AliTPCcalibTime::AliTPCcalibTime()
89 :AliTPCcalibBase(),
90 fHistDeDxTgl(0),
91 fHistDeDx(0),
92 fHistVdrift(0),
93 fIntegrationTimeDeDx(0),
94 fIntegrationTimeVdrift(0),
95 fCutMaxD(5), // maximal distance in rfi ditection
96 fCutTheta(0.03), // maximal distan theta
97 fCutMinDir(-0.99) // direction vector products
98
99{
100 AliInfo("Default Constructor");
101}
102
103
104AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, ULong64_t TriggerMask, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeDeDx, Int_t deltaIntegrationTimeVdrift)
105 :AliTPCcalibBase(),
106 fTriggerMask(0),
107 fHistDeDxTgl(0),
108 fHistDeDx(0),
109 fHistVdrift(0),
110 fIntegrationTimeDeDx(0),
111 fIntegrationTimeVdrift(0),
112 fCutMaxD(5), // maximal distance in rfi ditection
113 fCutTheta(0.03), // maximal distan theta
114 fCutMinDir(-0.99) // direction vector products
115{
116
117 SetName(name);
118 SetTitle(title);
119
120 AliInfo("Non Default Constructor");
121
122 fTriggerMask = TriggerMask;
123
124 fIntegrationTimeDeDx = deltaIntegrationTimeDeDx;
125 fIntegrationTimeVdrift = deltaIntegrationTimeVdrift;
126
127 Double_t deltaTime = EndTime - StartTime;
128
828383be 129 Int_t binsVdrift[2] = {TMath::Nint(deltaTime/deltaIntegrationTimeVdrift), 100};
c74c5f6c 130 Double_t xminVdrift[2] = {StartTime, -20};
131 Double_t xmaxVdrift[2] = {EndTime, 20};
132 fHistVdrift = new THnSparseF("HistVdrift","vDrift; time;#Delta z",2,binsVdrift,xminVdrift,xmaxVdrift);
133
828383be 134 Int_t binsDeDxTgl[3] = {TMath::Nint(deltaTime/deltaIntegrationTimeDeDx),30,100};
c74c5f6c 135 Double_t xminDeDxTgl[3] = {StartTime,-1,0.7};
136 Double_t xmaxDeDxTgl[3] = {EndTime,1,1.3};
137 fHistDeDxTgl = new THnSparseF("HistDeDxTgl","dEdx vs tgl;time;tgl;dEdxUp/dEdxLow",3,binsDeDxTgl,xminDeDxTgl,xmaxDeDxTgl);
138
0ffd2ae1 139 Int_t binsDeDx[2] = {TMath::Nint(deltaTime/deltaIntegrationTimeDeDx),100};
c74c5f6c 140 Double_t xminDeDx[2] = {StartTime,1};
141 Double_t xmaxDeDx[2] = {EndTime,100};
142 fHistDeDx = new THnSparseF("HistDeDx","dEdx l;time;dEdx",2,binsDeDx,xminDeDx,xmaxDeDx);
143
144}
145
146
147
148AliTPCcalibTime::~AliTPCcalibTime(){
149 //
150 //
151 //
152}
153
154void AliTPCcalibTime::Process(AliESDEvent *event) {
155 //
156 //
157 //
158
159 ProcessCosmic(event);
160
161}
162
163
164
165void AliTPCcalibTime::ProcessCosmic(AliESDEvent *event) {
166
167 if (!event) {
168 Printf("ERROR: ESD not available");
169 return;
170 }
171 if (event->GetTimeStamp() == 0 ) {
172 Printf("no time stamp!");
173 return;
174 }
175
176 if (fTriggerMask != 0 && event->GetTriggerMask() != fTriggerMask) return;
177
178 UInt_t time = event->GetTimeStamp();
179
180 //
181 // Find cosmic pairs
182 //
183 // Track0 is choosen in upper TPC part
184 // Track1 is choosen in lower TPC part
185 //
186 Int_t ntracks=event->GetNumberOfTracks();
187 if (ntracks==0) return;
188 if (ntracks > 10) return; // temporary debug to remove LASER events
189
190
191 if (GetDebugLevel()>1) printf("Hallo world: Im here\n");
192 AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
193
194 TObjArray tpcSeeds(ntracks);
195 Double_t vtxx[3]={0,0,0};
196 Double_t svtxx[3]={0.000001,0.000001,100.};
197 AliESDVertex vtx(vtxx,svtxx);
198 //
199 // track loop
200 //
201 for (Int_t i=0;i<ntracks;++i) {
202 AliESDtrack *track = event->GetTrack(i);
203
204 const AliExternalTrackParam * trackIn = track->GetInnerParam();
205 const AliExternalTrackParam * trackOut = track->GetOuterParam();
206 if (!trackIn) continue;
207 if (!trackOut) continue;
208
209 AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
210 TObject *calibObject;
211 AliTPCseed *seed = 0;
212 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
213 if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
214 }
215 if (seed) {
216 tpcSeeds.AddAt(seed,i);
217 if (track->GetTPCNcls() > 50) {
218 Double_t meanP = 0.5*(trackIn->GetP() + trackOut->GetP());
219 Double_t TPCsignal = seed->CookdEdxNorm(0.0,0.6,1,0,159,0x0,kTRUE,kTRUE);
220 Double_t vecDeDx[2] = {time, TPCsignal};
221 if (meanP > 12) fHistDeDx->Fill(vecDeDx);
222 }
223 }
224 }
225
226 if (ntracks<2) return;
227 //
228 // Find pairs
229 //
230 for (Int_t i=0;i<ntracks;++i) {
231 AliESDtrack *track0 = event->GetTrack(i);
232 // track0 - choosen upper part
233 if (!track0) continue;
234 if (!track0->GetOuterParam()) continue;
235 if (track0->GetOuterParam()->GetAlpha()<0) continue;
236 Double_t d1[3];
237 track0->GetDirection(d1);
238 for (Int_t j=0;j<ntracks;++j) {
239 if (i==j) continue;
240 AliESDtrack *track1 = event->GetTrack(j);
241 //track 1 lower part
242 if (!track1) continue;
243 if (!track1->GetOuterParam()) continue;
244 if (track1->GetOuterParam()->GetAlpha()>0) continue;
245 //
246 Double_t d2[3];
247 track1->GetDirection(d2);
248
249 AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
250 AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
251 if (! seed0) continue;
252 if (! seed1) continue;
253 Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
254 Float_t d0 = track0->GetLinearD(0,0);
255 Float_t d1 = track1->GetLinearD(0,0);
256 //
257 // conservative cuts - convergence to be guarantied
258 // applying before track propagation
259 if (TMath::Abs(d0+d1)>fCutMaxD) continue; // distance to the 0,0
260 if (dir>fCutMinDir) continue; // direction vector product
261 Float_t bz = AliTracker::GetBz();
262 Float_t dvertex0[2]; //distance to 0,0
263 Float_t dvertex1[2]; //distance to 0,0
264 track0->GetDZ(0,0,0,bz,dvertex0);
265 track1->GetDZ(0,0,0,bz,dvertex1);
266 if (TMath::Abs(dvertex0[1])>250) continue;
267 if (TMath::Abs(dvertex1[1])>250) continue;
268 //
269 //
270 //
271 Float_t dmax = TMath::Max(TMath::Abs(d0),TMath::Abs(d1));
272 AliExternalTrackParam param0(*track0);
273 AliExternalTrackParam param1(*track1);
274 //
275 // Propagate using Magnetic field and correct fo material budget
276 //
277 AliTracker::PropagateTrackTo(&param0,dmax+1,0.0005,3,kTRUE);
278 AliTracker::PropagateTrackTo(&param1,dmax+1,0.0005,3,kTRUE);
279 //
280 // Propagate rest to the 0,0 DCA - z should be ignored
281 //
828383be 282 //Bool_t b0 = ;
283 param0.PropagateToDCA(&vtx,bz,1000);
284 //Bool_t b1 =
285 param1.PropagateToDCA(&vtx,bz,1000);
c74c5f6c 286 //
287 param0.GetDZ(0,0,0,bz,dvertex0);
288 param1.GetDZ(0,0,0,bz,dvertex1);
289 //
290 Double_t xyz0[3];//,pxyz0[3];
291 Double_t xyz1[3];//,pxyz1[3];
292 param0.GetXYZ(xyz0);
293 param1.GetXYZ(xyz1);
294 Bool_t isPair = IsPair(&param0,&param1);
295
296 Double_t z0 = track0->GetOuterParam()->GetZ();
297 Double_t z1 = track1->GetOuterParam()->GetZ();
298
299 Double_t z0inner = track0->GetInnerParam()->GetZ();
300 Double_t z1inner = track1->GetInnerParam()->GetZ();
301
302 if (isPair && z0>0 && z0inner>0 && z1<0 && z1inner<0) {
303 Double_t vecVdrift[2] = {time, param0.GetZ() - param1.GetZ()};
304
305 if (track0->GetTPCNcls() > 80) fHistVdrift->Fill(vecVdrift);
306 }
307 if (isPair && z0 > 0 && z1 > 0) {
308 if (track1->GetTPCNcls()> 110 && track0->GetTPCNcls()> 110 && seed1->CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE) > 0) {
309 Double_t vecDeDxTgl[3] = {time, track0->GetTgl(), seed0->CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE)/seed1->CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE)};
310 fHistDeDxTgl->Fill(vecDeDxTgl);
311 }
312 }
313
314 } // end 2nd order loop
315 } // end 1st order loop
316
317}
318
319
320void AliTPCcalibTime::Analyze() {
828383be 321 //
322 //
323 //
c74c5f6c 324 TH2D * hVdrift = GetHistVdrift()->Projection(1,0);
0ffd2ae1 325 hVdrift->Draw();
c74c5f6c 326}
327
328
329Long64_t AliTPCcalibTime::Merge(TCollection *li) {
330
331 TIterator* iter = li->MakeIterator();
332 AliTPCcalibTime* cal = 0;
333
334 while ((cal = (AliTPCcalibTime*)iter->Next())) {
335 if (!cal->InheritsFrom(AliTPCcalibTime::Class())) {
336 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
337 return -1;
338 }
339
340 // add histograms here...
341 fHistDeDxTgl->Add(cal->GetHistDeDxVsTgl());
342 fHistVdrift->Add(cal->GetHistVdrift());
343 fHistDeDx->Add(cal->GetHistDeDx());
344
345 }
346
347 return 0;
348
349}
350
351
352
353Bool_t AliTPCcalibTime::IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1){
354 //
355 //
356 /*
357 // 0. Same direction - OPOSITE - cutDir +cutT
358 TCut cutDir("cutDir","dir<-0.99")
359 // 1.
360 TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03")
361 //
362 // 2. The same rphi
363 TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5")
364 //
365 //
366 //
367 TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
368 // 1/Pt diff cut
369 */
370 const Double_t *p0 = tr0->GetParameter();
371 const Double_t *p1 = tr1->GetParameter();
372 if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
373 if (TMath::Abs(p0[0]+p1[0])>fCutMaxD) return kFALSE;
374 Double_t d0[3], d1[3];
375 tr0->GetDirection(d0);
376 tr1->GetDirection(d1);
377 if (d0[0]*d1[0] + d0[1]*d1[1] + d0[2]*d1[2] >fCutMinDir) return kFALSE;
378 //
379 return kTRUE;
380}