Adding DDL corruption recovery option to dHLTdumpraw and making modifications require...
[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;
3e55050f 54 TChain * chainTime = tool.MakeChainRandom("time.txt","trackInfo",0,10000);
817766d5 55
3e55050f 56 AliXRDPROOFtoolkit::FilterList("timetpctpc.txt","* tpctpc",1)
57 AliXRDPROOFtoolkit::FilterList("timetoftpc.txt","* toftpc",1)
58 AliXRDPROOFtoolkit::FilterList("timeitstpc.txt","* itstpc",1)
59 AliXRDPROOFtoolkit::FilterList("timelaser.txt","* laserInfo",1)
60 TChain * chainTPCTPC = tool.MakeChainRandom("timetpctpc.txt.Good","tpctpc",0,10000);
817766d5 61 TChain * chainTPCITS = tool.MakeChainRandom("timeitstpc.txt.Good","itstpc",0,10000);
3e55050f 62 TChain * chainTPCTOF = tool.MakeChainRandom("timetoftpc.txt.Good","toftpc",0,10000);
817766d5 63 TChain * chainLaser = tool.MakeChainRandom("timelaser.txt.Good","laserInfo",0,10000);
74235403 64 chainTime->Lookup();
dde68d36 65 chainLaser->Lookup();
c74c5f6c 66*/
67
c74c5f6c 68#include "Riostream.h"
2811495d 69#include "TDatabasePDG.h"
70#include "TGraphErrors.h"
c74c5f6c 71#include "TH1F.h"
c74c5f6c 72#include "THnSparse.h"
73#include "TList.h"
74#include "TMath.h"
2811495d 75#include "TTimeStamp.h"
76#include "TTree.h"
c74c5f6c 77#include "TVectorD.h"
2811495d 78//#include "TChain.h"
79//#include "TFile.h"
80
81#include "AliDCSSensor.h"
82#include "AliDCSSensorArray.h"
c74c5f6c 83#include "AliESDEvent.h"
c74c5f6c 84#include "AliESDInputHandler.h"
2811495d 85#include "AliESDVertex.h"
86#include "AliESDfriend.h"
c74c5f6c 87#include "AliLog.h"
817766d5 88#include "AliRelAlignerKalman.h"
2811495d 89#include "AliTPCCalROC.h"
90#include "AliTPCParam.h"
c74c5f6c 91#include "AliTPCTracklet.h"
da6c0bc9 92#include "AliTPCcalibDB.h"
93#include "AliTPCcalibLaser.h"
2811495d 94#include "AliTPCcalibTime.h"
95#include "AliTPCclusterMI.h"
96#include "AliTPCseed.h"
817766d5 97#include "AliTrackPointArray.h"
2811495d 98#include "AliTracker.h"
817766d5 99
c74c5f6c 100ClassImp(AliTPCcalibTime)
101
102
103AliTPCcalibTime::AliTPCcalibTime()
da6c0bc9 104 :AliTPCcalibBase(),
da6c0bc9 105 fLaser(0), // pointer to laser calibration
106 fDz(0), // current delta z
d3ce44cb 107 fCutMaxD(3), // maximal distance in rfi ditection
108 fCutMaxDz(25), // maximal distance in rfi ditection
c74c5f6c 109 fCutTheta(0.03), // maximal distan theta
2be25cc9 110 fCutMinDir(-0.99), // direction vector products
b842d904 111 fCutTracks(100),
dde68d36 112 fArrayDz(0), //NEW! Tmap of V drifts for different triggers
817766d5 113 fAlignITSTPC(0), //alignemnt array ITS TPC match
114 fAlignTRDTPC(0), //alignemnt array TRD TPC match
115 fAlignTOFTPC(0), //alignemnt array TOF TPC match
2be25cc9 116 fTimeBins(0),
117 fTimeStart(0),
118 fTimeEnd(0),
119 fPtBins(0),
120 fPtStart(0),
121 fPtEnd(0),
122 fVdriftBins(0),
123 fVdriftStart(0),
124 fVdriftEnd(0),
125 fRunBins(0),
126 fRunStart(0),
127 fRunEnd(0)
c74c5f6c 128{
6ba68367 129 //
130 // default constructor
131 //
c74c5f6c 132 AliInfo("Default Constructor");
2be25cc9 133 for (Int_t i=0;i<3;i++) {
134 fHistVdriftLaserA[i]=0;
135 fHistVdriftLaserC[i]=0;
136 }
d3ce44cb 137 for (Int_t i=0;i<10;i++) {
74235403 138 fCosmiMatchingHisto[i]=0;
139 }
391ffdb2 140 //
141 for (Int_t i=0;i<5;i++) {
142 fResHistoTPCITS[i]=0;
143 fResHistoTPCTRD[i]=0;
144 fResHistoTPCvertex[i]=0;
145 }
146
c74c5f6c 147}
148
74235403 149AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeVdrift)
c74c5f6c 150 :AliTPCcalibBase(),
d3ce44cb 151 fLaser(0), // pointer to laser calibration
152 fDz(0), // current delta z
153 fCutMaxD(5*0.5356), // maximal distance in rfi ditection
dde68d36 154 fCutMaxDz(40), // maximal distance in rfi ditection
d3ce44cb 155 fCutTheta(5*0.004644),// maximal distan theta
156 fCutMinDir(-0.99), // direction vector products
b842d904 157 fCutTracks(100),
dde68d36 158 fArrayDz(0), //Tmap of V drifts for different triggers
817766d5 159 fAlignITSTPC(0), //alignemnt array ITS TPC match
160 fAlignTRDTPC(0), //alignemnt array TRD TPC match
161 fAlignTOFTPC(0), //alignemnt array TOF TPC match
2be25cc9 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{
6ba68367 175 //
176 // Non deafaul constructor - to be used in the Calibration setups
177 //
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
391ffdb2 186 for (Int_t i=0;i<5;i++) {
187 fResHistoTPCITS[i]=0;
188 fResHistoTPCTRD[i]=0;
189 fResHistoTPCvertex[i]=0;
190 }
191
192
c74c5f6c 193 AliInfo("Non Default Constructor");
74235403 194 fTimeBins =(EndTime-StartTime)/deltaIntegrationTimeVdrift;
2be25cc9 195 fTimeStart =StartTime; //(((TObjString*)(mapGRP->GetValue("fAliceStartTime")))->GetString()).Atoi();
196 fTimeEnd =EndTime; //(((TObjString*)(mapGRP->GetValue("fAliceStopTime")))->GetString()).Atoi();
dde68d36 197 fPtBins = 400;
2be25cc9 198 fPtStart = -0.04;
199 fPtEnd = 0.04;
dde68d36 200 fVdriftBins = 500;
201 fVdriftStart= -0.1;
202 fVdriftEnd = 0.1;
32100b1d 203 fRunBins = 1000001;
dde68d36 204 fRunStart = -1.5;
32100b1d 205 fRunEnd = 999999.5;
2be25cc9 206
207 Int_t binsVdriftLaser[4] = {fTimeBins , fPtBins , fVdriftBins*20, fRunBins };
208 Double_t xminVdriftLaser[4] = {fTimeStart, fPtStart, fVdriftStart , fRunStart};
209 Double_t xmaxVdriftLaser[4] = {fTimeEnd , fPtEnd , fVdriftEnd , fRunEnd };
dde68d36 210 TString axisTitle[4]={
211 "T",
212 "#delta_{P/T}",
213 "value",
214 "run"
215 };
216 TString histoName[3]={
217 "Loffset",
218 "Lcorr",
219 "Lgy"
220 };
2be25cc9 221
dde68d36 222
2be25cc9 223 for (Int_t i=0;i<3;i++) {
224 fHistVdriftLaserA[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
225 fHistVdriftLaserC[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
dde68d36 226 fHistVdriftLaserA[i]->SetName(histoName[i]);
227 fHistVdriftLaserC[i]->SetName(histoName[i]);
228 for (Int_t iaxis=0; iaxis<4;iaxis++){
229 fHistVdriftLaserA[i]->GetAxis(iaxis)->SetName(axisTitle[iaxis]);
230 fHistVdriftLaserC[i]->GetAxis(iaxis)->SetName(axisTitle[iaxis]);
231 }
2be25cc9 232 }
233 fBinsVdrift[0] = fTimeBins;
234 fBinsVdrift[1] = fPtBins;
235 fBinsVdrift[2] = fVdriftBins;
236 fBinsVdrift[3] = fRunBins;
237 fXminVdrift[0] = fTimeStart;
238 fXminVdrift[1] = fPtStart;
239 fXminVdrift[2] = fVdriftStart;
240 fXminVdrift[3] = fRunStart;
241 fXmaxVdrift[0] = fTimeEnd;
242 fXmaxVdrift[1] = fPtEnd;
243 fXmaxVdrift[2] = fVdriftEnd;
244 fXmaxVdrift[3] = fRunEnd;
245
dde68d36 246 fArrayDz=new TObjArray();
817766d5 247 fAlignITSTPC = new TObjArray; //alignemnt array ITS TPC match
248 fAlignTRDTPC = new TObjArray; //alignemnt array ITS TPC match
249 fAlignTOFTPC = new TObjArray; //alignemnt array ITS TPC match
3e55050f 250 fAlignITSTPC->SetOwner(kTRUE);
251 fAlignTRDTPC->SetOwner(kTRUE);
252 fAlignTOFTPC->SetOwner(kTRUE);
253
c74c5f6c 254
d3ce44cb 255 fCosmiMatchingHisto[0]=new TH1F("Cosmics matching","p0-all" ,100,-10*0.5356 ,10*0.5356 );
256 fCosmiMatchingHisto[1]=new TH1F("Cosmics matching","p1-all" ,100,-10*4.541 ,10*4.541 );
257 fCosmiMatchingHisto[2]=new TH1F("Cosmics matching","p2-all" ,100,-10*0.01134 ,10*0.01134 );
258 fCosmiMatchingHisto[3]=new TH1F("Cosmics matching","p3-all" ,100,-10*0.004644,10*0.004644);
259 fCosmiMatchingHisto[4]=new TH1F("Cosmics matching","p4-all" ,100,-10*0.03773 ,10*0.03773 );
260 fCosmiMatchingHisto[5]=new TH1F("Cosmics matching","p0-isPair",100,-10*0.5356 ,10*0.5356 );
261 fCosmiMatchingHisto[6]=new TH1F("Cosmics matching","p1-isPair",100,-10*4.541 ,10*4.541 );
262 fCosmiMatchingHisto[7]=new TH1F("Cosmics matching","p2-isPair",100,-10*0.01134 ,10*0.01134 );
263 fCosmiMatchingHisto[8]=new TH1F("Cosmics matching","p3-isPair",100,-10*0.004644,10*0.004644);
264 fCosmiMatchingHisto[9]=new TH1F("Cosmics matching","p4-isPair",100,-10*0.03773 ,10*0.03773 );
391ffdb2 265 BookDistortionMaps();
74235403 266}
c74c5f6c 267
268AliTPCcalibTime::~AliTPCcalibTime(){
269 //
6ba68367 270 // Virtual Destructor
c74c5f6c 271 //
d3ce44cb 272 for(Int_t i=0;i<3;i++){
273 if(fHistVdriftLaserA[i]){
274 delete fHistVdriftLaserA[i];
275 fHistVdriftLaserA[i]=NULL;
276 }
277 if(fHistVdriftLaserC[i]){
278 delete fHistVdriftLaserC[i];
279 fHistVdriftLaserC[i]=NULL;
280 }
281 }
dde68d36 282 if(fArrayDz){
283 fArrayDz->SetOwner();
284 fArrayDz->Delete();
285 delete fArrayDz;
286 fArrayDz=NULL;
2be25cc9 287 }
d3ce44cb 288 for(Int_t i=0;i<5;i++){
289 if(fCosmiMatchingHisto[i]){
290 delete fCosmiMatchingHisto[i];
291 fCosmiMatchingHisto[i]=NULL;
292 }
74235403 293 }
391ffdb2 294
295 for (Int_t i=0;i<5;i++) {
296 delete fResHistoTPCITS[i];
297 delete fResHistoTPCTRD[i];
298 delete fResHistoTPCvertex[i];
299 fResHistoTPCITS[i]=0;
300 fResHistoTPCTRD[i]=0;
301 fResHistoTPCvertex[i]=0;
302 }
303
304
3e55050f 305 fAlignITSTPC->SetOwner(kTRUE);
306 fAlignTRDTPC->SetOwner(kTRUE);
307 fAlignTOFTPC->SetOwner(kTRUE);
308
817766d5 309 fAlignITSTPC->Delete();
310 fAlignTRDTPC->Delete();
311 fAlignTOFTPC->Delete();
3e55050f 312 delete fAlignITSTPC;
313 delete fAlignTRDTPC;
314 delete fAlignTOFTPC;
c74c5f6c 315}
2be25cc9 316
2811495d 317Bool_t AliTPCcalibTime::IsLaser(const AliESDEvent *const /*event*/){
6ba68367 318 //
319 // Indicator is laser event not yet implemented - to be done using trigger info or event specie
320 //
dde68d36 321 return kTRUE; //More accurate creteria to be added
da6c0bc9 322}
2811495d 323Bool_t AliTPCcalibTime::IsCosmics(const AliESDEvent *const /*event*/){
6ba68367 324 //
325 // Indicator is cosmic event not yet implemented - to be done using trigger info or event specie
326 //
327
dde68d36 328 return kTRUE; //More accurate creteria to be added
329}
2811495d 330Bool_t AliTPCcalibTime::IsBeam(const AliESDEvent *const /*event*/){
6ba68367 331 //
332 // Indicator is physic event not yet implemented - to be done using trigger info or event specie
333 //
334
dde68d36 335 return kTRUE; //More accurate creteria to be added
336}
337void AliTPCcalibTime::ResetCurrent(){
338 fDz=0; //Reset current dz
2be25cc9 339}
6ba68367 340
341
342
2be25cc9 343void AliTPCcalibTime::Process(AliESDEvent *event){
6ba68367 344 //
345 // main function to make calibration
346 //
2be25cc9 347 if(!event) return;
348 if (event->GetNumberOfTracks()<2) return;
349 ResetCurrent();
dde68d36 350 if(IsLaser (event)) ProcessLaser (event);
351 if(IsCosmics(event)) ProcessCosmic(event);
352 if(IsBeam (event)) ProcessBeam (event);
2be25cc9 353}
354
355void AliTPCcalibTime::ProcessLaser(AliESDEvent *event){
c74c5f6c 356 //
2be25cc9 357 // Fit drift velocity using laser
358 //
359 // 0. cuts
dde68d36 360 const Int_t kMinTracks = 40; // minimal number of laser tracks
361 const Int_t kMinTracksSide = 20; // minimal number of tracks per side
362 const Float_t kMaxDeltaZ = 30.; // maximal trigger delay
363 const Float_t kMaxDeltaV = 0.05; // maximal deltaV
2be25cc9 364 const Float_t kMaxRMS = 0.1; // maximal RMS of tracks
dde68d36 365 //
2be25cc9 366 /*
367 TCut cutRMS("sqrt(laserA.fElements[4])<0.1&&sqrt(laserC.fElements[4])<0.1");
368 TCut cutZ("abs(laserA.fElements[0]-laserC.fElements[0])<3");
369 TCut cutV("abs(laserA.fElements[1]-laserC.fElements[1])<0.01");
370 TCut cutY("abs(laserA.fElements[2]-laserC.fElements[2])<2");
371 TCut cutAll = cutRMS+cutZ+cutV+cutY;
372 */
373 if (event->GetNumberOfTracks()<kMinTracks) return;
c74c5f6c 374 //
2be25cc9 375 if(!fLaser) fLaser = new AliTPCcalibLaser("laserTPC","laserTPC",kFALSE);
376 fLaser->Process(event);
377 if (fLaser->GetNtracks()<kMinTracks) return; // small amount of tracks cut
dde68d36 378 if (fLaser->fFitAside->GetNrows()==0 && fLaser->fFitCside->GetNrows()==0) return; // no fit neither a or C side
c74c5f6c 379 //
2be25cc9 380 // debug streamer - activate stream level
381 // Use it for tuning of the cuts
da6c0bc9 382 //
dde68d36 383 // cuts to be applied
384 //
385 Int_t isReject[2]={0,0};
386 //
387 // not enough tracks
388 if (TMath::Abs((*fLaser->fFitAside)[3]) < kMinTracksSide) isReject[0]|=1;
389 if (TMath::Abs((*fLaser->fFitCside)[3]) < kMinTracksSide) isReject[1]|=1;
390 // unreasonable z offset
391 if (TMath::Abs((*fLaser->fFitAside)[0])>kMaxDeltaZ) isReject[0]|=2;
392 if (TMath::Abs((*fLaser->fFitCside)[0])>kMaxDeltaZ) isReject[1]|=2;
393 // unreasonable drift velocity
394 if (TMath::Abs((*fLaser->fFitAside)[1]-1)>kMaxDeltaV) isReject[0]|=4;
395 if (TMath::Abs((*fLaser->fFitCside)[1]-1)>kMaxDeltaV) isReject[1]|=4;
396 // big chi2
397 if (TMath::Sqrt(TMath::Abs((*fLaser->fFitAside)[4]))>kMaxRMS ) isReject[0]|=8;
398 if (TMath::Sqrt(TMath::Abs((*fLaser->fFitCside)[4]))>kMaxRMS ) isReject[1]|=8;
399
400
401
402
2be25cc9 403 if (fStreamLevel>0){
404 printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
405
da6c0bc9 406 TTreeSRedirector *cstream = GetDebugStreamer();
407 if (cstream){
408 TTimeStamp tstamp(fTime);
2be25cc9 409 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
da6c0bc9 410 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
411 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
412 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
413 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
414 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
32100b1d 415 Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
31aa7c5c 416 TVectorD vecGoofie(20);
417 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
2be25cc9 418 if (goofieArray){
31aa7c5c 419 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
420 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
421 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
422 }
da6c0bc9 423 }
2be25cc9 424 (*cstream)<<"laserInfo"<<
da6c0bc9 425 "run="<<fRun<< // run number
426 "event="<<fEvent<< // event number
427 "time="<<fTime<< // time stamp of event
428 "trigger="<<fTrigger<< // trigger
429 "mag="<<fMagF<< // magnetic field
430 // Environment values
431 "press0="<<valuePressure0<<
432 "press1="<<valuePressure1<<
433 "pt0="<<ptrelative0<<
434 "pt1="<<ptrelative1<<
435 "temp0="<<temp0<<
436 "temp1="<<temp1<<
817766d5 437 "vecGoofie.="<<&vecGoofie<<
32100b1d 438 "vdcorr="<<vdcorr<<
da6c0bc9 439 //laser
dde68d36 440 "rejectA="<<isReject[0]<<
441 "rejectC="<<isReject[1]<<
2be25cc9 442 "laserA.="<<fLaser->fFitAside<<
443 "laserC.="<<fLaser->fFitCside<<
444 "laserAC.="<<fLaser->fFitACside<<
445 "trigger="<<event->GetFiredTriggerClasses()<<
da6c0bc9 446 "\n";
447 }
448 }
2be25cc9 449 //
2be25cc9 450 // fill histos
451 //
452 TVectorD vdriftA(5), vdriftC(5),vdriftAC(5);
453 vdriftA=*(fLaser->fFitAside);
454 vdriftC=*(fLaser->fFitCside);
455 vdriftAC=*(fLaser->fFitACside);
456 Int_t npointsA=0, npointsC=0;
457 Float_t chi2A=0, chi2C=0;
458 npointsA= TMath::Nint(vdriftA[3]);
459 chi2A= vdriftA[4];
460 npointsC= TMath::Nint(vdriftC[3]);
461 chi2C= vdriftC[4];
2be25cc9 462
2be25cc9 463 TTimeStamp tstamp(fTime);
464 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
465 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
dde68d36 466 Double_t driftA=0, driftC=0;
467 if (vdriftA[1]>1.-kMaxDeltaV) driftA = 1./vdriftA[1]-1.;
468 if (vdriftC[1]>1.-kMaxDeltaV) driftC = 1./vdriftC[1]-1.;
469 //
470 Double_t vecDriftLaserA[4]={fTime,(ptrelative0+ptrelative1)/2.0,driftA,event->GetRunNumber()};
471 Double_t vecDriftLaserC[4]={fTime,(ptrelative0+ptrelative1)/2.0,driftC,event->GetRunNumber()};
472 // Double_t vecDrift[4] ={fTime,(ptrelative0+ptrelative1)/2.0,1./((*(fLaser->fFitACside))[1])-1,event->GetRunNumber()};
473
474 for (Int_t icalib=0;icalib<3;icalib++){
475 if (icalib==0){ //z0 shift
476 vecDriftLaserA[2]=vdriftA[0]/250.;
477 vecDriftLaserC[2]=vdriftC[0]/250.;
2be25cc9 478 }
dde68d36 479 if (icalib==1){ //vdrel shift
480 vecDriftLaserA[2]=driftA;
481 vecDriftLaserC[2]=driftC;
2be25cc9 482 }
dde68d36 483 if (icalib==2){ //gy shift - full gy - full drift
484 vecDriftLaserA[2]=vdriftA[2]/250.;
485 vecDriftLaserC[2]=vdriftC[2]/250.;
2be25cc9 486 }
391ffdb2 487 //if (isReject[0]==0) fHistVdriftLaserA[icalib]->Fill(vecDriftLaserA);
488 //if (isReject[1]==0) fHistVdriftLaserC[icalib]->Fill(vecDriftLaserC);
489 fHistVdriftLaserA[icalib]->Fill(vecDriftLaserA);
490 fHistVdriftLaserC[icalib]->Fill(vecDriftLaserC);
2be25cc9 491 }
dde68d36 492
493// THnSparse* curHist=new THnSparseF("","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
494// TString shortName=curHist->ClassName();
495// shortName+="_MEAN_DRIFT_LASER_";
496// delete curHist;
497// curHist=NULL;
498// TString name="";
499
500// name=shortName;
501// name+=event->GetFiredTriggerClasses();
502// name.ToUpper();
503// curHist=(THnSparseF*)fArrayDz->FindObject(name);
504// if(!curHist){
505// curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
506// fArrayDz->AddLast(curHist);
507// }
508// curHist->Fill(vecDrift);
509
510// name=shortName;
511// name+="ALL";
512// name.ToUpper();
513// curHist=(THnSparseF*)fArrayDz->FindObject(name);
514// if(!curHist){
515// curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
516// fArrayDz->AddLast(curHist);
517// }
518// curHist->Fill(vecDrift);
2be25cc9 519}
c74c5f6c 520
2811495d 521void AliTPCcalibTime::ProcessCosmic(const AliESDEvent *const event){
6ba68367 522 //
523 // process Cosmic event - track matching A side C side
524 //
c74c5f6c 525 if (!event) {
526 Printf("ERROR: ESD not available");
527 return;
528 }
529 if (event->GetTimeStamp() == 0 ) {
530 Printf("no time stamp!");
531 return;
532 }
74235403 533
2be25cc9 534 //fd
c74c5f6c 535 // Find cosmic pairs
536 //
537 // Track0 is choosen in upper TPC part
538 // Track1 is choosen in lower TPC part
539 //
cc65e4f5 540 const Int_t kMinClustersCross =30;
541 const Int_t kMinClusters =80;
c74c5f6c 542 Int_t ntracks=event->GetNumberOfTracks();
543 if (ntracks==0) return;
74235403 544 if (ntracks > fCutTracks) return;
545
a390f11f 546 if (GetDebugLevel()>20) printf("Hallo world: Im here\n");
3e0d9569 547 AliESDfriend *esdFriend=(AliESDfriend*)(((AliESDEvent*)event)->FindListObject("AliESDfriend"));
c74c5f6c 548
549 TObjArray tpcSeeds(ntracks);
550 Double_t vtxx[3]={0,0,0};
551 Double_t svtxx[3]={0.000001,0.000001,100.};
552 AliESDVertex vtx(vtxx,svtxx);
553 //
554 // track loop
555 //
cc65e4f5 556 TArrayI clusterSideA(ntracks);
557 TArrayI clusterSideC(ntracks);
c74c5f6c 558 for (Int_t i=0;i<ntracks;++i) {
cc65e4f5 559 clusterSideA[i]=0;
560 clusterSideC[i]=0;
74235403 561 AliESDtrack *track = event->GetTrack(i);
562
563 const AliExternalTrackParam * trackIn = track->GetInnerParam();
564 const AliExternalTrackParam * trackOut = track->GetOuterParam();
565 if (!trackIn) continue;
566 if (!trackOut) continue;
567
6ba68367 568 AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
b1c27e4f 569 if (!friendTrack) continue;
3e55050f 570 if (friendTrack) ProcessSame(track,friendTrack,event);
6ba68367 571 if (friendTrack) ProcessAlignITS(track,friendTrack,event,esdFriend);
0dac7b3a 572 if (friendTrack) ProcessAlignTRD(track,friendTrack);
817766d5 573 if (friendTrack) ProcessAlignTOF(track,friendTrack);
74235403 574 TObject *calibObject;
575 AliTPCseed *seed = 0;
576 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
cc65e4f5 577 if (seed) {
578 tpcSeeds.AddAt(seed,i);
579 Int_t nA=0, nC=0;
580 for (Int_t irow=159;irow>0;irow--) {
581 AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
582 if (!cl) continue;
583 if ((cl->GetDetector()%36)<18) nA++;
584 if ((cl->GetDetector()%36)>=18) nC++;
585 }
586 clusterSideA[i]=nA;
587 clusterSideC[i]=nC;
588 }
c74c5f6c 589 }
c74c5f6c 590 if (ntracks<2) return;
591 //
592 // Find pairs
593 //
d3ce44cb 594
c74c5f6c 595 for (Int_t i=0;i<ntracks;++i) {
74235403 596 AliESDtrack *track0 = event->GetTrack(i);
c74c5f6c 597 // track0 - choosen upper part
598 if (!track0) continue;
599 if (!track0->GetOuterParam()) continue;
600 if (track0->GetOuterParam()->GetAlpha()<0) continue;
601 Double_t d1[3];
602 track0->GetDirection(d1);
603 for (Int_t j=0;j<ntracks;++j) {
604 if (i==j) continue;
74235403 605 AliESDtrack *track1 = event->GetTrack(j);
606 //track 1 lower part
607 if (!track1) continue;
608 if (!track1->GetOuterParam()) continue;
cc65e4f5 609 if (track0->GetTPCNcls()+ track1->GetTPCNcls()< kMinClusters) continue;
610 Int_t nAC = TMath::Max( TMath::Min(clusterSideA[i], clusterSideC[j]),
611 TMath::Min(clusterSideC[i], clusterSideA[j]));
612 if (nAC<kMinClustersCross) continue;
613 Int_t nA0=clusterSideA[i];
614 Int_t nC0=clusterSideC[i];
615 Int_t nA1=clusterSideA[j];
616 Int_t nC1=clusterSideC[j];
3e55050f 617 // if (track1->GetOuterParam()->GetAlpha()>0) continue;
74235403 618 //
619 Double_t d2[3];
620 track1->GetDirection(d2);
621
622 AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
623 AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
624 if (! seed0) continue;
625 if (! seed1) continue;
626 Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
627 Float_t dist0 = track0->GetLinearD(0,0);
628 Float_t dist1 = track1->GetLinearD(0,0);
629 //
630 // conservative cuts - convergence to be guarantied
631 // applying before track propagation
3e55050f 632 if (TMath::Abs(TMath::Abs(dist0)-TMath::Abs(dist1))>fCutMaxD) continue; // distance to the 0,0
633 if (TMath::Abs(dir)<TMath::Abs(fCutMinDir)) continue; // direction vector product
74235403 634 Float_t bz = AliTracker::GetBz();
635 Float_t dvertex0[2]; //distance to 0,0
636 Float_t dvertex1[2]; //distance to 0,0
637 track0->GetDZ(0,0,0,bz,dvertex0);
638 track1->GetDZ(0,0,0,bz,dvertex1);
639 if (TMath::Abs(dvertex0[1])>250) continue;
640 if (TMath::Abs(dvertex1[1])>250) continue;
641 //
642 //
643 //
644 Float_t dmax = TMath::Max(TMath::Abs(dist0),TMath::Abs(dist1));
645 AliExternalTrackParam param0(*track0);
646 AliExternalTrackParam param1(*track1);
647 //
648 // Propagate using Magnetic field and correct fo material budget
649 //
9963b5e2 650 AliTracker::PropagateTrackTo(&param0,dmax+1,TDatabasePDG::Instance()->GetParticle("e-")->Mass(),3,kTRUE);
651 AliTracker::PropagateTrackTo(&param1,dmax+1,TDatabasePDG::Instance()->GetParticle("e-")->Mass(),3,kTRUE);
74235403 652 //
653 // Propagate rest to the 0,0 DCA - z should be ignored
654 //
655 //Bool_t b0 = ;
656 param0.PropagateToDCA(&vtx,bz,1000);
657 //Bool_t b1 =
658 param1.PropagateToDCA(&vtx,bz,1000);
74235403 659 param0.GetDZ(0,0,0,bz,dvertex0);
660 param1.GetDZ(0,0,0,bz,dvertex1);
dde68d36 661 Double_t xyz0[3];
662 Double_t xyz1[3];
74235403 663 param0.GetXYZ(xyz0);
664 param1.GetXYZ(xyz1);
665 Bool_t isPair = IsPair(&param0,&param1);
d3ce44cb 666 Bool_t isCross = IsCross(track0, track1);
dde68d36 667 Bool_t isSame = IsSame(track0, track1);
d3ce44cb 668
dde68d36 669 THnSparse* hist=new THnSparseF("","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
670 TString shortName=hist->ClassName();
671 shortName+="_MEAN_VDRIFT_COSMICS_";
672 delete hist;
673 hist=NULL;
674
3e55050f 675 if((isSame) || (isCross && isPair)){
676 if (track0->GetTPCNcls()+ track1->GetTPCNcls()> 80) {
74235403 677 fDz = param0.GetZ() - param1.GetZ();
cc65e4f5 678 Double_t sign=(nA0>nA1)? 1:-1;
679 fDz*=sign;
74235403 680 TTimeStamp tstamp(fTime);
681 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
682 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
dde68d36 683 Double_t vecDrift[4]={fTime,(ptrelative0+ptrelative1)/2.0,fDz/500.0,event->GetRunNumber()};
684 THnSparse* curHist=NULL;
685 TString name="";
686
687 name=shortName;
688 name+=event->GetFiredTriggerClasses();
689 name.ToUpper();
690 curHist=(THnSparseF*)fArrayDz->FindObject(name);
74235403 691 if(!curHist){
dde68d36 692 curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
693 fArrayDz->AddLast(curHist);
74235403 694 }
dde68d36 695// curHist=(THnSparseF*)(fMapDz->GetValue(event->GetFiredTriggerClasses()));
696// if(!curHist){
697// curHist=new THnSparseF(event->GetFiredTriggerClasses(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
698// fMapDz->Add(new TObjString(event->GetFiredTriggerClasses()),curHist);
699// }
700 curHist->Fill(vecDrift);
74235403 701
dde68d36 702 name=shortName;
703 name+="ALL";
704 name.ToUpper();
705 curHist=(THnSparseF*)fArrayDz->FindObject(name);
74235403 706 if(!curHist){
dde68d36 707 curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
708 fArrayDz->AddLast(curHist);
74235403 709 }
dde68d36 710// curHist=(THnSparseF*)(fMapDz->GetValue("all"));
711// if(!curHist){
712// curHist=new THnSparseF("all","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
713// fMapDz->Add(new TObjString("all"),curHist);
714// }
715 curHist->Fill(vecDrift);
74235403 716 }
717 }
d3ce44cb 718 TTreeSRedirector *cstream = GetDebugStreamer();
719 if (fStreamLevel>0){
720 if (cstream){
721 (*cstream)<<"trackInfo"<<
722 "tr0.="<<track0<<
723 "tr1.="<<track1<<
724 "p0.="<<&param0<<
725 "p1.="<<&param1<<
cc65e4f5 726 "nAC="<<nAC<<
727 "nA0="<<nA0<<
728 "nA1="<<nA1<<
729 "nC0="<<nC0<<
730 "nC1="<<nC1<<
d3ce44cb 731 "isPair="<<isPair<<
732 "isCross="<<isCross<<
dde68d36 733 "isSame="<<isSame<<
d3ce44cb 734 "fDz="<<fDz<<
735 "fRun="<<fRun<<
736 "fTime="<<fTime<<
737 "\n";
738 }
739 }
c74c5f6c 740 } // end 2nd order loop
741 } // end 1st order loop
74235403 742
2be25cc9 743 if (fStreamLevel>0){
744 TTreeSRedirector *cstream = GetDebugStreamer();
745 if (cstream){
746 TTimeStamp tstamp(fTime);
747 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
748 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
749 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
750 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
751 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
752 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
32100b1d 753 Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
2be25cc9 754 TVectorD vecGoofie(20);
755 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
756 if (goofieArray){
757 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
758 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
759 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
760 }
761 }
762 (*cstream)<<"timeInfo"<<
763 "run="<<fRun<< // run number
764 "event="<<fEvent<< // event number
765 "time="<<fTime<< // time stamp of event
766 "trigger="<<fTrigger<< // trigger
767 "mag="<<fMagF<< // magnetic field
768 // Environment values
769 "press0="<<valuePressure0<<
770 "press1="<<valuePressure1<<
771 "pt0="<<ptrelative0<<
772 "pt1="<<ptrelative1<<
773 "temp0="<<temp0<<
774 "temp1="<<temp1<<
775 "vecGoofie.=<<"<<&vecGoofie<<
32100b1d 776 "vdcorr="<<vdcorr<<
2be25cc9 777 //
778 // accumulated values
779 //
780 "fDz="<<fDz<< //! current delta z
2be25cc9 781 "trigger="<<event->GetFiredTriggerClasses()<<
782 "\n";
783 }
784 }
a390f11f 785 if (GetDebugLevel()>20) printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
c74c5f6c 786}
787
2811495d 788void AliTPCcalibTime::ProcessBeam(const AliESDEvent *const /*event*/){
6ba68367 789 //
790 // Not special treatment yet - the same for cosmic and physic event
791 //
dde68d36 792}
793
6ba68367 794void AliTPCcalibTime::Analyze(){
795 //
796 // Special macro to analyze result of calibration and extract calibration entries
797 // Not yet ported to the Analyze function yet
798 //
799}
74235403 800
2811495d 801THnSparse* AliTPCcalibTime::GetHistoDrift(const char* name) const
802{
6ba68367 803 //
804 // Get histogram for given trigger mask
805 //
dde68d36 806 TIterator* iterator = fArrayDz->MakeIterator();
807 iterator->Reset();
808 TString newName=name;
809 newName.ToUpper();
810 THnSparse* newHist=new THnSparseF(newName,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
811 THnSparse* addHist=NULL;
812 while((addHist=(THnSparseF*)iterator->Next())){
813 if(!addHist) continue;
814 TString histName=addHist->GetName();
815 if(!histName.Contains(newName)) continue;
816 addHist->Print();
817 newHist->Add(addHist);
74235403 818 }
dde68d36 819 return newHist;
74235403 820}
821
2811495d 822TObjArray* AliTPCcalibTime::GetHistoDrift() const
823{
6ba68367 824 //
825 // return array of histograms
826 //
dde68d36 827 return fArrayDz;
74235403 828}
829
830TGraphErrors* AliTPCcalibTime::GetGraphDrift(const char* name){
6ba68367 831 //
832 // Make a drift velocity (delta Z) graph
833 //
dde68d36 834 THnSparse* histoDrift=GetHistoDrift(name);
835 TGraphErrors* graphDrift=NULL;
836 if(histoDrift){
837 graphDrift=FitSlices(histoDrift,2,0,400,100,0.05,0.95, kTRUE);
838 TString end=histoDrift->GetName();
839 Int_t pos=end.Index("_");
840 end=end(pos,end.Capacity()-pos);
841 TString graphName=graphDrift->ClassName();
842 graphName+=end;
843 graphName.ToUpper();
844 graphDrift->SetName(graphName);
74235403 845 }
846 return graphDrift;
847}
848
dde68d36 849TObjArray* AliTPCcalibTime::GetGraphDrift(){
6ba68367 850 //
851 // make a array of drift graphs
852 //
dde68d36 853 TObjArray* arrayGraphDrift=new TObjArray();
854 TIterator* iterator=fArrayDz->MakeIterator();
74235403 855 iterator->Reset();
dde68d36 856 THnSparse* addHist=NULL;
857 while((addHist=(THnSparseF*)iterator->Next())) arrayGraphDrift->AddLast(GetGraphDrift(addHist->GetName()));
858 return arrayGraphDrift;
74235403 859}
860
dde68d36 861AliSplineFit* AliTPCcalibTime::GetFitDrift(const char* name){
6ba68367 862 //
863 // Make a fit AliSplinefit of drift velocity
864 //
d3ce44cb 865 TGraph* graphDrift=GetGraphDrift(name);
dde68d36 866 AliSplineFit* fitDrift=NULL;
74235403 867 if(graphDrift && graphDrift->GetN()){
dde68d36 868 fitDrift=new AliSplineFit();
869 fitDrift->SetGraph(graphDrift);
870 fitDrift->SetMinPoints(graphDrift->GetN()+1);
871 fitDrift->InitKnots(graphDrift,2,0,0.001);
872 fitDrift->SplineFit(0);
873 TString end=graphDrift->GetName();
874 Int_t pos=end.Index("_");
875 end=end(pos,end.Capacity()-pos);
876 TString fitName=fitDrift->ClassName();
877 fitName+=end;
878 fitName.ToUpper();
879 //fitDrift->SetName(fitName);
74235403 880 delete graphDrift;
dde68d36 881 graphDrift=NULL;
74235403 882 }
883 return fitDrift;
884}
885
74235403 886
2811495d 887Long64_t AliTPCcalibTime::Merge(TCollection *const li) {
6ba68367 888 //
889 // Object specific merging procedure
890 //
c74c5f6c 891 TIterator* iter = li->MakeIterator();
892 AliTPCcalibTime* cal = 0;
893
894 while ((cal = (AliTPCcalibTime*)iter->Next())) {
895 if (!cal->InheritsFrom(AliTPCcalibTime::Class())) {
896 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
897 return -1;
898 }
2be25cc9 899 for (Int_t imeas=0; imeas<3; imeas++){
900 if (cal->GetHistVdriftLaserA(imeas) && cal->GetHistVdriftLaserA(imeas)){
901 fHistVdriftLaserA[imeas]->Add(cal->GetHistVdriftLaserA(imeas));
902 fHistVdriftLaserC[imeas]->Add(cal->GetHistVdriftLaserC(imeas));
903 }
904 }
9112627b 905 //
906 for (Int_t imeas=0; imeas<5; imeas++){
907 if (cal->GetResHistoTPCITS(imeas) && cal->GetResHistoTPCITS(imeas)){
908 fResHistoTPCITS[imeas]->Add(cal->fResHistoTPCITS[imeas]);
909 fResHistoTPCvertex[imeas]->Add(cal->fResHistoTPCvertex[imeas]);
910 fResHistoTPCTRD[imeas]->Add(cal->fResHistoTPCTRD[imeas]);
911 }
912 }
dde68d36 913 TObjArray* addArray=cal->GetHistoDrift();
914 if(!addArray) return 0;
915 TIterator* iterator = addArray->MakeIterator();
2be25cc9 916 iterator->Reset();
dde68d36 917 THnSparse* addHist=NULL;
918 while((addHist=(THnSparseF*)iterator->Next())){
919 if(!addHist) continue;
2be25cc9 920 addHist->Print();
dde68d36 921 THnSparse* localHist=(THnSparseF*)fArrayDz->FindObject(addHist->GetName());
2be25cc9 922 if(!localHist){
923 localHist=new THnSparseF(addHist->GetName(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
dde68d36 924 fArrayDz->AddLast(localHist);
2be25cc9 925 }
926 localHist->Add(addHist);
927 }
391ffdb2 928
d3ce44cb 929 for(Int_t i=0;i<10;i++) if (cal->GetCosmiMatchingHisto(i)) fCosmiMatchingHisto[i]->Add(cal->GetCosmiMatchingHisto(i));
3e55050f 930 //
931 // Merge alignment
932 //
933 for (Int_t itype=0; itype<3; itype++){
934 //
935 //
936 TObjArray *arr0= 0;
937 TObjArray *arr1= 0;
938 if (itype==0) {arr0=fAlignITSTPC; arr1=cal->fAlignITSTPC;}
939 if (itype==1) {arr0=fAlignTRDTPC; arr1=cal->fAlignTRDTPC;}
940 if (itype==2) {arr0=fAlignTOFTPC; arr1=cal->fAlignTOFTPC;}
941 if (!arr1) continue;
942 if (!arr0) arr0=new TObjArray(arr1->GetEntriesFast());
943 if (arr1->GetEntriesFast()>arr0->GetEntriesFast()){
944 arr0->Expand(arr1->GetEntriesFast());
945 }
946 for (Int_t i=0;i<arr1->GetEntriesFast(); i++){
947 AliRelAlignerKalman *kalman1 = (AliRelAlignerKalman *)arr1->UncheckedAt(i);
948 AliRelAlignerKalman *kalman0 = (AliRelAlignerKalman *)arr0->UncheckedAt(i);
949 if (!kalman1) continue;
950 if (!kalman0) {arr0->AddAt(new AliRelAlignerKalman(*kalman1),i); continue;}
951 kalman0->SetRejectOutliers(kFALSE);
952 kalman0->Merge(kalman1);
953 }
954 }
955
c74c5f6c 956 }
c74c5f6c 957 return 0;
c74c5f6c 958}
959
c74c5f6c 960Bool_t AliTPCcalibTime::IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1){
c74c5f6c 961 /*
962 // 0. Same direction - OPOSITE - cutDir +cutT
963 TCut cutDir("cutDir","dir<-0.99")
964 // 1.
965 TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03")
966 //
967 // 2. The same rphi
968 TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5")
969 //
c74c5f6c 970 TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
971 // 1/Pt diff cut
972 */
973 const Double_t *p0 = tr0->GetParameter();
974 const Double_t *p1 = tr1->GetParameter();
74235403 975 fCosmiMatchingHisto[0]->Fill(p0[0]+p1[0]);
976 fCosmiMatchingHisto[1]->Fill(p0[1]-p1[1]);
d3ce44cb 977 fCosmiMatchingHisto[2]->Fill(tr1->GetAlpha()-tr0->GetAlpha()+TMath::Pi());
74235403 978 fCosmiMatchingHisto[3]->Fill(p0[3]+p1[3]);
979 fCosmiMatchingHisto[4]->Fill(p0[4]+p1[4]);
980
c74c5f6c 981 if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
982 if (TMath::Abs(p0[0]+p1[0])>fCutMaxD) return kFALSE;
da6c0bc9 983 if (TMath::Abs(p0[1]-p1[1])>fCutMaxDz) return kFALSE;
c74c5f6c 984 Double_t d0[3], d1[3];
985 tr0->GetDirection(d0);
986 tr1->GetDirection(d1);
987 if (d0[0]*d1[0] + d0[1]*d1[1] + d0[2]*d1[2] >fCutMinDir) return kFALSE;
74235403 988
d3ce44cb 989 fCosmiMatchingHisto[5]->Fill(p0[0]+p1[0]);
990 fCosmiMatchingHisto[6]->Fill(p0[1]-p1[1]);
991 fCosmiMatchingHisto[7]->Fill(tr1->GetAlpha()-tr0->GetAlpha()+TMath::Pi());
992 fCosmiMatchingHisto[8]->Fill(p0[3]+p1[3]);
993 fCosmiMatchingHisto[9]->Fill(p0[4]+p1[4]);
994
c74c5f6c 995 return kTRUE;
996}
2811495d 997Bool_t AliTPCcalibTime::IsCross(AliESDtrack *const tr0, AliESDtrack *const tr1){
cc65e4f5 998 //
999 // check if the cosmic pair of tracks crossed A/C side
1000 //
1001 Bool_t result= tr0->GetOuterParam()->GetZ()*tr1->GetOuterParam()->GetZ()<0;
1002 if (result==kFALSE) return result;
1003 result=kTRUE;
1004 return result;
d3ce44cb 1005}
1006
2811495d 1007Bool_t AliTPCcalibTime::IsSame(AliESDtrack *const tr0, AliESDtrack *const tr1){
3e55050f 1008 //
1009 // track crossing the CE
1010 // 0. minimal number of clusters
1011 // 1. Same sector +-1
1012 // 2. Inner and outer track param on opposite side
1013 // 3. Outer and inner track parameter close each to other
1014 // 3.
1015 Bool_t result=kTRUE;
1016 //
1017 // inner and outer on opposite sides in z
1018 //
1019 const Int_t knclCut0 = 30;
1020 const Double_t kalphaCut = 0.4;
1021 //
1022 // 0. minimal number of clusters
1023 //
1024 if (tr0->GetTPCNcls()<knclCut0) return kFALSE;
1025 if (tr1->GetTPCNcls()<knclCut0) return kFALSE;
1026 //
1027 // 1. alpha cut - sector+-1
1028 //
1029 if (TMath::Abs(tr0->GetOuterParam()->GetAlpha()-tr1->GetOuterParam()->GetAlpha())>kalphaCut) return kFALSE;
1030 //
1031 // 2. Z crossing
1032 //
1033 if (tr0->GetOuterParam()->GetZ()*tr0->GetInnerParam()->GetZ()>0) result&=kFALSE;
1034 if (tr1->GetOuterParam()->GetZ()*tr1->GetInnerParam()->GetZ()>0) result&=kFALSE;
1035 if (result==kFALSE){
1036 return result;
1037 }
1038 //
1039 //
1040 const Double_t *p0I = tr0->GetInnerParam()->GetParameter();
1041 const Double_t *p1I = tr1->GetInnerParam()->GetParameter();
1042 const Double_t *p0O = tr0->GetOuterParam()->GetParameter();
1043 const Double_t *p1O = tr1->GetOuterParam()->GetParameter();
1044 //
1045 if (TMath::Abs(p0I[0]-p1I[0])>fCutMaxD) result&=kFALSE;
1046 if (TMath::Abs(p0I[1]-p1I[1])>fCutMaxDz) result&=kFALSE;
1047 if (TMath::Abs(p0I[2]-p1I[2])>fCutTheta) result&=kFALSE;
1048 if (TMath::Abs(p0I[3]-p1I[3])>fCutTheta) result&=kFALSE;
1049 if (TMath::Abs(p0O[0]-p1O[0])>fCutMaxD) result&=kFALSE;
1050 if (TMath::Abs(p0O[1]-p1O[1])>fCutMaxDz) result&=kFALSE;
1051 if (TMath::Abs(p0O[2]-p1O[2])>fCutTheta) result&=kFALSE;
1052 if (TMath::Abs(p0O[3]-p1O[3])>fCutTheta) result&=kFALSE;
1053 if (result==kTRUE){
1054 result=kTRUE; // just to put break point here
1055 }
1056 return result;
dde68d36 1057}
1058
31aa7c5c 1059
2811495d 1060void AliTPCcalibTime::ProcessSame(AliESDtrack *const track, AliESDfriendTrack *const friendTrack, const AliESDEvent *const event){
3e55050f 1061 //
1062 // Process TPC tracks crossing CE
1063 //
1064 // 0. Select only track crossing the CE
1065 // 1. Cut on the track length
1066 // 2. Refit the terack on A and C side separatelly
1067 // 3. Fill time histograms
1068 const Int_t kMinNcl=100;
1069 const Int_t kMinNclS=25; // minimul number of clusters on the sides
1070 if (!friendTrack->GetTPCOut()) return;
1071 //
1072 // 0. Select only track crossing the CE
1073 //
1074 if (track->GetInnerParam()->GetZ()*friendTrack->GetTPCOut()->GetZ()>0) return;
1075 //
1076 // 1. cut on track length
1077 //
1078 if (track->GetTPCNcls()<kMinNcl) return;
1079 //
1080 // 2. Refit track sepparatel on A and C side
1081 //
1082 TObject *calibObject;
1083 AliTPCseed *seed = 0;
1084 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
1085 if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
1086 }
1087 if (!seed) return;
1088 //
1089 AliExternalTrackParam trackIn(*track->GetInnerParam());
1090 AliExternalTrackParam trackOut(*track->GetOuterParam());
1091 Double_t cov[3]={0.01,0.,0.01}; //use the same errors
1092 Double_t xyz[3]={0,0.,0.0};
1093 Double_t bz =0;
1094 Int_t nclIn=0,nclOut=0;
1095 trackIn.ResetCovariance(30.);
1096 trackOut.ResetCovariance(30.);
1097 //
1098 //2.a Refit inner
1099 //
1100 for (Int_t irow=0;irow<159;irow++) {
1101 AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
1102 if (!cl) continue;
1103 if (cl->GetX()<80) continue;
1104 if (track->GetInnerParam()->GetZ()<0 &&(cl->GetDetector()%36)<18) break;
1105 if (track->GetInnerParam()->GetZ()>0 &&(cl->GetDetector()%36)>=18) break;
1106 Int_t sector = cl->GetDetector();
1107 Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
1108 if (TMath::Abs(dalpha)>0.01){
1109 if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
1110 }
1111 Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
1112 trackIn.GetXYZ(xyz);
1113 bz = AliTracker::GetBz(xyz);
1114 if (!trackIn.PropagateTo(r[0],bz)) break;
1115 nclIn++;
1116 trackIn.Update(&r[1],cov);
1117 }
1118 //
1119 //2.b Refit outer
1120 //
1121 for (Int_t irow=159;irow>0;irow--) {
1122 AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
1123 if (!cl) continue;
1124 if (cl->GetX()<80) continue;
1125 if (cl->GetZ()*track->GetOuterParam()->GetZ()<0) break;
1126 if (friendTrack->GetTPCOut()->GetZ()<0 &&(cl->GetDetector()%36)<18) break;
1127 if (friendTrack->GetTPCOut()->GetZ()>0 &&(cl->GetDetector()%36)>=18) break;
1128 Int_t sector = cl->GetDetector();
1129 Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();
1130 if (TMath::Abs(dalpha)>0.01){
1131 if (!trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
1132 }
1133 Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
1134 trackOut.GetXYZ(xyz);
1135 bz = AliTracker::GetBz(xyz);
1136 if (!trackOut.PropagateTo(r[0],bz)) break;
1137 nclOut++;
1138 trackOut.Update(&r[1],cov);
1139 }
1140 trackOut.Rotate(trackIn.GetAlpha());
1141 Double_t meanX = (trackIn.GetX()+trackOut.GetX())*0.5;
1142 trackIn.PropagateTo(meanX,bz);
1143 trackOut.PropagateTo(meanX,bz);
1144 TTreeSRedirector *cstream = GetDebugStreamer();
1145 if (cstream){
1146 TVectorD gxyz(3);
1147 trackIn.GetXYZ(gxyz.GetMatrixArray());
1148 TTimeStamp tstamp(fTime);
1149 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
1150 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
32100b1d 1151 Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
3e55050f 1152 (*cstream)<<"tpctpc"<<
1153 "run="<<fRun<< // run number
1154 "event="<<fEvent<< // event number
1155 "time="<<fTime<< // time stamp of event
1156 "trigger="<<fTrigger<< // trigger
1157 "mag="<<fMagF<< // magnetic field
1158 "ptrel0.="<<ptrelative0<<
1159 "ptrel1.="<<ptrelative1<<
32100b1d 1160 "vdcorr="<<vdcorr<< // drift correction applied
3e55050f 1161 //
1162 "xyz.="<<&gxyz<< // global position
1163 "tIn.="<<&trackIn<< // refitterd track in
1164 "tOut.="<<&trackOut<< // refitter track out
1165 "nclIn="<<nclIn<< //
1166 "nclOut="<<nclOut<< //
1167 "\n";
1168 }
1169 //
1170 // 3. Fill time histograms
1171 // Debug stremaer expression
1172 // chainTPCTPC->Draw("(tIn.fP[1]-tOut.fP[1])*sign(-tIn.fP[3]):tIn.fP[3]","min(nclIn,nclOut)>30","")
1173 if (TMath::Min(nclIn,nclOut)>kMinNclS){
1174 fDz = trackOut.GetZ()-trackIn.GetZ();
1175 if (trackOut.GetTgl()<0) fDz*=-1.;
1176 TTimeStamp tstamp(fTime);
1177 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
1178 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
1179 Double_t vecDrift[4]={fTime,(ptrelative0+ptrelative1)/2.0,fDz/500.0,event->GetRunNumber()};
1180 //
1181 // fill histograms per trigger class and itegrated
1182 //
1183 THnSparse* curHist=NULL;
1184 for (Int_t itype=0; itype<2; itype++){
1185 TString name="MEAN_VDRIFT_CROSS_";
1186 if (itype==0){
1187 name+=event->GetFiredTriggerClasses();
1188 name.ToUpper();
1189 }else{
1190 name+="ALL";
1191 }
1192 curHist=(THnSparseF*)fArrayDz->FindObject(name);
1193 if(!curHist){
1194 curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
1195 fArrayDz->AddLast(curHist);
1196 }
1197 curHist->Fill(vecDrift);
1198 }
1199 }
1200
1201}
817766d5 1202
2811495d 1203void AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTrack *const friendTrack, const AliESDEvent *const event, AliESDfriend *const esdFriend){
817766d5 1204 //
3e55050f 1205 // Process track - Update TPC-ITS alignment
1206 // Updates:
1207 // 0. Apply standartd cuts
1208 // 1. Recalucluate the current statistic median/RMS
1209 // 2. Apply median+-rms cut
1210 // 3. Update kalman filter
1211 //
1212 const Int_t kMinTPC = 80; // minimal number of TPC cluster
1213 const Int_t kMinITS = 3; // minimal number of ITS cluster
1214 const Double_t kMinZ = 10; // maximal dz distance
b96c3aef 1215 const Double_t kMaxDy = 2.; // maximal dy distance
1216 const Double_t kMaxAngle= 0.015; // maximal angular distance
3e55050f 1217 const Double_t kSigmaCut= 5; // maximal sigma distance to median
1218 const Double_t kVdErr = 0.1; // initial uncertainty of the vd correction
1219 const Double_t kVdYErr = 0.05; // initial uncertainty of the vd correction
1220 const Double_t kOutCut = 1.0; // outlyer cut in AliRelAlgnmentKalman
b96c3aef 1221 const Double_t kMinPt = 0.3; // minimal pt
be67055b 1222 const Int_t kN=50; // deepnes of history
3e55050f 1223 static Int_t kglast=0;
1224 static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
1225 /*
1226 0. Standrd cuts:
1227 TCut cut="abs(pTPC.fP[2]-pITS.fP[2])<0.01&&abs(pTPC.fP[3]-pITS.fP[3])<0.01&&abs(pTPC.fP[2]-pITS.fP[2])<1";
1228 */
817766d5 1229 //
3e55050f 1230 // 0. Apply standard cuts
817766d5 1231 //
1232 Int_t dummycl[1000];
817766d5 1233 if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
b842d904 1234 if (track->GetITSclusters(dummycl)<kMinITS) return; // minimal amount of clusters
b96c3aef 1235 if (!track->IsOn(AliESDtrack::kTPCrefit)) return;
3e55050f 1236 if (!friendTrack->GetITSOut()) return;
817766d5 1237 if (!track->GetInnerParam()) return;
1238 if (!track->GetOuterParam()) return;
b842d904 1239 if (track->GetInnerParam()->Pt()<kMinPt) return;
817766d5 1240 // exclude crossing track
1241 if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
1242 if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ) return;
b96c3aef 1243 if (track->GetInnerParam()->GetX()>90) return;
817766d5 1244 //
1245 AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(track->GetInnerParam()));
b842d904 1246 AliExternalTrackParam pITS(*(friendTrack->GetITSOut())); // ITS standalone if possible
1247 AliExternalTrackParam pITS2(*(friendTrack->GetITSOut())); //TPC-ITS track
1248 pITS2.Rotate(pTPC.GetAlpha());
be67055b 1249 // pITS2.PropagateTo(pTPC.GetX(),fMagF);
1250 AliTracker::PropagateTrackToBxByBz(&pITS2,pTPC.GetX(),0.1,0.1,kFALSE);
1251
b842d904 1252 AliESDfriendTrack *itsfriendTrack=0;
1253 //
1254 // try to find standalone ITS track corresponing to the TPC if possible
1255 //
1256 Bool_t hasAlone=kFALSE;
1257 Int_t ntracks=event->GetNumberOfTracks();
1258 for (Int_t i=0; i<ntracks; i++){
1259 AliESDtrack *trackS = event->GetTrack(i);
1260 if (trackS->GetTPCNcls()>0) continue; //continue if has TPC info
6ba68367 1261 itsfriendTrack = esdFriend->GetTrack(i);
b842d904 1262 if (!itsfriendTrack) continue;
1263 if (!itsfriendTrack->GetITSOut()) continue;
1264 if (TMath::Abs(pITS2.GetTgl()-itsfriendTrack->GetITSOut()->GetTgl())> kMaxAngle) continue;
1265 pITS=(*(itsfriendTrack->GetITSOut()));
1266 //
1267 pITS.Rotate(pTPC.GetAlpha());
be67055b 1268 //pITS.PropagateTo(pTPC.GetX(),fMagF);
1269 AliTracker::PropagateTrackToBxByBz(&pITS,pTPC.GetX(),0.1,0.1,kFALSE);
b842d904 1270 if (TMath::Abs(pITS2.GetY()-pITS.GetY())> kMaxDy) continue;
1271 hasAlone=kTRUE;
1272 }
1273 if (!hasAlone) pITS=pITS2;
1274 //
3e55050f 1275 if (TMath::Abs(pITS.GetY()-pTPC.GetY()) >kMaxDy) return;
1276 if (TMath::Abs(pITS.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
1277 if (TMath::Abs(pITS.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
1278 //
1279 // 1. Update median and RMS info
1280 //
b322e06a 1281 TVectorD vecDelta(5),vecMedian(5), vecRMS(5);
3e55050f 1282 TVectorD vecDeltaN(5);
1283 Double_t sign=(pITS.GetParameter()[1]>0)? 1.:-1.;
1284 vecDelta[4]=0;
1285 for (Int_t i=0;i<4;i++){
1286 vecDelta[i]=(pITS.GetParameter()[i]-pTPC.GetParameter()[i])*sign;
1287 kgdP[i][kglast%kN]=vecDelta[i];
1288 }
1289 kglast=(kglast+1);
1290 Int_t entries=(kglast<kN)?kglast:kN;
1291 for (Int_t i=0;i<4;i++){
1292 vecMedian[i] = TMath::Median(entries,kgdP[i]);
1293 vecRMS[i] = TMath::RMS(entries,kgdP[i]);
1294 vecDeltaN[i] = 0;
1295 if (vecRMS[i]>0.){
1296 vecDeltaN[i] = (vecDelta[i]-vecMedian[i])/vecRMS[i];
1297 vecDeltaN[4]+= TMath::Abs(vecDeltaN[i]); //sum of abs residuals
1298 }
1299 }
817766d5 1300 //
3e55050f 1301 // 2. Apply median+-rms cut
817766d5 1302 //
3e55050f 1303 if (kglast<3) return; //median and RMS to be defined
1304 if ( vecDeltaN[4]/4.>kSigmaCut) return;
1305 //
1306 // 3. Update alignment
817766d5 1307 //
1308 Int_t htime = fTime/3600; //time in hours
1309 if (fAlignITSTPC->GetEntries()<htime){
1310 fAlignITSTPC->Expand(htime*2+20);
1311 }
1312 AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignITSTPC->At(htime);
1313 if (!align){
3e55050f 1314 // make Alignment object if doesn't exist
817766d5 1315 align=new AliRelAlignerKalman();
3e55050f 1316 align->SetRunNumber(fRun);
1317 (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
1318 (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
1319 align->SetOutRejSigma(kOutCut+kOutCut*kN);
1320 align->SetRejectOutliers(kFALSE);
1321
1322 align->SetTPCvd(AliTPCcalibDB::Instance()->GetParameters()->GetDriftV()/1000000.);
817766d5 1323 align->SetMagField(fMagF);
1324 fAlignITSTPC->AddAt(align,htime);
1325 }
817766d5 1326 align->AddTrackParams(&pITS,&pTPC);
1327 align->SetTimeStamp(fTime);
3e55050f 1328 align->SetRunNumber(fRun );
391ffdb2 1329 Float_t dca[2],cov[3];
1330 track->GetImpactParameters(dca,cov);
882f5c06 1331 if (TMath::Abs(dca[0])<kMaxDy){
391ffdb2 1332 FillResHistoTPCITS(&pTPC,&pITS);
1333 FillResHistoTPC(track);
1334 }
3e55050f 1335 //
1336 Int_t nupdates=align->GetNUpdates();
1337 align->SetOutRejSigma(kOutCut+kOutCut*kN/Double_t(nupdates));
1338 align->SetRejectOutliers(kFALSE);
817766d5 1339 TTreeSRedirector *cstream = GetDebugStreamer();
1340 if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
1341 TTimeStamp tstamp(fTime);
1342 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
1343 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
1344 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
1345 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
1346 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
1347 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
1348 TVectorD vecGoofie(20);
1349 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
1350 if (goofieArray){
1351 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
1352 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
1353 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
1354 }
1355 }
1356 TVectorD gpTPC(3), gdTPC(3);
1357 TVectorD gpITS(3), gdITS(3);
1358 pTPC.GetXYZ(gpTPC.GetMatrixArray());
1359 pTPC.GetDirection(gdTPC.GetMatrixArray());
1360 pITS.GetXYZ(gpITS.GetMatrixArray());
1361 pITS.GetDirection(gdITS.GetMatrixArray());
32100b1d 1362 Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
817766d5 1363 (*cstream)<<"itstpc"<<
1364 "run="<<fRun<< // run number
1365 "event="<<fEvent<< // event number
1366 "time="<<fTime<< // time stamp of event
1367 "trigger="<<fTrigger<< // trigger
1368 "mag="<<fMagF<< // magnetic field
1369 // Environment values
1370 "press0="<<valuePressure0<<
1371 "press1="<<valuePressure1<<
1372 "pt0="<<ptrelative0<<
1373 "pt1="<<ptrelative1<<
1374 "temp0="<<temp0<<
1375 "temp1="<<temp1<<
1376 "vecGoofie.="<<&vecGoofie<<
32100b1d 1377 "vdcorr="<<vdcorr<< // drift correction applied
817766d5 1378 //
b842d904 1379 "hasAlone="<<hasAlone<< // has ITS standalone ?
1380 "track.="<<track<< // track info
3e55050f 1381 "nmed="<<kglast<< // number of entries to define median and RMS
1382 "vMed.="<<&vecMedian<< // median of deltas
1383 "vRMS.="<<&vecRMS<< // rms of deltas
1384 "vDelta.="<<&vecDelta<< // delta in respect to median
1385 "vDeltaN.="<<&vecDeltaN<< // normalized delta in respect to median
1386 "t.="<<track<< // ful track - find proper cuts
1387 "a.="<<align<< // current alignment
b842d904 1388 "pITS.="<<&pITS<< // track param ITS
1389 "pITS2.="<<&pITS2<< // track param ITS+TPC
3e55050f 1390 "pTPC.="<<&pTPC<< // track param TPC
1391 "gpTPC.="<<&gpTPC<< // global position TPC
1392 "gdTPC.="<<&gdTPC<< // global direction TPC
1393 "gpITS.="<<&gpITS<< // global position ITS
1394 "gdITS.="<<&gdITS<< // global position ITS
817766d5 1395 "\n";
1396 }
817766d5 1397}
3e55050f 1398
1399
1400
1401
2811495d 1402void AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTrack *const friendTrack){
0dac7b3a 1403 //
3e55050f 1404 // Process track - Update TPC-TRD alignment
1405 // Updates:
1406 // 0. Apply standartd cuts
1407 // 1. Recalucluate the current statistic median/RMS
1408 // 2. Apply median+-rms cut
1409 // 3. Update kalman filter
1410 //
1411 const Int_t kMinTPC = 80; // minimal number of TPC cluster
1412 const Int_t kMinTRD = 50; // minimal number of TRD cluster
1413 const Double_t kMinZ = 20; // maximal dz distance
9112627b 1414 const Double_t kMaxDy = 5.; // maximal dy distance
1415 const Double_t kMaxAngle= 0.1; // maximal angular distance
1416 const Double_t kSigmaCut= 10; // maximal sigma distance to median
3e55050f 1417 const Double_t kVdErr = 0.1; // initial uncertainty of the vd correction
1418 const Double_t kVdYErr = 0.05; // initial uncertainty of the vd correction
1419 const Double_t kOutCut = 1.0; // outlyer cut in AliRelAlgnmentKalman
be67055b 1420 const Double_t kRefX = 275; // reference X
1421 const Int_t kN=50; // deepnes of history
3e55050f 1422 static Int_t kglast=0;
1423 static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
0dac7b3a 1424 //
3e55050f 1425 // 0. Apply standard cuts
0dac7b3a 1426 //
1427 Int_t dummycl[1000];
1428 if (track->GetTRDclusters(dummycl)<kMinTRD) return; // minimal amount of clusters
1429 if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
3e55050f 1430 if (!friendTrack->GetTRDIn()) return;
be67055b 1431 if (!track->IsOn(AliESDtrack::kTRDrefit)) return;
1432 if (!track->IsOn(AliESDtrack::kTRDout)) return;
0dac7b3a 1433 if (!track->GetInnerParam()) return;
be67055b 1434 if (!friendTrack->GetTPCOut()) return;
0dac7b3a 1435 // exclude crossing track
be67055b 1436 if (friendTrack->GetTPCOut()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
3e55050f 1437 if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ) return;
0dac7b3a 1438 //
be67055b 1439 AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(friendTrack->GetTPCOut()));
1440 AliTracker::PropagateTrackToBxByBz(&pTPC,kRefX,0.1,0.1,kFALSE);
0dac7b3a 1441 AliExternalTrackParam pTRD(*(friendTrack->GetTRDIn()));
3e55050f 1442 pTRD.Rotate(pTPC.GetAlpha());
be67055b 1443 // pTRD.PropagateTo(pTPC.GetX(),fMagF);
1444 AliTracker::PropagateTrackToBxByBz(&pTRD,pTPC.GetX(),0.1,0.1,kFALSE);
1445
3e55050f 1446 ((Double_t*)pTRD.GetCovariance())[2]+=3.*3.; // increas sys errors
1447 ((Double_t*)pTRD.GetCovariance())[9]+=0.1*0.1; // increse sys errors
1448
1449 if (TMath::Abs(pTRD.GetY()-pTPC.GetY()) >kMaxDy) return;
1450 if (TMath::Abs(pTRD.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
9112627b 1451 // if (TMath::Abs(pTRD.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
0dac7b3a 1452 //
3e55050f 1453 // 1. Update median and RMS info
0dac7b3a 1454 //
b322e06a 1455 TVectorD vecDelta(5),vecMedian(5), vecRMS(5);
3e55050f 1456 TVectorD vecDeltaN(5);
1457 Double_t sign=(pTRD.GetParameter()[1]>0)? 1.:-1.;
1458 vecDelta[4]=0;
1459 for (Int_t i=0;i<4;i++){
1460 vecDelta[i]=(pTRD.GetParameter()[i]-pTPC.GetParameter()[i])*sign;
1461 kgdP[i][kglast%kN]=vecDelta[i];
1462 }
1463 kglast=(kglast+1);
1464 Int_t entries=(kglast<kN)?kglast:kN;
1465 for (Int_t i=0;i<4;i++){
1466 vecMedian[i] = TMath::Median(entries,kgdP[i]);
be67055b 1467
3e55050f 1468 vecRMS[i] = TMath::RMS(entries,kgdP[i]);
1469 vecDeltaN[i] = 0;
1470 if (vecRMS[i]>0.){
1471 vecDeltaN[i] = (vecDelta[i]-vecMedian[i])/vecRMS[i];
1472 vecDeltaN[4]+= TMath::Abs(vecDeltaN[i]); //sum of abs residuals
1473 }
1474 }
1475 //
1476 // 2. Apply median+-rms cut
1477 //
1478 if (kglast<3) return; //median and RMS to be defined
1479 if ( vecDeltaN[4]/4.>kSigmaCut) return;
1480 //
1481 // 3. Update alignment
0dac7b3a 1482 //
1483 Int_t htime = fTime/3600; //time in hours
1484 if (fAlignTRDTPC->GetEntries()<htime){
1485 fAlignTRDTPC->Expand(htime*2+20);
1486 }
1487 AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTRDTPC->At(htime);
1488 if (!align){
3e55050f 1489 // make Alignment object if doesn't exist
0dac7b3a 1490 align=new AliRelAlignerKalman();
3e55050f 1491 align->SetRunNumber(fRun);
1492 (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
1493 (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
1494 align->SetOutRejSigma(kOutCut+kOutCut*kN);
1495 align->SetRejectOutliers(kFALSE);
1496 align->SetTPCvd(AliTPCcalibDB::Instance()->GetParameters()->GetDriftV()/1000000.);
0dac7b3a 1497 align->SetMagField(fMagF);
1498 fAlignTRDTPC->AddAt(align,htime);
1499 }
0dac7b3a 1500 align->AddTrackParams(&pTRD,&pTPC);
1501 align->SetTimeStamp(fTime);
3e55050f 1502 align->SetRunNumber(fRun );
9112627b 1503 FillResHistoTPCTRD(&pTPC,&pTRD);
3e55050f 1504 //
1505 Int_t nupdates=align->GetNUpdates();
1506 align->SetOutRejSigma(kOutCut+kOutCut*kN/Double_t(nupdates));
1507 align->SetRejectOutliers(kFALSE);
0dac7b3a 1508 TTreeSRedirector *cstream = GetDebugStreamer();
1509 if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
1510 TTimeStamp tstamp(fTime);
1511 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
1512 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
1513 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
1514 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
1515 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
1516 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
1517 TVectorD vecGoofie(20);
1518 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
1519 if (goofieArray){
1520 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
1521 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
1522 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
1523 }
1524 }
1525 TVectorD gpTPC(3), gdTPC(3);
1526 TVectorD gpTRD(3), gdTRD(3);
1527 pTPC.GetXYZ(gpTPC.GetMatrixArray());
1528 pTPC.GetDirection(gdTPC.GetMatrixArray());
1529 pTRD.GetXYZ(gpTRD.GetMatrixArray());
1530 pTRD.GetDirection(gdTRD.GetMatrixArray());
32100b1d 1531 Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
3e55050f 1532 (*cstream)<<"trdtpc"<<
0dac7b3a 1533 "run="<<fRun<< // run number
1534 "event="<<fEvent<< // event number
1535 "time="<<fTime<< // time stamp of event
1536 "trigger="<<fTrigger<< // trigger
1537 "mag="<<fMagF<< // magnetic field
1538 // Environment values
1539 "press0="<<valuePressure0<<
1540 "press1="<<valuePressure1<<
1541 "pt0="<<ptrelative0<<
1542 "pt1="<<ptrelative1<<
1543 "temp0="<<temp0<<
1544 "temp1="<<temp1<<
1545 "vecGoofie.="<<&vecGoofie<<
32100b1d 1546 "vdcorr="<<vdcorr<< // drift correction applied
0dac7b3a 1547 //
3e55050f 1548 "nmed="<<kglast<< // number of entries to define median and RMS
1549 "vMed.="<<&vecMedian<< // median of deltas
1550 "vRMS.="<<&vecRMS<< // rms of deltas
1551 "vDelta.="<<&vecDelta<< // delta in respect to median
1552 "vDeltaN.="<<&vecDeltaN<< // normalized delta in respect to median
1553 "t.="<<track<< // ful track - find proper cuts
1554 "a.="<<align<< // current alignment
1555 "pTRD.="<<&pTRD<< // track param TRD
1556 "pTPC.="<<&pTPC<< // track param TPC
1557 "gpTPC.="<<&gpTPC<< // global position TPC
1558 "gdTPC.="<<&gdTPC<< // global direction TPC
1559 "gpTRD.="<<&gpTRD<< // global position TRD
1560 "gdTRD.="<<&gdTRD<< // global position TRD
0dac7b3a 1561 "\n";
1562 }
0dac7b3a 1563}
817766d5 1564
1565
2811495d 1566void AliTPCcalibTime::ProcessAlignTOF(AliESDtrack *const track, AliESDfriendTrack *const friendTrack){
817766d5 1567 //
817766d5 1568 //
3e55050f 1569 // Process track - Update TPC-TOF alignment
1570 // Updates:
1571 // -1. Make a TOF "track"
1572 // 0. Apply standartd cuts
1573 // 1. Recalucluate the current statistic median/RMS
1574 // 2. Apply median+-rms cut
1575 // 3. Update kalman filter
1576 //
1577 const Int_t kMinTPC = 80; // minimal number of TPC cluster
b842d904 1578 // const Double_t kMinZ = 10; // maximal dz distance
3e55050f 1579 const Double_t kMaxDy = 5.; // maximal dy distance
b96c3aef 1580 const Double_t kMaxAngle= 0.015; // maximal angular distance
3e55050f 1581 const Double_t kSigmaCut= 5; // maximal sigma distance to median
1582 const Double_t kVdErr = 0.1; // initial uncertainty of the vd correction
1583 const Double_t kVdYErr = 0.05; // initial uncertainty of the vd correction
1584
1585 const Double_t kOutCut = 1.0; // outlyer cut in AliRelAlgnmentKalman
be67055b 1586 const Int_t kN=50; // deepnes of history
3e55050f 1587 static Int_t kglast=0;
1588 static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
817766d5 1589 //
3e55050f 1590 // -1. Make a TOF track-
1591 // Clusters are not in friends - use alingment points
1592 //
1593 if (track->GetTOFsignal()<=0) return;
1594 if (!friendTrack->GetTPCOut()) return;
1595 if (!track->GetInnerParam()) return;
be67055b 1596 if (!friendTrack->GetTPCOut()) return;
817766d5 1597 const AliTrackPointArray *points=friendTrack->GetTrackPointArray();
1598 if (!points) return;
be67055b 1599 AliExternalTrackParam pTPC(*(friendTrack->GetTPCOut()));
3e55050f 1600 AliExternalTrackParam pTOF(pTPC);
1601 Double_t mass = TDatabasePDG::Instance()->GetParticle("mu+")->Mass();
817766d5 1602 Int_t npoints = points->GetNPoints();
1603 AliTrackPoint point;
3e55050f 1604 Int_t naccept=0;
817766d5 1605 //
1606 for (Int_t ipoint=0;ipoint<npoints;ipoint++){
817766d5 1607 points->GetPoint(point,ipoint);
1608 Float_t xyz[3];
1609 point.GetXYZ(xyz);
1610 Double_t r=TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
3e55050f 1611 if (r<350) continue;
817766d5 1612 if (r>400) continue;
3e55050f 1613 AliTracker::PropagateTrackToBxByBz(&pTPC,r,mass,2.,kTRUE);
1614 AliTracker::PropagateTrackToBxByBz(&pTPC,r,mass,0.1,kTRUE);
1615 AliTrackPoint lpoint = point.Rotate(pTPC.GetAlpha());
1616 pTPC.PropagateTo(lpoint.GetX(),fMagF);
1617 pTOF=pTPC;
817766d5 1618 ((Double_t*)pTOF.GetParameter())[0] =lpoint.GetY();
1619 ((Double_t*)pTOF.GetParameter())[1] =lpoint.GetZ();
3e55050f 1620 ((Double_t*)pTOF.GetCovariance())[0]+=3.*3./12.;
1621 ((Double_t*)pTOF.GetCovariance())[2]+=3.*3./12.;
0dac7b3a 1622 ((Double_t*)pTOF.GetCovariance())[5]+=0.1*0.1;
1623 ((Double_t*)pTOF.GetCovariance())[9]+=0.1*0.1;
3e55050f 1624 naccept++;
1625 }
1626 if (naccept==0) return; // no tof match clusters
1627 //
1628 // 0. Apply standard cuts
1629 //
1630 if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
1631 // exclude crossing track
be67055b 1632 if (friendTrack->GetTPCOut()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
3e55050f 1633 //
1634 if (TMath::Abs(pTOF.GetY()-pTPC.GetY()) >kMaxDy) return;
1635 if (TMath::Abs(pTOF.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
1636 if (TMath::Abs(pTOF.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
1637 //
1638 // 1. Update median and RMS info
1639 //
b322e06a 1640 TVectorD vecDelta(5),vecMedian(5), vecRMS(5);
3e55050f 1641 TVectorD vecDeltaN(5);
1642 Double_t sign=(pTOF.GetParameter()[1]>0)? 1.:-1.;
1643 vecDelta[4]=0;
1644 for (Int_t i=0;i<4;i++){
1645 vecDelta[i]=(pTOF.GetParameter()[i]-pTPC.GetParameter()[i])*sign;
1646 kgdP[i][kglast%kN]=vecDelta[i];
1647 }
1648 kglast=(kglast+1);
1649 Int_t entries=(kglast<kN)?kglast:kN;
1650 Bool_t isOK=kTRUE;
1651 for (Int_t i=0;i<4;i++){
1652 vecMedian[i] = TMath::Median(entries,kgdP[i]);
1653 vecRMS[i] = TMath::RMS(entries,kgdP[i]);
1654 vecDeltaN[i] = 0;
1655 if (vecRMS[i]>0.){
1656 vecDeltaN[i] = (vecDelta[i]-vecMedian[i])/(vecRMS[i]+1.);
1657 vecDeltaN[4]+= TMath::Abs(vecDeltaN[i]); //sum of abs residuals
1658 if (TMath::Abs(vecDeltaN[i])>kSigmaCut) isOK=kFALSE;
817766d5 1659 }
3e55050f 1660 }
1661 //
1662 // 2. Apply median+-rms cut
1663 //
1664 if (kglast<10) return; //median and RMS to be defined
1665 if (!isOK) return;
1666 //
1667 // 3. Update alignment
1668 //
1669 Int_t htime = fTime/3600; //time in hours
1670 if (fAlignTOFTPC->GetEntries()<htime){
1671 fAlignTOFTPC->Expand(htime*2+20);
1672 }
1673 AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTOFTPC->At(htime);
1674 if (!align){
1675 // make Alignment object if doesn't exist
1676 align=new AliRelAlignerKalman();
1677 align->SetRunNumber(fRun);
1678 (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
1679 (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
1680 align->SetOutRejSigma(kOutCut+kOutCut*kN);
1681 align->SetRejectOutliers(kFALSE);
1682 align->SetTPCvd(AliTPCcalibDB::Instance()->GetParameters()->GetDriftV()/1000000.);
1683 align->SetMagField(fMagF);
1684 fAlignTOFTPC->AddAt(align,htime);
1685 }
1686 align->AddTrackParams(&pTOF,&pTPC);
1687 align->SetTimeStamp(fTime);
1688 align->SetRunNumber(fRun );
1689 //
1690 Int_t nupdates=align->GetNUpdates();
1691 align->SetOutRejSigma(kOutCut+kOutCut*kN/Double_t(nupdates));
1692 align->SetRejectOutliers(kFALSE);
1693 TTreeSRedirector *cstream = GetDebugStreamer();
1694 if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
1695 TTimeStamp tstamp(fTime);
1696 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
1697 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
1698 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
1699 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
1700 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
1701 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
1702 TVectorD vecGoofie(20);
1703 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
1704 if (goofieArray){
1705 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
1706 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
1707 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
1708 }
817766d5 1709 }
3e55050f 1710 TVectorD gpTPC(3), gdTPC(3);
1711 TVectorD gpTOF(3), gdTOF(3);
1712 pTPC.GetXYZ(gpTPC.GetMatrixArray());
1713 pTPC.GetDirection(gdTPC.GetMatrixArray());
1714 pTOF.GetXYZ(gpTOF.GetMatrixArray());
1715 pTOF.GetDirection(gdTOF.GetMatrixArray());
32100b1d 1716 Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
3e55050f 1717 (*cstream)<<"toftpc"<<
1718 "run="<<fRun<< // run number
1719 "event="<<fEvent<< // event number
1720 "time="<<fTime<< // time stamp of event
1721 "trigger="<<fTrigger<< // trigger
1722 "mag="<<fMagF<< // magnetic field
1723 // Environment values
1724 "press0="<<valuePressure0<<
1725 "press1="<<valuePressure1<<
1726 "pt0="<<ptrelative0<<
1727 "pt1="<<ptrelative1<<
1728 "temp0="<<temp0<<
1729 "temp1="<<temp1<<
1730 "vecGoofie.="<<&vecGoofie<<
32100b1d 1731 "vdcorr="<<vdcorr<< // drift correction applied
3e55050f 1732 //
1733 "nmed="<<kglast<< // number of entries to define median and RMS
1734 "vMed.="<<&vecMedian<< // median of deltas
1735 "vRMS.="<<&vecRMS<< // rms of deltas
1736 "vDelta.="<<&vecDelta<< // delta in respect to median
1737 "vDeltaN.="<<&vecDeltaN<< // normalized delta in respect to median
1738 "t.="<<track<< // ful track - find proper cuts
1739 "a.="<<align<< // current alignment
1740 "pTOF.="<<&pTOF<< // track param TOF
1741 "pTPC.="<<&pTPC<< // track param TPC
1742 "gpTPC.="<<&gpTPC<< // global position TPC
1743 "gdTPC.="<<&gdTPC<< // global direction TPC
1744 "gpTOF.="<<&gpTOF<< // global position TOF
1745 "gdTOF.="<<&gdTOF<< // global position TOF
1746 "\n";
817766d5 1747 }
817766d5 1748}
3e55050f 1749
1750
391ffdb2 1751void AliTPCcalibTime::BookDistortionMaps(){
1752 //
1753 // Book ndimensional histograms of distortions/residuals
1754 // Only primary tracks are selected for analysis
1755 //
1756
1757 Double_t xminTrack[4], xmaxTrack[4];
1758 Int_t binsTrack[4];
1759 TString axisName[4];
1760 TString axisTitle[4];
1761 //
882f5c06 1762 binsTrack[0] =50;
391ffdb2 1763 axisName[0] ="#Delta";
1764 axisTitle[0] ="#Delta";
1765 //
1766 binsTrack[1] =20;
1767 xminTrack[1] =-1.5; xmaxTrack[1]=1.5;
1768 axisName[1] ="tanTheta";
1769 axisTitle[1] ="tan(#Theta)";
1770 //
1771 binsTrack[2] =90;
1772 xminTrack[2] =-TMath::Pi(); xmaxTrack[2]=TMath::Pi();
1773 axisName[2] ="phi";
1774 axisTitle[2] ="#phi";
1775 //
1776 binsTrack[3] =20;
1777 xminTrack[3] =-1.; xmaxTrack[3]=1.; // 0.33 GeV cut
1778 axisName[3] ="snp";
1779 //
1780 // delta y
be67055b 1781 xminTrack[0] =-1.5; xmaxTrack[0]=1.5; //
391ffdb2 1782 fResHistoTPCITS[0] = new THnSparseS("TPCITS#Delta_{Y} (cm)","#Delta_{Y} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
1783 fResHistoTPCvertex[0] = new THnSparseS("TPCVertex#Delta_{Y} (cm)","#Delta_{Y} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
be67055b 1784 xminTrack[0] =-1.5; xmaxTrack[0]=1.5; //
391ffdb2 1785 fResHistoTPCTRD[0] = new THnSparseS("TPCTRD#Delta_{Y} (cm)","#Delta_{Y} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
1786 //
1787 // delta z
9112627b 1788 xminTrack[0] =-3.; xmaxTrack[0]=3.; //
391ffdb2 1789 fResHistoTPCITS[1] = new THnSparseS("TPCITS#Delta_{Z} (cm)","#Delta_{Z} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
1790 fResHistoTPCvertex[1] = new THnSparseS("TPCVertex#Delta_{Z} (cm)","#Delta_{Z} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
1791 fResHistoTPCTRD[1] = new THnSparseS("TPCTRD#Delta_{Z} (cm)","#Delta_{Z} (cm)", 4, binsTrack,xminTrack, xmaxTrack);
1792 //
1793 // delta snp-P2
882f5c06 1794 xminTrack[0] =-0.015; xmaxTrack[0]=0.015; //
391ffdb2 1795 fResHistoTPCITS[2] = new THnSparseS("TPCITS#Delta_{#phi}","#Delta_{#phi}", 4, binsTrack,xminTrack, xmaxTrack);
9112627b 1796 fResHistoTPCvertex[2] = new THnSparseS("TPCITSv#Delta_{#phi}","#Delta_{#phi}", 4, binsTrack,xminTrack, xmaxTrack);
391ffdb2 1797 fResHistoTPCTRD[2] = new THnSparseS("TPCTRD#Delta_{#phi}","#Delta_{#phi}", 4, binsTrack,xminTrack, xmaxTrack);
1798 //
1799 // delta theta-P3
9112627b 1800 xminTrack[0] =-0.025; xmaxTrack[0]=0.025; //
391ffdb2 1801 fResHistoTPCITS[3] = new THnSparseS("TPCITS#Delta_{#theta}","#Delta_{#theta}", 4, binsTrack,xminTrack, xmaxTrack);
9112627b 1802 fResHistoTPCvertex[3] = new THnSparseS("TPCITSv#Delta_{#theta}","#Delta_{#theta}", 4, binsTrack,xminTrack, xmaxTrack);
391ffdb2 1803 fResHistoTPCTRD[3] = new THnSparseS("TPCTRD#Delta_{#theta}","#Delta_{#theta}", 4, binsTrack,xminTrack, xmaxTrack);
1804 //
1805 // delta theta-P4
9112627b 1806 xminTrack[0] =-0.2; xmaxTrack[0]=0.2; //
391ffdb2 1807 fResHistoTPCITS[4] = new THnSparseS("TPCITS#Delta_{1/pt}","#Delta_{1/pt}", 4, binsTrack,xminTrack, xmaxTrack);
9112627b 1808 fResHistoTPCvertex[4] = new THnSparseS("TPCITSv#Delta_{1/pt}","#Delta_{1/pt}", 4, binsTrack,xminTrack, xmaxTrack);
391ffdb2 1809 fResHistoTPCTRD[4] = new THnSparseS("TPCTRD#Delta_{1/pt}","#Delta_{1/pt}", 4, binsTrack,xminTrack, xmaxTrack);
1810 //
1811 for (Int_t ivar=0;ivar<4;ivar++){
1812 for (Int_t ivar2=0;ivar2<4;ivar2++){
1813 fResHistoTPCITS[ivar]->GetAxis(ivar2)->SetName(axisName[ivar2].Data());
1814 fResHistoTPCITS[ivar]->GetAxis(ivar2)->SetTitle(axisTitle[ivar2].Data());
1815 fResHistoTPCTRD[ivar]->GetAxis(ivar2)->SetName(axisName[ivar2].Data());
1816 fResHistoTPCTRD[ivar]->GetAxis(ivar2)->SetTitle(axisTitle[ivar2].Data());
1817 fResHistoTPCvertex[ivar]->GetAxis(ivar2)->SetName(axisName[ivar2].Data());
1818 fResHistoTPCvertex[ivar]->GetAxis(ivar2)->SetTitle(axisTitle[ivar2].Data());
1819 }
1820 }
1821}
1822
1823
1824void AliTPCcalibTime::FillResHistoTPCITS(const AliExternalTrackParam * pTPCIn, const AliExternalTrackParam * pITSOut ){
1825 //
1826 // fill residual histograms pTPCIn-pITSOut
1827 // Histogram is filled only for primary tracks
1828 //
1829 Double_t histoX[4];
1830 Double_t xyz[3];
1831 pTPCIn->GetXYZ(xyz);
1832 Double_t phi= TMath::ATan2(xyz[1],xyz[0]);
1833 histoX[1]= pTPCIn->GetTgl();
1834 histoX[2]= phi;
1835 histoX[3]= pTPCIn->GetSnp();
9112627b 1836 AliExternalTrackParam lits(*pITSOut);
1837 lits.Rotate(pTPCIn->GetAlpha());
1838 lits.PropagateTo(pTPCIn->GetX(),fMagF);
391ffdb2 1839 //
1840 for (Int_t ihisto=0; ihisto<5; ihisto++){
9112627b 1841 histoX[0]=pTPCIn->GetParameter()[ihisto]-lits.GetParameter()[ihisto];
391ffdb2 1842 fResHistoTPCITS[ihisto]->Fill(histoX);
1843 }
1844}
1845
1846
1847void AliTPCcalibTime::FillResHistoTPC(const AliESDtrack * pTrack){
1848 //
1849 // fill residual histograms pTPC - vertex
1850 // Histogram is filled only for primary tracks
1851 //
1852 Double_t histoX[4];
1853 const AliExternalTrackParam * pTPCIn = pTrack->GetInnerParam();
be67055b 1854 AliExternalTrackParam pTPCvertex(*(pTrack->GetInnerParam()));
9112627b 1855 //
1856 AliExternalTrackParam lits(*pTrack);
be67055b 1857 if (TMath::Abs(pTrack->GetY())>3) return; // beam pipe
1858 pTPCvertex.Rotate(lits.GetAlpha());
1859 //pTPCvertex.PropagateTo(pTPCvertex->GetX(),fMagF);
1860 AliTracker::PropagateTrackToBxByBz(&pTPCvertex,lits.GetX(),0.1,2,kFALSE);
1861 AliTracker::PropagateTrackToBxByBz(&pTPCvertex,lits.GetX(),0.1,0.1,kFALSE);
391ffdb2 1862 Double_t xyz[3];
1863 pTPCIn->GetXYZ(xyz);
1864 Double_t phi= TMath::ATan2(xyz[1],xyz[0]);
1865 histoX[1]= pTPCIn->GetTgl();
1866 histoX[2]= phi;
1867 histoX[3]= pTPCIn->GetSnp();
1868 //
1869 Float_t dca[2], cov[3];
1870 pTrack->GetImpactParametersTPC(dca,cov);
882f5c06 1871 for (Int_t ihisto=0; ihisto<5; ihisto++){
be67055b 1872 histoX[0]=pTPCvertex.GetParameter()[ihisto]-lits.GetParameter()[ihisto];
1873 // if (ihisto<2) histoX[0]=dca[ihisto];
9112627b 1874 fResHistoTPCvertex[ihisto]->Fill(histoX);
391ffdb2 1875 }
1876}
1877
1878
1879void AliTPCcalibTime::FillResHistoTPCTRD(const AliExternalTrackParam * pTPCOut, const AliExternalTrackParam * pTRDIn ){
1880 //
1881 // fill resuidual histogram TPCout-TRDin
1882 //
1883 Double_t histoX[4];
1884 Double_t xyz[3];
1885 pTPCOut->GetXYZ(xyz);
1886 Double_t phi= TMath::ATan2(xyz[1],xyz[0]);
1887 histoX[1]= pTPCOut->GetTgl();
1888 histoX[2]= phi;
1889 histoX[3]= pTPCOut->GetSnp();
1890 //
882f5c06 1891 AliExternalTrackParam ltrd(*pTRDIn);
1892 ltrd.Rotate(pTPCOut->GetAlpha());
be67055b 1893 // ltrd.PropagateTo(pTPCOut->GetX(),fMagF);
1894 AliTracker::PropagateTrackToBxByBz(&ltrd,pTPCOut->GetX(),0.1,0.1,kFALSE);
882f5c06 1895
391ffdb2 1896 for (Int_t ihisto=0; ihisto<5; ihisto++){
882f5c06 1897 histoX[0]=pTPCOut->GetParameter()[ihisto]-ltrd.GetParameter()[ihisto];
391ffdb2 1898 fResHistoTPCTRD[ihisto]->Fill(histoX);
1899 }
1900
1901}