Possibilty to use median instead of the sigma of gaus fit
[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);
57dc06f2 37s
c74c5f6c 38c) plot results
2be25cc9 39 .x ~/NimStyle.C
40 gSystem->Load("libANALYSIS");
41 gSystem->Load("libTPCcalib");
c74c5f6c 42
43TFile f("CalibObjects.root");
2be25cc9 44AliTPCcalibTime *cal = (AliTPCcalibTime *)f->Get("TPCCalib")->FindObject("calibTime");
c74c5f6c 45cal->GetHistVdrift()->Projection(1,0)->Draw()
46
47 4. Analysis using debug streamers.
48
da6c0bc9 49 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
50 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
51 AliXRDPROOFtoolkit tool;
52 TChain * chainTime = tool.MakeChain("time.txt","timeInfo",0,10200);
53 chainTime->Lookup();
54
55
c74c5f6c 56*/
57
58
59#include "Riostream.h"
60#include "TChain.h"
61#include "TTree.h"
62#include "TH1F.h"
63#include "TH2F.h"
64#include "TH3F.h"
65#include "THnSparse.h"
66#include "TList.h"
67#include "TMath.h"
68#include "TCanvas.h"
69#include "TFile.h"
70#include "TF1.h"
71#include "TVectorD.h"
72#include "TProfile.h"
73#include "TGraphErrors.h"
74#include "TCanvas.h"
75
76#include "AliTPCclusterMI.h"
77#include "AliTPCseed.h"
78#include "AliESDVertex.h"
79#include "AliESDEvent.h"
80#include "AliESDfriend.h"
81#include "AliESDInputHandler.h"
82#include "AliAnalysisManager.h"
83
84#include "AliTracker.h"
f7a1cc68 85#include "AliMagF.h"
c74c5f6c 86#include "AliTPCCalROC.h"
87
88#include "AliLog.h"
89
90#include "AliTPCcalibTime.h"
91
92#include "TTreeStream.h"
93#include "AliTPCTracklet.h"
da6c0bc9 94#include "TTimeStamp.h"
95#include "AliTPCcalibDB.h"
96#include "AliTPCcalibLaser.h"
31aa7c5c 97#include "AliDCSSensorArray.h"
98#include "AliDCSSensor.h"
c74c5f6c 99
100ClassImp(AliTPCcalibTime)
101
102
103AliTPCcalibTime::AliTPCcalibTime()
da6c0bc9 104 :AliTPCcalibBase(),
c74c5f6c 105 fHistDeDxTgl(0),
106 fHistDeDx(0),
107 fHistVdrift(0),
108 fIntegrationTimeDeDx(0),
da6c0bc9 109 fIntegrationTimeVdrift(0),
110 fLaser(0), // pointer to laser calibration
111 fDz(0), // current delta z
112 fdEdx(0), // current dEdx
113 fdEdxRatio(0), // current dEdx ratio
114 fTl(0), // current tan(lambda)
c74c5f6c 115 fCutMaxD(5), // maximal distance in rfi ditection
2be25cc9 116 fCutMaxDz(20), // maximal distance in rfi ditection
c74c5f6c 117 fCutTheta(0.03), // maximal distan theta
2be25cc9 118 fCutMinDir(-0.99), // direction vector products
119 fMapDz(0), //NEW! Tmap of V drifts for different triggers
120 fTimeBins(0),
121 fTimeStart(0),
122 fTimeEnd(0),
123 fPtBins(0),
124 fPtStart(0),
125 fPtEnd(0),
126 fVdriftBins(0),
127 fVdriftStart(0),
128 fVdriftEnd(0),
129 fRunBins(0),
130 fRunStart(0),
131 fRunEnd(0)
132// fBinsVdrift(fTimeBins,fPtBins,fVdriftBins),
133// fXminVdrift(fTimeStart,fPtStart,fVdriftStart),
134// fXmaxVdrift(fTimeEnd,fPtEnd,fVdriftEnd)
c74c5f6c 135
136{
137 AliInfo("Default Constructor");
2be25cc9 138 for (Int_t i=0;i<3;i++) {
139 fHistVdriftLaserA[i]=0;
140 fHistVdriftLaserC[i]=0;
141 }
c74c5f6c 142}
143
144
2be25cc9 145AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeDeDx, Int_t deltaIntegrationTimeVdrift)
c74c5f6c 146 :AliTPCcalibBase(),
c74c5f6c 147 fHistDeDxTgl(0),
148 fHistDeDx(0),
149 fHistVdrift(0),
150 fIntegrationTimeDeDx(0),
da6c0bc9 151 fIntegrationTimeVdrift(0),
152 fLaser(0), // pointer to laser calibration
153 fDz(0), // current delta z
154 fdEdx(0), // current dEdx
155 fdEdxRatio(0), // current dEdx ratio
156 fTl(0), // current tan(lambda)
c74c5f6c 157 fCutMaxD(5), // maximal distance in rfi ditection
2be25cc9 158 fCutMaxDz(20), // maximal distance in rfi ditection
c74c5f6c 159 fCutTheta(0.03), // maximal distan theta
2be25cc9 160 fCutMinDir(-0.99), // direction vector products
161 fMapDz(0), //NEW! Tmap of V drifts for different triggers
162 fTimeBins(0),
163 fTimeStart(0),
164 fTimeEnd(0),
165 fPtBins(0),
166 fPtStart(0),
167 fPtEnd(0),
168 fVdriftBins(0),
169 fVdriftStart(0),
170 fVdriftEnd(0),
171 fRunBins(0),
172 fRunStart(0),
173 fRunEnd(0)
c74c5f6c 174{
175
176 SetName(name);
177 SetTitle(title);
2be25cc9 178 for (Int_t i=0;i<3;i++) {
179 fHistVdriftLaserA[i]=0;
180 fHistVdriftLaserC[i]=0;
181 }
c74c5f6c 182
183 AliInfo("Non Default Constructor");
184
c74c5f6c 185
186 fIntegrationTimeDeDx = deltaIntegrationTimeDeDx;
187 fIntegrationTimeVdrift = deltaIntegrationTimeVdrift;
188
189 Double_t deltaTime = EndTime - StartTime;
190
828383be 191 Int_t binsVdrift[2] = {TMath::Nint(deltaTime/deltaIntegrationTimeVdrift), 100};
c74c5f6c 192 Double_t xminVdrift[2] = {StartTime, -20};
193 Double_t xmaxVdrift[2] = {EndTime, 20};
194 fHistVdrift = new THnSparseF("HistVdrift","vDrift; time;#Delta z",2,binsVdrift,xminVdrift,xmaxVdrift);
195
2be25cc9 196 Int_t binsDeDxTgl[4] = {TMath::Nint(deltaTime/deltaIntegrationTimeDeDx),30,100,10000};
197 Double_t xminDeDxTgl[4] = {StartTime,-1,0.7,-0.5};
198 Double_t xmaxDeDxTgl[4] = {EndTime,1,1.3,999.5};
c74c5f6c 199 fHistDeDxTgl = new THnSparseF("HistDeDxTgl","dEdx vs tgl;time;tgl;dEdxUp/dEdxLow",3,binsDeDxTgl,xminDeDxTgl,xmaxDeDxTgl);
200
0ffd2ae1 201 Int_t binsDeDx[2] = {TMath::Nint(deltaTime/deltaIntegrationTimeDeDx),100};
c74c5f6c 202 Double_t xminDeDx[2] = {StartTime,1};
203 Double_t xmaxDeDx[2] = {EndTime,100};
204 fHistDeDx = new THnSparseF("HistDeDx","dEdx l;time;dEdx",2,binsDeDx,xminDeDx,xmaxDeDx);
205
2be25cc9 206// fLaser = new AliTPCcalibLaser("laserTPC","laserTPC",kFALSE);
207
208 fTimeBins =TMath::Nint(deltaTime/deltaIntegrationTimeVdrift);
209 fTimeStart =StartTime; //(((TObjString*)(mapGRP->GetValue("fAliceStartTime")))->GetString()).Atoi();
210 fTimeEnd =EndTime; //(((TObjString*)(mapGRP->GetValue("fAliceStopTime")))->GetString()).Atoi();
211 fPtBins = 200;
212 fPtStart = -0.04;
213 fPtEnd = 0.04;
214 fVdriftBins = 200;
215 fVdriftStart= -20.0/500.0;
216 fVdriftEnd = 20.0/500.0;
217 fRunBins = 10000;
218 fRunStart = -0.5;
219 fRunEnd = 0.5;
220
221 Int_t binsVdriftLaser[4] = {fTimeBins , fPtBins , fVdriftBins*20, fRunBins };
222 Double_t xminVdriftLaser[4] = {fTimeStart, fPtStart, fVdriftStart , fRunStart};
223 Double_t xmaxVdriftLaser[4] = {fTimeEnd , fPtEnd , fVdriftEnd , fRunEnd };
224
225 for (Int_t i=0;i<3;i++) {
226 fHistVdriftLaserA[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
227 fHistVdriftLaserC[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
228 }
229 fBinsVdrift[0] = fTimeBins;
230 fBinsVdrift[1] = fPtBins;
231 fBinsVdrift[2] = fVdriftBins;
232 fBinsVdrift[3] = fRunBins;
233 fXminVdrift[0] = fTimeStart;
234 fXminVdrift[1] = fPtStart;
235 fXminVdrift[2] = fVdriftStart;
236 fXminVdrift[3] = fRunStart;
237 fXmaxVdrift[0] = fTimeEnd;
238 fXmaxVdrift[1] = fPtEnd;
239 fXmaxVdrift[2] = fVdriftEnd;
240 fXmaxVdrift[3] = fRunEnd;
241
242 fMapDz=new TMap();
c74c5f6c 243}
244
245
246
247AliTPCcalibTime::~AliTPCcalibTime(){
248 //
2be25cc9 249 // Destructor
c74c5f6c 250 //
2be25cc9 251 delete fHistDeDxTgl;
252 delete fHistDeDx;
253 delete fHistVdrift;
254 for (Int_t i=0;i<3;i++){
255 delete fHistVdriftLaserA[i];
256 delete fHistVdriftLaserC[i];
257 }
258 fMapDz->SetOwner();
259 fMapDz->Delete();
260 delete fMapDz;
c74c5f6c 261}
2be25cc9 262
263
da6c0bc9 264void AliTPCcalibTime::ResetCurrent(){
265 //
266 // reset current values
267 //
268 fDz=0; // current delta z
269 fdEdx=0; // current dEdx
270 fdEdxRatio=0; // current dEdx ratio
271 fTl=0; // current tan(lambda)
272
273}
274
2be25cc9 275Bool_t AliTPCcalibTime::IsLaser(AliESDEvent *event){
276 return ((event->GetFiredTriggerClasses()).Contains("0LSR")==1);
277}
c74c5f6c 278
2be25cc9 279void AliTPCcalibTime::Process(AliESDEvent *event){
280 if(!event) return;
281 if (event->GetNumberOfTracks()<2) return;
282 ResetCurrent();
283
284// if(IsLaser(event))
285 ProcessLaser (event);
286// else
287 ProcessCosmic(event);
288}
289
290void AliTPCcalibTime::ProcessLaser(AliESDEvent *event){
c74c5f6c 291 //
2be25cc9 292 // Fit drift velocity using laser
293 //
294 // 0. cuts
295 const Int_t kMinTracks = 20; // minimal number of laser tracks
296 const Int_t kMinTracksSide = 10; // minimal number of tracks per side
297 const Float_t kMaxRMS = 0.1; // maximal RMS of tracks
298 const Float_t kMaxDeltaZ = 3.; // maximal deltaZ A-C side
299 const Float_t kMaxDeltaV = 0.01; // maximal deltaV A-C side
300 const Float_t kMaxDeltaY = 2.; // maximal deltaY A-C side
301 /*
302 TCut cutRMS("sqrt(laserA.fElements[4])<0.1&&sqrt(laserC.fElements[4])<0.1");
303 TCut cutZ("abs(laserA.fElements[0]-laserC.fElements[0])<3");
304 TCut cutV("abs(laserA.fElements[1]-laserC.fElements[1])<0.01");
305 TCut cutY("abs(laserA.fElements[2]-laserC.fElements[2])<2");
306 TCut cutAll = cutRMS+cutZ+cutV+cutY;
307 */
308 if (event->GetNumberOfTracks()<kMinTracks) return;
c74c5f6c 309 //
2be25cc9 310 if(!fLaser) fLaser = new AliTPCcalibLaser("laserTPC","laserTPC",kFALSE);
311 fLaser->Process(event);
312 if (fLaser->GetNtracks()<kMinTracks) return; // small amount of tracks cut
313 if (fLaser->fFitAside->GetNrows()==0) return; // no fit A side
314 if (fLaser->fFitCside->GetNrows()==0) return; // no fit C side
c74c5f6c 315 //
2be25cc9 316 // debug streamer - activate stream level
317 // Use it for tuning of the cuts
da6c0bc9 318 //
2be25cc9 319 if (fStreamLevel>0){
320 printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
321
da6c0bc9 322 TTreeSRedirector *cstream = GetDebugStreamer();
323 if (cstream){
324 TTimeStamp tstamp(fTime);
2be25cc9 325 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
da6c0bc9 326 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
327 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
328 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
329 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
330 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
31aa7c5c 331 TVectorD vecGoofie(20);
332 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
2be25cc9 333 if (goofieArray){
31aa7c5c 334 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
335 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
336 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
337 }
da6c0bc9 338 }
2be25cc9 339 (*cstream)<<"laserInfo"<<
da6c0bc9 340 "run="<<fRun<< // run number
341 "event="<<fEvent<< // event number
342 "time="<<fTime<< // time stamp of event
343 "trigger="<<fTrigger<< // trigger
344 "mag="<<fMagF<< // magnetic field
345 // Environment values
346 "press0="<<valuePressure0<<
347 "press1="<<valuePressure1<<
348 "pt0="<<ptrelative0<<
349 "pt1="<<ptrelative1<<
350 "temp0="<<temp0<<
351 "temp1="<<temp1<<
31aa7c5c 352 "vecGoofie.=<<"<<&vecGoofie<<
da6c0bc9 353 //laser
2be25cc9 354 "laserA.="<<fLaser->fFitAside<<
355 "laserC.="<<fLaser->fFitCside<<
356 "laserAC.="<<fLaser->fFitACside<<
357 "trigger="<<event->GetFiredTriggerClasses()<<
da6c0bc9 358 "\n";
359 }
360 }
2be25cc9 361 //
362 // Apply custs
363 //
364 if ((*fLaser->fFitAside)[3] <kMinTracksSide) return; // enough tracks A side
365 if ((*fLaser->fFitCside)[3]<kMinTracksSide) return; // enough tracks C side
366 //
367 if (TMath::Abs((*fLaser->fFitAside)[0]-(*fLaser->fFitCside)[0])>kMaxDeltaZ) return;
368 if (TMath::Abs((*fLaser->fFitAside)[2]-(*fLaser->fFitCside)[2])>kMaxDeltaY) return;
369 if (TMath::Abs((*fLaser->fFitAside)[1]-(*fLaser->fFitCside)[1])>kMaxDeltaV) return;
370 if (TMath::Sqrt(TMath::Abs((*fLaser->fFitAside)[4]))>kMaxRMS) return;
371 if (TMath::Sqrt(TMath::Abs((*fLaser->fFitCside)[4]))>kMaxRMS) return;
372 //
373 // fill histos
374 //
375 TVectorD vdriftA(5), vdriftC(5),vdriftAC(5);
376 vdriftA=*(fLaser->fFitAside);
377 vdriftC=*(fLaser->fFitCside);
378 vdriftAC=*(fLaser->fFitACside);
379 Int_t npointsA=0, npointsC=0;
380 Float_t chi2A=0, chi2C=0;
381 npointsA= TMath::Nint(vdriftA[3]);
382 chi2A= vdriftA[4];
383 npointsC= TMath::Nint(vdriftC[3]);
384 chi2C= vdriftC[4];
385 //
386 //
387
388 //
389 TTimeStamp tstamp(fTime);
390 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
391 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
392 //
393 Double_t vecVdriftLaserA[4]={fTime,(ptrelative0+ptrelative1)/2.0,1./((*(fLaser->fFitAside))[1])-1,event->GetRunNumber()};
394 Double_t vecVdriftLaserC[4]={fTime,(ptrelative0+ptrelative1)/2.0,1./((*(fLaser->fFitCside))[1])-1,event->GetRunNumber()};
395
396 for (Int_t i=0;i<3;i++){
397 if (i==0){ //z0 shift
398 vecVdriftLaserA[3]=(*(fLaser->fFitAside))[0]/250.;
399 vecVdriftLaserA[3]=(*(fLaser->fFitCside))[0]/250.;
400 }
401 if (i==1){ //vdrel shift
402 vecVdriftLaserA[3]=1./(*(fLaser->fFitAside))[1]-1.;
403 vecVdriftLaserA[3]=1./(*(fLaser->fFitCside))[1]-1.;
404 }
405 if (i==2){ //gy shift - full gy - full drift
406 vecVdriftLaserA[3]=(*(fLaser->fFitAside))[2]/250.;
407 vecVdriftLaserA[3]=(*(fLaser->fFitCside))[2]/250.;
408 }
409 fHistVdriftLaserA[i]->Fill(vecVdriftLaserA);
410 fHistVdriftLaserC[i]->Fill(vecVdriftLaserC);
411 }
412 //
413 //
414 //
c74c5f6c 415
416
2be25cc9 417}
c74c5f6c 418
2be25cc9 419void AliTPCcalibTime::ProcessCosmic(AliESDEvent *event){
c74c5f6c 420
421 if (!event) {
422 Printf("ERROR: ESD not available");
423 return;
424 }
425 if (event->GetTimeStamp() == 0 ) {
426 Printf("no time stamp!");
427 return;
428 }
429
c74c5f6c 430
431 UInt_t time = event->GetTimeStamp();
432
2be25cc9 433 //fd
c74c5f6c 434 // Find cosmic pairs
435 //
436 // Track0 is choosen in upper TPC part
437 // Track1 is choosen in lower TPC part
438 //
439 Int_t ntracks=event->GetNumberOfTracks();
440 if (ntracks==0) return;
441 if (ntracks > 10) return; // temporary debug to remove LASER events
442
443
444 if (GetDebugLevel()>1) printf("Hallo world: Im here\n");
445 AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
446
447 TObjArray tpcSeeds(ntracks);
448 Double_t vtxx[3]={0,0,0};
449 Double_t svtxx[3]={0.000001,0.000001,100.};
450 AliESDVertex vtx(vtxx,svtxx);
451 //
452 // track loop
453 //
454 for (Int_t i=0;i<ntracks;++i) {
455 AliESDtrack *track = event->GetTrack(i);
456
457 const AliExternalTrackParam * trackIn = track->GetInnerParam();
458 const AliExternalTrackParam * trackOut = track->GetOuterParam();
459 if (!trackIn) continue;
460 if (!trackOut) continue;
461
462 AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
463 TObject *calibObject;
464 AliTPCseed *seed = 0;
465 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
466 if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
467 }
468 if (seed) {
469 tpcSeeds.AddAt(seed,i);
470 if (track->GetTPCNcls() > 50) {
471 Double_t meanP = 0.5*(trackIn->GetP() + trackOut->GetP());
472 Double_t TPCsignal = seed->CookdEdxNorm(0.0,0.6,1,0,159,0x0,kTRUE,kTRUE);
473 Double_t vecDeDx[2] = {time, TPCsignal};
da6c0bc9 474 if (meanP > 12) {
475 fdEdx = TPCsignal;
476 fHistDeDx->Fill(vecDeDx);
477 }
c74c5f6c 478 }
479 }
480 }
481
482 if (ntracks<2) return;
483 //
484 // Find pairs
485 //
486 for (Int_t i=0;i<ntracks;++i) {
487 AliESDtrack *track0 = event->GetTrack(i);
488 // track0 - choosen upper part
489 if (!track0) continue;
490 if (!track0->GetOuterParam()) continue;
491 if (track0->GetOuterParam()->GetAlpha()<0) continue;
492 Double_t d1[3];
493 track0->GetDirection(d1);
494 for (Int_t j=0;j<ntracks;++j) {
495 if (i==j) continue;
496 AliESDtrack *track1 = event->GetTrack(j);
497 //track 1 lower part
498 if (!track1) continue;
499 if (!track1->GetOuterParam()) continue;
500 if (track1->GetOuterParam()->GetAlpha()>0) continue;
501 //
502 Double_t d2[3];
503 track1->GetDirection(d2);
504
505 AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
506 AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
507 if (! seed0) continue;
508 if (! seed1) continue;
509 Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
57dc06f2 510 Float_t dist0 = track0->GetLinearD(0,0);
511 Float_t dist1 = track1->GetLinearD(0,0);
c74c5f6c 512 //
513 // conservative cuts - convergence to be guarantied
514 // applying before track propagation
57dc06f2 515 if (TMath::Abs(dist0+dist1)>fCutMaxD) continue; // distance to the 0,0
c74c5f6c 516 if (dir>fCutMinDir) continue; // direction vector product
517 Float_t bz = AliTracker::GetBz();
518 Float_t dvertex0[2]; //distance to 0,0
519 Float_t dvertex1[2]; //distance to 0,0
520 track0->GetDZ(0,0,0,bz,dvertex0);
521 track1->GetDZ(0,0,0,bz,dvertex1);
522 if (TMath::Abs(dvertex0[1])>250) continue;
523 if (TMath::Abs(dvertex1[1])>250) continue;
524 //
525 //
526 //
57dc06f2 527 Float_t dmax = TMath::Max(TMath::Abs(dist0),TMath::Abs(dist1));
c74c5f6c 528 AliExternalTrackParam param0(*track0);
529 AliExternalTrackParam param1(*track1);
530 //
531 // Propagate using Magnetic field and correct fo material budget
532 //
533 AliTracker::PropagateTrackTo(&param0,dmax+1,0.0005,3,kTRUE);
534 AliTracker::PropagateTrackTo(&param1,dmax+1,0.0005,3,kTRUE);
535 //
536 // Propagate rest to the 0,0 DCA - z should be ignored
537 //
828383be 538 //Bool_t b0 = ;
539 param0.PropagateToDCA(&vtx,bz,1000);
540 //Bool_t b1 =
541 param1.PropagateToDCA(&vtx,bz,1000);
c74c5f6c 542 //
543 param0.GetDZ(0,0,0,bz,dvertex0);
544 param1.GetDZ(0,0,0,bz,dvertex1);
545 //
546 Double_t xyz0[3];//,pxyz0[3];
547 Double_t xyz1[3];//,pxyz1[3];
548 param0.GetXYZ(xyz0);
549 param1.GetXYZ(xyz1);
550 Bool_t isPair = IsPair(&param0,&param1);
551
552 Double_t z0 = track0->GetOuterParam()->GetZ();
553 Double_t z1 = track1->GetOuterParam()->GetZ();
554
555 Double_t z0inner = track0->GetInnerParam()->GetZ();
556 Double_t z1inner = track1->GetInnerParam()->GetZ();
557
558 if (isPair && z0>0 && z0inner>0 && z1<0 && z1inner<0) {
559 Double_t vecVdrift[2] = {time, param0.GetZ() - param1.GetZ()};
560
da6c0bc9 561 if (track0->GetTPCNcls() > 80) {
da6c0bc9 562 fDz = param0.GetZ() - param1.GetZ();
2be25cc9 563 TTimeStamp tstamp(fTime);
564 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
565 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
566 THnSparse* curHist =(THnSparseF*)(fMapDz->GetValue(event->GetFiredTriggerClasses()));
567 if(!curHist){
568 curHist=new THnSparseF(event->GetFiredTriggerClasses(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
569 fMapDz->Add(new TObjString(event->GetFiredTriggerClasses()),curHist);
570 }
571 Double_t vecVdrift2[4]={fTime,(ptrelative0+ptrelative1)/2.0,fDz/500.0,event->GetRunNumber()};
572 curHist->Fill(vecVdrift2);
573
574 fHistVdrift->Fill(vecVdrift);
da6c0bc9 575 }
c74c5f6c 576 }
577 if (isPair && z0 > 0 && z1 > 0) {
578 if (track1->GetTPCNcls()> 110 && track0->GetTPCNcls()> 110 && seed1->CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE) > 0) {
da6c0bc9 579 Float_t dedxratio = seed0->CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE)/seed1->CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE);
c74c5f6c 580 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)};
581 fHistDeDxTgl->Fill(vecDeDxTgl);
da6c0bc9 582 fdEdxRatio = dedxratio;
583 fTl = track0->GetTgl() ;
c74c5f6c 584 }
585 }
586
587 } // end 2nd order loop
588 } // end 1st order loop
589
2be25cc9 590 if (fStreamLevel>0){
591 TTreeSRedirector *cstream = GetDebugStreamer();
592 if (cstream){
593 TTimeStamp tstamp(fTime);
594 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
595 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
596 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
597 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
598 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
599 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
600 TVectorD vecGoofie(20);
601 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
602 if (goofieArray){
603 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
604 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
605 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
606 }
607 }
608 (*cstream)<<"timeInfo"<<
609 "run="<<fRun<< // run number
610 "event="<<fEvent<< // event number
611 "time="<<fTime<< // time stamp of event
612 "trigger="<<fTrigger<< // trigger
613 "mag="<<fMagF<< // magnetic field
614 // Environment values
615 "press0="<<valuePressure0<<
616 "press1="<<valuePressure1<<
617 "pt0="<<ptrelative0<<
618 "pt1="<<ptrelative1<<
619 "temp0="<<temp0<<
620 "temp1="<<temp1<<
621 "vecGoofie.=<<"<<&vecGoofie<<
622 //
623 // accumulated values
624 //
625 "fDz="<<fDz<< //! current delta z
626 "fdEdx="<<fdEdx<< //! current dEdx
627 "fdEdxRatio="<<fdEdxRatio<< //! current dEdx ratio
628 "fTl="<<fTl<< //! current tan(lambda)
629 "trigger="<<event->GetFiredTriggerClasses()<<
630 "\n";
631 }
632 }
633 printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
c74c5f6c 634}
635
c74c5f6c 636void AliTPCcalibTime::Analyze() {
828383be 637 //
638 //
639 //
c74c5f6c 640 TH2D * hVdrift = GetHistVdrift()->Projection(1,0);
0ffd2ae1 641 hVdrift->Draw();
c74c5f6c 642}
643
c74c5f6c 644Long64_t AliTPCcalibTime::Merge(TCollection *li) {
645
646 TIterator* iter = li->MakeIterator();
647 AliTPCcalibTime* cal = 0;
648
649 while ((cal = (AliTPCcalibTime*)iter->Next())) {
650 if (!cal->InheritsFrom(AliTPCcalibTime::Class())) {
651 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
652 return -1;
653 }
654
655 // add histograms here...
656 fHistDeDxTgl->Add(cal->GetHistDeDxVsTgl());
657 fHistVdrift->Add(cal->GetHistVdrift());
658 fHistDeDx->Add(cal->GetHistDeDx());
659
2be25cc9 660 for (Int_t imeas=0; imeas<3; imeas++){
661 if (cal->GetHistVdriftLaserA(imeas) && cal->GetHistVdriftLaserA(imeas)){
662 fHistVdriftLaserA[imeas]->Add(cal->GetHistVdriftLaserA(imeas));
663 fHistVdriftLaserC[imeas]->Add(cal->GetHistVdriftLaserC(imeas));
664 }
665 }
666 TMap * addMap=cal->GetMapDz();
667 if(!addMap) return 0;
668 TIterator* iterator = addMap->MakeIterator();
669 iterator->Reset();
670 TPair* addPair=0;
671 while((addPair=(TPair *)(addMap->FindObject(iterator->Next())))){
672 THnSparse* addHist=dynamic_cast<THnSparseF*>(addPair->Value());
673 if (!addHist) continue;
674 addHist->Print();
675 THnSparse* localHist=dynamic_cast<THnSparseF*>(fMapDz->GetValue(addHist->GetName()));
676 if(!localHist){
677 localHist=new THnSparseF(addHist->GetName(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
678 fMapDz->Add(new TObjString(addHist->GetName()),localHist);
679 }
680 localHist->Add(addHist);
681 }
c74c5f6c 682 }
c74c5f6c 683 return 0;
c74c5f6c 684}
685
c74c5f6c 686Bool_t AliTPCcalibTime::IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1){
687 //
688 //
689 /*
690 // 0. Same direction - OPOSITE - cutDir +cutT
691 TCut cutDir("cutDir","dir<-0.99")
692 // 1.
693 TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03")
694 //
695 // 2. The same rphi
696 TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5")
697 //
698 //
699 //
700 TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
701 // 1/Pt diff cut
702 */
703 const Double_t *p0 = tr0->GetParameter();
704 const Double_t *p1 = tr1->GetParameter();
705 if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
706 if (TMath::Abs(p0[0]+p1[0])>fCutMaxD) return kFALSE;
da6c0bc9 707 if (TMath::Abs(p0[1]-p1[1])>fCutMaxDz) return kFALSE;
c74c5f6c 708 Double_t d0[3], d1[3];
709 tr0->GetDirection(d0);
710 tr1->GetDirection(d1);
711 if (d0[0]*d1[0] + d0[1]*d1[1] + d0[2]*d1[2] >fCutMinDir) return kFALSE;
712 //
713 return kTRUE;
714}
31aa7c5c 715
716