Adding visualization of the tree aliases
[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/*
74235403 17Comments to be written here:
c74c5f6c 18
74235403 191. What do we calibrate.
20
21 Time dependence of gain and drift velocity in order to account for changes in: temperature, pressure, gas composition.
c74c5f6c 22
23 AliTPCcalibTime *calibTime = new AliTPCcalibTime("cosmicTime","cosmicTime",0, 1213.9e+06, 1213.96e+06, 0.04e+04, 0.04e+04);
24
74235403 252. How to interpret results
26
273. Simple example
c74c5f6c 28
74235403 29 a) determine the required time range:
c74c5f6c 30
74235403 31 AliXRDPROOFtoolkit tool;
32 TChain * chain = tool.MakeChain("pass2.txt","esdTree",0,6000);
33 chain->Draw("GetTimeStamp()")
c74c5f6c 34
74235403 35 b) analyse calibration object on Proof in calibration train
c74c5f6c 36
74235403 37 AliTPCcalibTime *calibTime = new AliTPCcalibTime("cosmicTime","cosmicTime", StartTimeStamp, EndTimeStamp, IntegrationTimeVdrift);
c74c5f6c 38
74235403 39 c) plot results
40 .x ~/NimStyle.C
2be25cc9 41 gSystem->Load("libANALYSIS");
42 gSystem->Load("libTPCcalib");
c74c5f6c 43
dde68d36 44 TFile f("CalibObjectsTrain1.root");
45 AliTPCcalibTime *calib = (AliTPCcalibTime *)f->Get("calibTime");
46 calib->GetHistoDrift("all")->Projection(2,0)->Draw()
47 calib->GetFitDrift("all")->Draw("lp")
da6c0bc9 48
74235403 494. Analysis using debug streamers.
da6c0bc9 50
74235403 51 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
52 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
53 AliXRDPROOFtoolkit tool;
817766d5 54 AliXRDPROOFtoolkit::FilterList("timeitstpc.txt","* itstpc",1)
a66da2d4 55 AliXRDPROOFtoolkit::FilterList("timetoftpc.txt","* pointMatch",1)
817766d5 56 AliXRDPROOFtoolkit::FilterList("time.txt","* trackInfo",1)
57 AliXRDPROOFtoolkit::FilterList("timelaser.txt","* laserInfo",1)
58
59 TChain * chainTPCITS = tool.MakeChainRandom("timeitstpc.txt.Good","itstpc",0,10000);
a66da2d4 60 TChain * chainTPCTOF = tool.MakeChainRandom("timetoftpc.txt","pointMatch",0,500);
817766d5 61 TChain * chainTime = tool.MakeChainRandom("time.txt.Good","trackInfo",0,10000);
62 TChain * chainLaser = tool.MakeChainRandom("timelaser.txt.Good","laserInfo",0,10000);
74235403 63 chainTime->Lookup();
dde68d36 64 chainLaser->Lookup();
c74c5f6c 65*/
66
c74c5f6c 67#include "Riostream.h"
68#include "TChain.h"
69#include "TTree.h"
70#include "TH1F.h"
71#include "TH2F.h"
72#include "TH3F.h"
73#include "THnSparse.h"
74#include "TList.h"
75#include "TMath.h"
76#include "TCanvas.h"
77#include "TFile.h"
78#include "TF1.h"
79#include "TVectorD.h"
80#include "TProfile.h"
81#include "TGraphErrors.h"
82#include "TCanvas.h"
c74c5f6c 83#include "AliTPCclusterMI.h"
84#include "AliTPCseed.h"
85#include "AliESDVertex.h"
86#include "AliESDEvent.h"
87#include "AliESDfriend.h"
88#include "AliESDInputHandler.h"
89#include "AliAnalysisManager.h"
90
91#include "AliTracker.h"
f7a1cc68 92#include "AliMagF.h"
c74c5f6c 93#include "AliTPCCalROC.h"
94
95#include "AliLog.h"
96
97#include "AliTPCcalibTime.h"
817766d5 98#include "AliRelAlignerKalman.h"
c74c5f6c 99
100#include "TTreeStream.h"
101#include "AliTPCTracklet.h"
da6c0bc9 102#include "TTimeStamp.h"
103#include "AliTPCcalibDB.h"
104#include "AliTPCcalibLaser.h"
31aa7c5c 105#include "AliDCSSensorArray.h"
106#include "AliDCSSensor.h"
c74c5f6c 107
817766d5 108#include "TDatabasePDG.h"
109#include "AliTrackPointArray.h"
110
c74c5f6c 111ClassImp(AliTPCcalibTime)
112
113
114AliTPCcalibTime::AliTPCcalibTime()
da6c0bc9 115 :AliTPCcalibBase(),
da6c0bc9 116 fLaser(0), // pointer to laser calibration
117 fDz(0), // current delta z
d3ce44cb 118 fCutMaxD(3), // maximal distance in rfi ditection
119 fCutMaxDz(25), // maximal distance in rfi ditection
c74c5f6c 120 fCutTheta(0.03), // maximal distan theta
2be25cc9 121 fCutMinDir(-0.99), // direction vector products
74235403 122 fCutTracks(10),
dde68d36 123 fArrayDz(0), //NEW! Tmap of V drifts for different triggers
817766d5 124 fAlignITSTPC(0), //alignemnt array ITS TPC match
125 fAlignTRDTPC(0), //alignemnt array TRD TPC match
126 fAlignTOFTPC(0), //alignemnt array TOF TPC match
2be25cc9 127 fTimeBins(0),
128 fTimeStart(0),
129 fTimeEnd(0),
130 fPtBins(0),
131 fPtStart(0),
132 fPtEnd(0),
133 fVdriftBins(0),
134 fVdriftStart(0),
135 fVdriftEnd(0),
136 fRunBins(0),
137 fRunStart(0),
138 fRunEnd(0)
139// fBinsVdrift(fTimeBins,fPtBins,fVdriftBins),
140// fXminVdrift(fTimeStart,fPtStart,fVdriftStart),
141// fXmaxVdrift(fTimeEnd,fPtEnd,fVdriftEnd)
c74c5f6c 142{
143 AliInfo("Default Constructor");
2be25cc9 144 for (Int_t i=0;i<3;i++) {
145 fHistVdriftLaserA[i]=0;
146 fHistVdriftLaserC[i]=0;
147 }
d3ce44cb 148 for (Int_t i=0;i<10;i++) {
74235403 149 fCosmiMatchingHisto[i]=0;
150 }
c74c5f6c 151}
152
74235403 153AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeVdrift)
c74c5f6c 154 :AliTPCcalibBase(),
d3ce44cb 155 fLaser(0), // pointer to laser calibration
156 fDz(0), // current delta z
157 fCutMaxD(5*0.5356), // maximal distance in rfi ditection
dde68d36 158 fCutMaxDz(40), // maximal distance in rfi ditection
d3ce44cb 159 fCutTheta(5*0.004644),// maximal distan theta
160 fCutMinDir(-0.99), // direction vector products
74235403 161 fCutTracks(10),
dde68d36 162 fArrayDz(0), //Tmap of V drifts for different triggers
817766d5 163 fAlignITSTPC(0), //alignemnt array ITS TPC match
164 fAlignTRDTPC(0), //alignemnt array TRD TPC match
165 fAlignTOFTPC(0), //alignemnt array TOF TPC match
2be25cc9 166 fTimeBins(0),
167 fTimeStart(0),
168 fTimeEnd(0),
169 fPtBins(0),
170 fPtStart(0),
171 fPtEnd(0),
172 fVdriftBins(0),
173 fVdriftStart(0),
174 fVdriftEnd(0),
175 fRunBins(0),
176 fRunStart(0),
177 fRunEnd(0)
c74c5f6c 178{
c74c5f6c 179 SetName(name);
180 SetTitle(title);
2be25cc9 181 for (Int_t i=0;i<3;i++) {
182 fHistVdriftLaserA[i]=0;
183 fHistVdriftLaserC[i]=0;
184 }
c74c5f6c 185
186 AliInfo("Non Default Constructor");
74235403 187 fTimeBins =(EndTime-StartTime)/deltaIntegrationTimeVdrift;
2be25cc9 188 fTimeStart =StartTime; //(((TObjString*)(mapGRP->GetValue("fAliceStartTime")))->GetString()).Atoi();
189 fTimeEnd =EndTime; //(((TObjString*)(mapGRP->GetValue("fAliceStopTime")))->GetString()).Atoi();
dde68d36 190 fPtBins = 400;
2be25cc9 191 fPtStart = -0.04;
192 fPtEnd = 0.04;
dde68d36 193 fVdriftBins = 500;
194 fVdriftStart= -0.1;
195 fVdriftEnd = 0.1;
196 fRunBins = 100001;
197 fRunStart = -1.5;
d3ce44cb 198 fRunEnd = 99999.5;
2be25cc9 199
200 Int_t binsVdriftLaser[4] = {fTimeBins , fPtBins , fVdriftBins*20, fRunBins };
201 Double_t xminVdriftLaser[4] = {fTimeStart, fPtStart, fVdriftStart , fRunStart};
202 Double_t xmaxVdriftLaser[4] = {fTimeEnd , fPtEnd , fVdriftEnd , fRunEnd };
dde68d36 203 TString axisTitle[4]={
204 "T",
205 "#delta_{P/T}",
206 "value",
207 "run"
208 };
209 TString histoName[3]={
210 "Loffset",
211 "Lcorr",
212 "Lgy"
213 };
2be25cc9 214
dde68d36 215
2be25cc9 216 for (Int_t i=0;i<3;i++) {
217 fHistVdriftLaserA[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
218 fHistVdriftLaserC[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
dde68d36 219 fHistVdriftLaserA[i]->SetName(histoName[i]);
220 fHistVdriftLaserC[i]->SetName(histoName[i]);
221 for (Int_t iaxis=0; iaxis<4;iaxis++){
222 fHistVdriftLaserA[i]->GetAxis(iaxis)->SetName(axisTitle[iaxis]);
223 fHistVdriftLaserC[i]->GetAxis(iaxis)->SetName(axisTitle[iaxis]);
224 }
2be25cc9 225 }
226 fBinsVdrift[0] = fTimeBins;
227 fBinsVdrift[1] = fPtBins;
228 fBinsVdrift[2] = fVdriftBins;
229 fBinsVdrift[3] = fRunBins;
230 fXminVdrift[0] = fTimeStart;
231 fXminVdrift[1] = fPtStart;
232 fXminVdrift[2] = fVdriftStart;
233 fXminVdrift[3] = fRunStart;
234 fXmaxVdrift[0] = fTimeEnd;
235 fXmaxVdrift[1] = fPtEnd;
236 fXmaxVdrift[2] = fVdriftEnd;
237 fXmaxVdrift[3] = fRunEnd;
238
dde68d36 239 fArrayDz=new TObjArray();
817766d5 240 fAlignITSTPC = new TObjArray; //alignemnt array ITS TPC match
241 fAlignTRDTPC = new TObjArray; //alignemnt array ITS TPC match
242 fAlignTOFTPC = new TObjArray; //alignemnt array ITS TPC match
243
abae1b29 244 // fArrayDz->AddLast(fHistVdriftLaserA[0]);
245// fArrayDz->AddLast(fHistVdriftLaserA[1]);
246// fArrayDz->AddLast(fHistVdriftLaserA[2]);
247// fArrayDz->AddLast(fHistVdriftLaserC[0]);
248// fArrayDz->AddLast(fHistVdriftLaserC[1]);
249// fArrayDz->AddLast(fHistVdriftLaserC[2]);
c74c5f6c 250
d3ce44cb 251 fCosmiMatchingHisto[0]=new TH1F("Cosmics matching","p0-all" ,100,-10*0.5356 ,10*0.5356 );
252 fCosmiMatchingHisto[1]=new TH1F("Cosmics matching","p1-all" ,100,-10*4.541 ,10*4.541 );
253 fCosmiMatchingHisto[2]=new TH1F("Cosmics matching","p2-all" ,100,-10*0.01134 ,10*0.01134 );
254 fCosmiMatchingHisto[3]=new TH1F("Cosmics matching","p3-all" ,100,-10*0.004644,10*0.004644);
255 fCosmiMatchingHisto[4]=new TH1F("Cosmics matching","p4-all" ,100,-10*0.03773 ,10*0.03773 );
256 fCosmiMatchingHisto[5]=new TH1F("Cosmics matching","p0-isPair",100,-10*0.5356 ,10*0.5356 );
257 fCosmiMatchingHisto[6]=new TH1F("Cosmics matching","p1-isPair",100,-10*4.541 ,10*4.541 );
258 fCosmiMatchingHisto[7]=new TH1F("Cosmics matching","p2-isPair",100,-10*0.01134 ,10*0.01134 );
259 fCosmiMatchingHisto[8]=new TH1F("Cosmics matching","p3-isPair",100,-10*0.004644,10*0.004644);
260 fCosmiMatchingHisto[9]=new TH1F("Cosmics matching","p4-isPair",100,-10*0.03773 ,10*0.03773 );
261// Char_t nameHisto[3]={'p','0','\n'};
262// for (Int_t i=0;i<10;i++){
263// fCosmiMatchingHisto[i]=new TH1F("Cosmics matching",nameHisto,8192,0,0);
264// nameHisto[1]++;
265// if(i==4) nameHisto[1]='0';
266// }
74235403 267}
c74c5f6c 268
269AliTPCcalibTime::~AliTPCcalibTime(){
270 //
2be25cc9 271 // Destructor
c74c5f6c 272 //
d3ce44cb 273 for(Int_t i=0;i<3;i++){
274 if(fHistVdriftLaserA[i]){
275 delete fHistVdriftLaserA[i];
276 fHistVdriftLaserA[i]=NULL;
277 }
278 if(fHistVdriftLaserC[i]){
279 delete fHistVdriftLaserC[i];
280 fHistVdriftLaserC[i]=NULL;
281 }
282 }
dde68d36 283 if(fArrayDz){
284 fArrayDz->SetOwner();
285 fArrayDz->Delete();
286 delete fArrayDz;
287 fArrayDz=NULL;
2be25cc9 288 }
d3ce44cb 289 for(Int_t i=0;i<5;i++){
290 if(fCosmiMatchingHisto[i]){
291 delete fCosmiMatchingHisto[i];
292 fCosmiMatchingHisto[i]=NULL;
293 }
74235403 294 }
817766d5 295 fAlignITSTPC->Delete();
296 fAlignTRDTPC->Delete();
297 fAlignTOFTPC->Delete();
c74c5f6c 298}
2be25cc9 299
dde68d36 300Bool_t AliTPCcalibTime::IsLaser(AliESDEvent */*event*/){
301 return kTRUE; //More accurate creteria to be added
da6c0bc9 302}
dde68d36 303Bool_t AliTPCcalibTime::IsCosmics(AliESDEvent */*event*/){
304 return kTRUE; //More accurate creteria to be added
305}
306Bool_t AliTPCcalibTime::IsBeam(AliESDEvent */*event*/){
307 return kTRUE; //More accurate creteria to be added
308}
309void AliTPCcalibTime::ResetCurrent(){
310 fDz=0; //Reset current dz
2be25cc9 311}
2be25cc9 312void AliTPCcalibTime::Process(AliESDEvent *event){
313 if(!event) return;
314 if (event->GetNumberOfTracks()<2) return;
315 ResetCurrent();
dde68d36 316 if(IsLaser (event)) ProcessLaser (event);
317 if(IsCosmics(event)) ProcessCosmic(event);
318 if(IsBeam (event)) ProcessBeam (event);
2be25cc9 319}
320
321void AliTPCcalibTime::ProcessLaser(AliESDEvent *event){
c74c5f6c 322 //
2be25cc9 323 // Fit drift velocity using laser
324 //
325 // 0. cuts
dde68d36 326 const Int_t kMinTracks = 40; // minimal number of laser tracks
327 const Int_t kMinTracksSide = 20; // minimal number of tracks per side
328 const Float_t kMaxDeltaZ = 30.; // maximal trigger delay
329 const Float_t kMaxDeltaV = 0.05; // maximal deltaV
2be25cc9 330 const Float_t kMaxRMS = 0.1; // maximal RMS of tracks
dde68d36 331 //
2be25cc9 332 /*
333 TCut cutRMS("sqrt(laserA.fElements[4])<0.1&&sqrt(laserC.fElements[4])<0.1");
334 TCut cutZ("abs(laserA.fElements[0]-laserC.fElements[0])<3");
335 TCut cutV("abs(laserA.fElements[1]-laserC.fElements[1])<0.01");
336 TCut cutY("abs(laserA.fElements[2]-laserC.fElements[2])<2");
337 TCut cutAll = cutRMS+cutZ+cutV+cutY;
338 */
339 if (event->GetNumberOfTracks()<kMinTracks) return;
c74c5f6c 340 //
2be25cc9 341 if(!fLaser) fLaser = new AliTPCcalibLaser("laserTPC","laserTPC",kFALSE);
342 fLaser->Process(event);
343 if (fLaser->GetNtracks()<kMinTracks) return; // small amount of tracks cut
dde68d36 344 if (fLaser->fFitAside->GetNrows()==0 && fLaser->fFitCside->GetNrows()==0) return; // no fit neither a or C side
c74c5f6c 345 //
2be25cc9 346 // debug streamer - activate stream level
347 // Use it for tuning of the cuts
da6c0bc9 348 //
dde68d36 349 // cuts to be applied
350 //
351 Int_t isReject[2]={0,0};
352 //
353 // not enough tracks
354 if (TMath::Abs((*fLaser->fFitAside)[3]) < kMinTracksSide) isReject[0]|=1;
355 if (TMath::Abs((*fLaser->fFitCside)[3]) < kMinTracksSide) isReject[1]|=1;
356 // unreasonable z offset
357 if (TMath::Abs((*fLaser->fFitAside)[0])>kMaxDeltaZ) isReject[0]|=2;
358 if (TMath::Abs((*fLaser->fFitCside)[0])>kMaxDeltaZ) isReject[1]|=2;
359 // unreasonable drift velocity
360 if (TMath::Abs((*fLaser->fFitAside)[1]-1)>kMaxDeltaV) isReject[0]|=4;
361 if (TMath::Abs((*fLaser->fFitCside)[1]-1)>kMaxDeltaV) isReject[1]|=4;
362 // big chi2
363 if (TMath::Sqrt(TMath::Abs((*fLaser->fFitAside)[4]))>kMaxRMS ) isReject[0]|=8;
364 if (TMath::Sqrt(TMath::Abs((*fLaser->fFitCside)[4]))>kMaxRMS ) isReject[1]|=8;
365
366
367
368
2be25cc9 369 if (fStreamLevel>0){
370 printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
371
da6c0bc9 372 TTreeSRedirector *cstream = GetDebugStreamer();
373 if (cstream){
374 TTimeStamp tstamp(fTime);
2be25cc9 375 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
da6c0bc9 376 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
377 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
378 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
379 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
380 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
31aa7c5c 381 TVectorD vecGoofie(20);
382 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
2be25cc9 383 if (goofieArray){
31aa7c5c 384 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
385 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
386 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
387 }
da6c0bc9 388 }
2be25cc9 389 (*cstream)<<"laserInfo"<<
da6c0bc9 390 "run="<<fRun<< // run number
391 "event="<<fEvent<< // event number
392 "time="<<fTime<< // time stamp of event
393 "trigger="<<fTrigger<< // trigger
394 "mag="<<fMagF<< // magnetic field
395 // Environment values
396 "press0="<<valuePressure0<<
397 "press1="<<valuePressure1<<
398 "pt0="<<ptrelative0<<
399 "pt1="<<ptrelative1<<
400 "temp0="<<temp0<<
401 "temp1="<<temp1<<
817766d5 402 "vecGoofie.="<<&vecGoofie<<
da6c0bc9 403 //laser
dde68d36 404 "rejectA="<<isReject[0]<<
405 "rejectC="<<isReject[1]<<
2be25cc9 406 "laserA.="<<fLaser->fFitAside<<
407 "laserC.="<<fLaser->fFitCside<<
408 "laserAC.="<<fLaser->fFitACside<<
409 "trigger="<<event->GetFiredTriggerClasses()<<
da6c0bc9 410 "\n";
411 }
412 }
2be25cc9 413 //
2be25cc9 414 // fill histos
415 //
416 TVectorD vdriftA(5), vdriftC(5),vdriftAC(5);
417 vdriftA=*(fLaser->fFitAside);
418 vdriftC=*(fLaser->fFitCside);
419 vdriftAC=*(fLaser->fFitACside);
420 Int_t npointsA=0, npointsC=0;
421 Float_t chi2A=0, chi2C=0;
422 npointsA= TMath::Nint(vdriftA[3]);
423 chi2A= vdriftA[4];
424 npointsC= TMath::Nint(vdriftC[3]);
425 chi2C= vdriftC[4];
2be25cc9 426
2be25cc9 427 TTimeStamp tstamp(fTime);
428 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
429 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
dde68d36 430 Double_t driftA=0, driftC=0;
431 if (vdriftA[1]>1.-kMaxDeltaV) driftA = 1./vdriftA[1]-1.;
432 if (vdriftC[1]>1.-kMaxDeltaV) driftC = 1./vdriftC[1]-1.;
433 //
434 Double_t vecDriftLaserA[4]={fTime,(ptrelative0+ptrelative1)/2.0,driftA,event->GetRunNumber()};
435 Double_t vecDriftLaserC[4]={fTime,(ptrelative0+ptrelative1)/2.0,driftC,event->GetRunNumber()};
436 // Double_t vecDrift[4] ={fTime,(ptrelative0+ptrelative1)/2.0,1./((*(fLaser->fFitACside))[1])-1,event->GetRunNumber()};
437
438 for (Int_t icalib=0;icalib<3;icalib++){
439 if (icalib==0){ //z0 shift
440 vecDriftLaserA[2]=vdriftA[0]/250.;
441 vecDriftLaserC[2]=vdriftC[0]/250.;
2be25cc9 442 }
dde68d36 443 if (icalib==1){ //vdrel shift
444 vecDriftLaserA[2]=driftA;
445 vecDriftLaserC[2]=driftC;
2be25cc9 446 }
dde68d36 447 if (icalib==2){ //gy shift - full gy - full drift
448 vecDriftLaserA[2]=vdriftA[2]/250.;
449 vecDriftLaserC[2]=vdriftC[2]/250.;
2be25cc9 450 }
abae1b29 451 if (isReject[0]==0) fHistVdriftLaserA[icalib]->Fill(vecDriftLaserA);
452 if (isReject[1]==0) fHistVdriftLaserC[icalib]->Fill(vecDriftLaserC);
2be25cc9 453 }
dde68d36 454
455// THnSparse* curHist=new THnSparseF("","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
456// TString shortName=curHist->ClassName();
457// shortName+="_MEAN_DRIFT_LASER_";
458// delete curHist;
459// curHist=NULL;
460// TString name="";
461
462// name=shortName;
463// name+=event->GetFiredTriggerClasses();
464// name.ToUpper();
465// curHist=(THnSparseF*)fArrayDz->FindObject(name);
466// if(!curHist){
467// curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
468// fArrayDz->AddLast(curHist);
469// }
470// curHist->Fill(vecDrift);
471
472// name=shortName;
473// name+="ALL";
474// name.ToUpper();
475// curHist=(THnSparseF*)fArrayDz->FindObject(name);
476// if(!curHist){
477// curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
478// fArrayDz->AddLast(curHist);
479// }
480// curHist->Fill(vecDrift);
2be25cc9 481}
c74c5f6c 482
2be25cc9 483void AliTPCcalibTime::ProcessCosmic(AliESDEvent *event){
c74c5f6c 484 if (!event) {
485 Printf("ERROR: ESD not available");
486 return;
487 }
488 if (event->GetTimeStamp() == 0 ) {
489 Printf("no time stamp!");
490 return;
491 }
74235403 492
2be25cc9 493 //fd
c74c5f6c 494 // Find cosmic pairs
495 //
496 // Track0 is choosen in upper TPC part
497 // Track1 is choosen in lower TPC part
498 //
499 Int_t ntracks=event->GetNumberOfTracks();
500 if (ntracks==0) return;
74235403 501 if (ntracks > fCutTracks) return;
502
c74c5f6c 503 if (GetDebugLevel()>1) printf("Hallo world: Im here\n");
504 AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
505
506 TObjArray tpcSeeds(ntracks);
507 Double_t vtxx[3]={0,0,0};
508 Double_t svtxx[3]={0.000001,0.000001,100.};
509 AliESDVertex vtx(vtxx,svtxx);
510 //
511 // track loop
512 //
513 for (Int_t i=0;i<ntracks;++i) {
74235403 514 AliESDtrack *track = event->GetTrack(i);
515
516 const AliExternalTrackParam * trackIn = track->GetInnerParam();
517 const AliExternalTrackParam * trackOut = track->GetOuterParam();
518 if (!trackIn) continue;
519 if (!trackOut) continue;
520
521 AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
817766d5 522 if (friendTrack) ProcessAlignITS(track,friendTrack);
523 if (friendTrack) ProcessAlignTOF(track,friendTrack);
74235403 524 TObject *calibObject;
525 AliTPCseed *seed = 0;
526 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
527 if (seed) tpcSeeds.AddAt(seed,i);
c74c5f6c 528 }
c74c5f6c 529 if (ntracks<2) return;
530 //
531 // Find pairs
532 //
d3ce44cb 533
c74c5f6c 534 for (Int_t i=0;i<ntracks;++i) {
74235403 535 AliESDtrack *track0 = event->GetTrack(i);
c74c5f6c 536 // track0 - choosen upper part
537 if (!track0) continue;
538 if (!track0->GetOuterParam()) continue;
539 if (track0->GetOuterParam()->GetAlpha()<0) continue;
540 Double_t d1[3];
541 track0->GetDirection(d1);
542 for (Int_t j=0;j<ntracks;++j) {
543 if (i==j) continue;
74235403 544 AliESDtrack *track1 = event->GetTrack(j);
545 //track 1 lower part
546 if (!track1) continue;
547 if (!track1->GetOuterParam()) continue;
548 if (track1->GetOuterParam()->GetAlpha()>0) continue;
549 //
550 Double_t d2[3];
551 track1->GetDirection(d2);
552
553 AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
554 AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
555 if (! seed0) continue;
556 if (! seed1) continue;
557 Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
558 Float_t dist0 = track0->GetLinearD(0,0);
559 Float_t dist1 = track1->GetLinearD(0,0);
560 //
561 // conservative cuts - convergence to be guarantied
562 // applying before track propagation
563 if (TMath::Abs(dist0+dist1)>fCutMaxD) continue; // distance to the 0,0
564 if (dir>fCutMinDir) continue; // direction vector product
565 Float_t bz = AliTracker::GetBz();
566 Float_t dvertex0[2]; //distance to 0,0
567 Float_t dvertex1[2]; //distance to 0,0
568 track0->GetDZ(0,0,0,bz,dvertex0);
569 track1->GetDZ(0,0,0,bz,dvertex1);
570 if (TMath::Abs(dvertex0[1])>250) continue;
571 if (TMath::Abs(dvertex1[1])>250) continue;
572 //
573 //
574 //
575 Float_t dmax = TMath::Max(TMath::Abs(dist0),TMath::Abs(dist1));
576 AliExternalTrackParam param0(*track0);
577 AliExternalTrackParam param1(*track1);
578 //
579 // Propagate using Magnetic field and correct fo material budget
580 //
581 AliTracker::PropagateTrackTo(&param0,dmax+1,0.0005,3,kTRUE);
582 AliTracker::PropagateTrackTo(&param1,dmax+1,0.0005,3,kTRUE);
583 //
584 // Propagate rest to the 0,0 DCA - z should be ignored
585 //
586 //Bool_t b0 = ;
587 param0.PropagateToDCA(&vtx,bz,1000);
588 //Bool_t b1 =
589 param1.PropagateToDCA(&vtx,bz,1000);
74235403 590 param0.GetDZ(0,0,0,bz,dvertex0);
591 param1.GetDZ(0,0,0,bz,dvertex1);
dde68d36 592 Double_t xyz0[3];
593 Double_t xyz1[3];
74235403 594 param0.GetXYZ(xyz0);
595 param1.GetXYZ(xyz1);
596 Bool_t isPair = IsPair(&param0,&param1);
d3ce44cb 597 Bool_t isCross = IsCross(track0, track1);
dde68d36 598 Bool_t isSame = IsSame(track0, track1);
d3ce44cb 599
dde68d36 600 THnSparse* hist=new THnSparseF("","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
601 TString shortName=hist->ClassName();
602 shortName+="_MEAN_VDRIFT_COSMICS_";
603 delete hist;
604 hist=NULL;
605
606 if(isSame || (isCross && isPair)){
74235403 607 if (track0->GetTPCNcls() > 80) {
608 fDz = param0.GetZ() - param1.GetZ();
d3ce44cb 609 if(track0->GetOuterParam()->GetZ()<0) fDz=-fDz;
74235403 610 TTimeStamp tstamp(fTime);
611 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
612 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
dde68d36 613 Double_t vecDrift[4]={fTime,(ptrelative0+ptrelative1)/2.0,fDz/500.0,event->GetRunNumber()};
614 THnSparse* curHist=NULL;
615 TString name="";
616
617 name=shortName;
618 name+=event->GetFiredTriggerClasses();
619 name.ToUpper();
620 curHist=(THnSparseF*)fArrayDz->FindObject(name);
74235403 621 if(!curHist){
dde68d36 622 curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
623 fArrayDz->AddLast(curHist);
74235403 624 }
dde68d36 625// curHist=(THnSparseF*)(fMapDz->GetValue(event->GetFiredTriggerClasses()));
626// if(!curHist){
627// curHist=new THnSparseF(event->GetFiredTriggerClasses(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
628// fMapDz->Add(new TObjString(event->GetFiredTriggerClasses()),curHist);
629// }
630 curHist->Fill(vecDrift);
74235403 631
dde68d36 632 name=shortName;
633 name+="ALL";
634 name.ToUpper();
635 curHist=(THnSparseF*)fArrayDz->FindObject(name);
74235403 636 if(!curHist){
dde68d36 637 curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
638 fArrayDz->AddLast(curHist);
74235403 639 }
dde68d36 640// curHist=(THnSparseF*)(fMapDz->GetValue("all"));
641// if(!curHist){
642// curHist=new THnSparseF("all","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
643// fMapDz->Add(new TObjString("all"),curHist);
644// }
645 curHist->Fill(vecDrift);
74235403 646 }
647 }
d3ce44cb 648 TTreeSRedirector *cstream = GetDebugStreamer();
649 if (fStreamLevel>0){
650 if (cstream){
651 (*cstream)<<"trackInfo"<<
652 "tr0.="<<track0<<
653 "tr1.="<<track1<<
654 "p0.="<<&param0<<
655 "p1.="<<&param1<<
656 "isPair="<<isPair<<
657 "isCross="<<isCross<<
dde68d36 658 "isSame="<<isSame<<
d3ce44cb 659 "fDz="<<fDz<<
660 "fRun="<<fRun<<
661 "fTime="<<fTime<<
662 "\n";
663 }
664 }
c74c5f6c 665 } // end 2nd order loop
666 } // end 1st order loop
74235403 667
2be25cc9 668 if (fStreamLevel>0){
669 TTreeSRedirector *cstream = GetDebugStreamer();
670 if (cstream){
671 TTimeStamp tstamp(fTime);
672 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
673 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
674 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
675 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
676 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
677 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
678 TVectorD vecGoofie(20);
679 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
680 if (goofieArray){
681 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
682 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
683 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
684 }
685 }
686 (*cstream)<<"timeInfo"<<
687 "run="<<fRun<< // run number
688 "event="<<fEvent<< // event number
689 "time="<<fTime<< // time stamp of event
690 "trigger="<<fTrigger<< // trigger
691 "mag="<<fMagF<< // magnetic field
692 // Environment values
693 "press0="<<valuePressure0<<
694 "press1="<<valuePressure1<<
695 "pt0="<<ptrelative0<<
696 "pt1="<<ptrelative1<<
697 "temp0="<<temp0<<
698 "temp1="<<temp1<<
699 "vecGoofie.=<<"<<&vecGoofie<<
700 //
701 // accumulated values
702 //
703 "fDz="<<fDz<< //! current delta z
2be25cc9 704 "trigger="<<event->GetFiredTriggerClasses()<<
705 "\n";
706 }
707 }
708 printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
c74c5f6c 709}
710
dde68d36 711void AliTPCcalibTime::ProcessBeam(AliESDEvent */*event*/){
712}
713
74235403 714void AliTPCcalibTime::Analyze(){}
715
74235403 716THnSparse* AliTPCcalibTime::GetHistoDrift(const char* name){
dde68d36 717 TIterator* iterator = fArrayDz->MakeIterator();
718 iterator->Reset();
719 TString newName=name;
720 newName.ToUpper();
721 THnSparse* newHist=new THnSparseF(newName,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
722 THnSparse* addHist=NULL;
723 while((addHist=(THnSparseF*)iterator->Next())){
724 if(!addHist) continue;
725 TString histName=addHist->GetName();
726 if(!histName.Contains(newName)) continue;
727 addHist->Print();
728 newHist->Add(addHist);
74235403 729 }
dde68d36 730 return newHist;
74235403 731}
732
dde68d36 733TObjArray* AliTPCcalibTime::GetHistoDrift(){
734 return fArrayDz;
74235403 735}
736
737TGraphErrors* AliTPCcalibTime::GetGraphDrift(const char* name){
dde68d36 738 THnSparse* histoDrift=GetHistoDrift(name);
739 TGraphErrors* graphDrift=NULL;
740 if(histoDrift){
741 graphDrift=FitSlices(histoDrift,2,0,400,100,0.05,0.95, kTRUE);
742 TString end=histoDrift->GetName();
743 Int_t pos=end.Index("_");
744 end=end(pos,end.Capacity()-pos);
745 TString graphName=graphDrift->ClassName();
746 graphName+=end;
747 graphName.ToUpper();
748 graphDrift->SetName(graphName);
74235403 749 }
750 return graphDrift;
751}
752
dde68d36 753TObjArray* AliTPCcalibTime::GetGraphDrift(){
754 TObjArray* arrayGraphDrift=new TObjArray();
755 TIterator* iterator=fArrayDz->MakeIterator();
74235403 756 iterator->Reset();
dde68d36 757 THnSparse* addHist=NULL;
758 while((addHist=(THnSparseF*)iterator->Next())) arrayGraphDrift->AddLast(GetGraphDrift(addHist->GetName()));
759 return arrayGraphDrift;
74235403 760}
761
dde68d36 762AliSplineFit* AliTPCcalibTime::GetFitDrift(const char* name){
d3ce44cb 763 TGraph* graphDrift=GetGraphDrift(name);
dde68d36 764 AliSplineFit* fitDrift=NULL;
74235403 765 if(graphDrift && graphDrift->GetN()){
dde68d36 766 fitDrift=new AliSplineFit();
767 fitDrift->SetGraph(graphDrift);
768 fitDrift->SetMinPoints(graphDrift->GetN()+1);
769 fitDrift->InitKnots(graphDrift,2,0,0.001);
770 fitDrift->SplineFit(0);
771 TString end=graphDrift->GetName();
772 Int_t pos=end.Index("_");
773 end=end(pos,end.Capacity()-pos);
774 TString fitName=fitDrift->ClassName();
775 fitName+=end;
776 fitName.ToUpper();
777 //fitDrift->SetName(fitName);
74235403 778 delete graphDrift;
dde68d36 779 graphDrift=NULL;
74235403 780 }
781 return fitDrift;
782}
783
dde68d36 784//TObjArray* AliTPCcalibTime::GetFitDrift(){
785// TObjArray* arrayFitDrift=new TObjArray();
786// TIterator* iterator = fArrayDz->MakeIterator();
787// iterator->Reset();
788// THnSparse* addHist=NULL;
789// while((addHist=(THnSparseF*)iterator->Next())) arrayFitDrift->AddLast(GetFitDrift(addHist->GetName()));
790// return arrayFitDrift;
791//}
74235403 792
c74c5f6c 793Long64_t AliTPCcalibTime::Merge(TCollection *li) {
c74c5f6c 794 TIterator* iter = li->MakeIterator();
795 AliTPCcalibTime* cal = 0;
796
797 while ((cal = (AliTPCcalibTime*)iter->Next())) {
798 if (!cal->InheritsFrom(AliTPCcalibTime::Class())) {
799 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
800 return -1;
801 }
2be25cc9 802 for (Int_t imeas=0; imeas<3; imeas++){
803 if (cal->GetHistVdriftLaserA(imeas) && cal->GetHistVdriftLaserA(imeas)){
804 fHistVdriftLaserA[imeas]->Add(cal->GetHistVdriftLaserA(imeas));
805 fHistVdriftLaserC[imeas]->Add(cal->GetHistVdriftLaserC(imeas));
806 }
807 }
dde68d36 808 TObjArray* addArray=cal->GetHistoDrift();
809 if(!addArray) return 0;
810 TIterator* iterator = addArray->MakeIterator();
2be25cc9 811 iterator->Reset();
dde68d36 812 THnSparse* addHist=NULL;
813 while((addHist=(THnSparseF*)iterator->Next())){
814 if(!addHist) continue;
2be25cc9 815 addHist->Print();
dde68d36 816 THnSparse* localHist=(THnSparseF*)fArrayDz->FindObject(addHist->GetName());
2be25cc9 817 if(!localHist){
818 localHist=new THnSparseF(addHist->GetName(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
dde68d36 819 fArrayDz->AddLast(localHist);
2be25cc9 820 }
821 localHist->Add(addHist);
822 }
dde68d36 823// TMap * addMap=cal->GetHistoDrift();
824// if(!addMap) return 0;
825// TIterator* iterator = addMap->MakeIterator();
826// iterator->Reset();
827// TPair* addPair=0;
828// while((addPair=(TPair *)(addMap->FindObject(iterator->Next())))){
829// THnSparse* addHist=dynamic_cast<THnSparseF*>(addPair->Value());
830// if (!addHist) continue;
831// addHist->Print();
832// THnSparse* localHist=dynamic_cast<THnSparseF*>(fMapDz->GetValue(addHist->GetName()));
833// if(!localHist){
834// localHist=new THnSparseF(addHist->GetName(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
835// fMapDz->Add(new TObjString(addHist->GetName()),localHist);
836// }
837// localHist->Add(addHist);
838// }
d3ce44cb 839 for(Int_t i=0;i<10;i++) if (cal->GetCosmiMatchingHisto(i)) fCosmiMatchingHisto[i]->Add(cal->GetCosmiMatchingHisto(i));
c74c5f6c 840 }
c74c5f6c 841 return 0;
c74c5f6c 842}
843
c74c5f6c 844Bool_t AliTPCcalibTime::IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1){
c74c5f6c 845 /*
846 // 0. Same direction - OPOSITE - cutDir +cutT
847 TCut cutDir("cutDir","dir<-0.99")
848 // 1.
849 TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03")
850 //
851 // 2. The same rphi
852 TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5")
853 //
c74c5f6c 854 TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
855 // 1/Pt diff cut
856 */
857 const Double_t *p0 = tr0->GetParameter();
858 const Double_t *p1 = tr1->GetParameter();
74235403 859 fCosmiMatchingHisto[0]->Fill(p0[0]+p1[0]);
860 fCosmiMatchingHisto[1]->Fill(p0[1]-p1[1]);
d3ce44cb 861 fCosmiMatchingHisto[2]->Fill(tr1->GetAlpha()-tr0->GetAlpha()+TMath::Pi());
74235403 862 fCosmiMatchingHisto[3]->Fill(p0[3]+p1[3]);
863 fCosmiMatchingHisto[4]->Fill(p0[4]+p1[4]);
864
c74c5f6c 865 if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
866 if (TMath::Abs(p0[0]+p1[0])>fCutMaxD) return kFALSE;
da6c0bc9 867 if (TMath::Abs(p0[1]-p1[1])>fCutMaxDz) return kFALSE;
c74c5f6c 868 Double_t d0[3], d1[3];
869 tr0->GetDirection(d0);
870 tr1->GetDirection(d1);
871 if (d0[0]*d1[0] + d0[1]*d1[1] + d0[2]*d1[2] >fCutMinDir) return kFALSE;
74235403 872
d3ce44cb 873 fCosmiMatchingHisto[5]->Fill(p0[0]+p1[0]);
874 fCosmiMatchingHisto[6]->Fill(p0[1]-p1[1]);
875 fCosmiMatchingHisto[7]->Fill(tr1->GetAlpha()-tr0->GetAlpha()+TMath::Pi());
876 fCosmiMatchingHisto[8]->Fill(p0[3]+p1[3]);
877 fCosmiMatchingHisto[9]->Fill(p0[4]+p1[4]);
878
c74c5f6c 879 return kTRUE;
880}
d3ce44cb 881Bool_t AliTPCcalibTime::IsCross(AliESDtrack *tr0, AliESDtrack *tr1){
882 return tr0->GetOuterParam()->GetZ()*tr1->GetOuterParam()->GetZ()<0 && tr0->GetInnerParam()->GetZ()*tr1->GetInnerParam()->GetZ()<0 && tr0->GetOuterParam()->GetZ()*tr0->GetInnerParam()->GetZ()>0 && tr1->GetOuterParam()->GetZ()*tr1->GetInnerParam()->GetZ()>0;
883}
884
dde68d36 885Bool_t AliTPCcalibTime::IsSame(AliESDtrack */*tr0*/, AliESDtrack */*tr1*/){
886 // To be implemented
887 return kFALSE;
888}
889
d3ce44cb 890/*
891chainDrift->Draw("p0.fP[0]+p1.fP[0]","isPair");
892 mean ~-0.02 ~-0.03913
893 RMS ~ 0.5 ~ 0.5356 --> 3 (fCutMaxD)
894
895chainDrift->Draw("p0.fP[1]-p1.fP[1]","isPair");
896 mean ~ 0.1855
897 RMS ~ 4.541 -->25 (fCutMaxDz)
898
899chainDrift->Draw("p1.fAlpha-p0.fAlpha+pi","isPair");
900//chainDrift->Draw("p1.fAlpha+p0.fAlpha","isPair");
901//chainDrift->Draw("p1.fP[2]-p0.fP[2]+pi","isPair");
902//chainDrift->Draw("p1.fP[2]+p0.fP[2]","isPair");
903 mean ~ 0 ~ 0.001898
904 RMS ~ 0.009 ~ 0.01134 --> 0.06
905
906chainDrift->Draw("p0.fP[3]+p1.fP[3]","isPair");
907 mean ~ 0.0013 ~ 0.001539
908 RMS ~ 0.003 ~ 0.004644 --> 0.03 (fCutTheta)
909
910chainDrift->Draw("p1.fP[4]+p0.fP[4]>>his(100,-0.2,0.2)","isPair")
911 mean ~ 0.012 ~-0.0009729
912 RMS ~ 0.036 ~ 0.03773 --> 0.2
913*/
31aa7c5c 914
817766d5 915
916void AliTPCcalibTime::ProcessAlignITS(AliESDtrack* track, AliESDfriendTrack *friendTrack){
917 //
918 // Process track
919 // Update TPC-ITS alignment
920 //
921 const Int_t kMinTPC = 80;
922 const Int_t kMinITS = 3;
923 const Double_t kMinZ = 10;
924 const Double_t kMaxDy = 2;
925 const Double_t kMaxAngle= 0.02;
926 //
927 Int_t dummycl[1000];
928 if (track->GetITSclusters(dummycl)<kMinITS) return; // minimal amount of clusters
929 if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
930 //
931 if (!friendTrack->GetITSOut()) return;
932 if (!track->GetInnerParam()) return;
933 if (!track->GetOuterParam()) return;
934 // exclude crossing track
935 if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
936 if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ) return;
937 //
938 AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(track->GetInnerParam()));
939 AliExternalTrackParam pITS(*(friendTrack->GetITSOut()));
940 //
941 //
942 //
943 Int_t htime = fTime/3600; //time in hours
944 if (fAlignITSTPC->GetEntries()<htime){
945 fAlignITSTPC->Expand(htime*2+20);
946 }
947 AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignITSTPC->At(htime);
948 if (!align){
949 align=new AliRelAlignerKalman();
950 align->SetOutRejSigma(2.);
951 //align->SetRejectOutliers(kFALSE);
952 align->SetMagField(fMagF);
953 fAlignITSTPC->AddAt(align,htime);
954 }
955 pITS.Rotate(pTPC.GetAlpha());
956 pITS.PropagateTo(pTPC.GetX(),fMagF);
957 if (TMath::Abs(pITS.GetY()-pTPC.GetY())>kMaxDy) return;
958 if (TMath::Abs(pITS.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
959 if (TMath::Abs(pITS.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
960 align->AddTrackParams(&pITS,&pTPC);
961 align->SetTimeStamp(fTime);
817766d5 962 // align->SetRunNumber(fRun );
963 static Int_t entry=-1;
964 entry++;
a66da2d4 965 // Int_t nupdates=align->GetNUpdates();
966 Int_t nupdates=entry;
817766d5 967 align->SetOutRejSigma(1.+1./Double_t(nupdates));
968 TTreeSRedirector *cstream = GetDebugStreamer();
969 if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
970 TTimeStamp tstamp(fTime);
971 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
972 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
973 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
974 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
975 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
976 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
977 TVectorD vecGoofie(20);
978 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
979 if (goofieArray){
980 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
981 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
982 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
983 }
984 }
985 TVectorD gpTPC(3), gdTPC(3);
986 TVectorD gpITS(3), gdITS(3);
987 pTPC.GetXYZ(gpTPC.GetMatrixArray());
988 pTPC.GetDirection(gdTPC.GetMatrixArray());
989 pITS.GetXYZ(gpITS.GetMatrixArray());
990 pITS.GetDirection(gdITS.GetMatrixArray());
991 (*cstream)<<"itstpc"<<
992 "run="<<fRun<< // run number
993 "event="<<fEvent<< // event number
994 "time="<<fTime<< // time stamp of event
995 "trigger="<<fTrigger<< // trigger
996 "mag="<<fMagF<< // magnetic field
997 // Environment values
998 "press0="<<valuePressure0<<
999 "press1="<<valuePressure1<<
1000 "pt0="<<ptrelative0<<
1001 "pt1="<<ptrelative1<<
1002 "temp0="<<temp0<<
1003 "temp1="<<temp1<<
1004 "vecGoofie.="<<&vecGoofie<<
1005 "entry="<<entry<< // current entry
1006 //
1007 "a.="<<align<< // current alignment
1008 "pITS.="<<&pITS<< // track param ITS
1009 "pTPC.="<<&pTPC<< // track param TPC
1010 "gpTPC.="<<&gpTPC<<
1011 "gdTPC.="<<&gdTPC<<
1012 "gpITS.="<<&gpITS<<
1013 "gdITS.="<<&gdITS<<
1014 "\n";
1015 }
1016
1017}
1018
1019
1020void AliTPCcalibTime::ProcessAlignTOF(AliESDtrack* track, AliESDfriendTrack *friendTrack){
1021 //
1022 //process TOF-TPC alignment
1023 //
1024 Int_t kminNcl=80;
1025 Float_t kMaxDy=6;
1026 Float_t kMaxDz=10;
1027 if (track->GetTPCNcls()<kminNcl) return;
1028 if (track->GetOuterParam()==0) return;
1029 if (track->GetInnerParam()==0) return;
1030 if (track->GetTOFsignal()<=0) return;
1031 //
1032 AliExternalTrackParam *paramOut = new AliExternalTrackParam(*(track->GetOuterParam()));
1033 AliExternalTrackParam *param=0;
1034 const AliTrackPointArray *points=friendTrack->GetTrackPointArray();
1035 if (!points) return;
1036 Int_t npoints = points->GetNPoints();
1037 AliTrackPoint point;
1038 //Double_t alpha=
1039 Double_t mass = TDatabasePDG::Instance()->GetParticle("mu+")->Mass();
1040 TTreeSRedirector * cstream = GetDebugStreamer();
1041 //
1042 //
1043 //
1044 for (Int_t ipoint=0;ipoint<npoints;ipoint++){
1045 //
1046 points->GetPoint(point,ipoint);
1047 Float_t xyz[3];
1048 point.GetXYZ(xyz);
1049 Double_t r=TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
1050 if (r<300) continue;
1051 if (r>400) continue;
1052 param=paramOut;
1053 if (!param) continue;
1054 AliTracker::PropagateTrackToBxByBz(param,r,mass,2.,kTRUE);
1055 AliTracker::PropagateTrackToBxByBz(param,r,mass,0.1,kTRUE);
1056 AliTrackPoint lpoint = point.Rotate(param->GetAlpha());
1057 param->PropagateTo(lpoint.GetX(),fMagF);
1058 //
1059 //
1060 // this is ugly - we need AliCluster constructor
1061 //
1062 AliExternalTrackParam &pTPC=*param;
1063 AliExternalTrackParam pTOF(*param);
1064 ((Double_t*)pTOF.GetParameter())[0] =lpoint.GetY();
1065 ((Double_t*)pTOF.GetParameter())[1] =lpoint.GetZ();
1066 pTOF.ResetCovariance(20);
1067 ((Double_t*)pTOF.GetCovariance())[0]+=3.*3.;
1068 ((Double_t*)pTOF.GetCovariance())[2]+=3.*3.;
1069 if (TMath::Abs(pTOF.GetY()-pTPC.GetY())>kMaxDy) continue;
1070 if (TMath::Abs(pTOF.GetZ()-pTPC.GetZ())>kMaxDz) continue;
1071 //
1072 Int_t htime = fTime/3600; //time in hours
1073
1074 if (fAlignTOFTPC->GetEntries()<htime){
1075 fAlignTOFTPC->Expand(htime*2+20);
1076 }
1077 AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTOFTPC->At(htime);
1078 if (!align){
1079 align=new AliRelAlignerKalman();
1080 align->SetOutRejSigma(2.);
1081 //align->SetRejectOutliers(kFALSE);
1082 align->SetMagField(fMagF);
1083 fAlignTOFTPC->AddAt(align,htime);
1084 }
1085 pTOF.Rotate(pTPC.GetAlpha());
1086 pTOF.PropagateTo(pTPC.GetX(),fMagF);
1087 align->AddTrackParams(&pTOF,&pTPC);
1088 align->SetTimeStamp(fTime);
817766d5 1089 static Int_t entry=-1;
1090 entry++;
a66da2d4 1091 // Int_t nupdates=align->GetNUpdates();
1092 Int_t nupdates=entry;
817766d5 1093 align->SetOutRejSigma(1.+1./Double_t(nupdates));
1094
1095 //
1096 //
1097 if (cstream) {
1098 (*cstream) << "pointMatch" <<
1099 "run="<<fRun<< // run number
1100 "event="<<fEvent<< // event number
1101 "time="<<fTime<< // time stamp of event
1102 "trigger="<<fTrigger<< // trigger
1103 "mag="<<fMagF<< // magnetic field
1104 //
1105 "a.="<<align<< // current alignment
1106 "p.="<<&point<<
1107 "lp.="<<&lpoint<<
1108 "pTPC.="<<&pTPC<<
1109 "pTOF.="<<&pTOF<<
1110 "\n";
1111 }
1112
1113
1114
1115 }
1116 delete paramOut;
1117}