]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TOF/AliTOFCalibTask.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
New default values for baselines (F.Prino)
[u/mrichter/AliRoot.git] / TOF / AliTOFCalibTask.cxx
CommitLineData
1ed77580 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/*
45af7330 16$Log$
ea23eb9f 17Revision 1.6 2007/10/05 10:44:11 zampolli
18Oversight for debugging purpose fixed
19
b8047c99 20Revision 1.5 2007/10/05 10:38:10 zampolli
21Oversight fixed
22
1e76f1ec 23Revision 1.4 2007/10/04 15:36:44 zampolli
24Updates to new TOF offline calibration schema
25
e26353e8 26Revision 1.3 2007/07/31 07:26:16 zampolli
27Bug fixed in the $ field
28
1ed77580 29*/
30
31// task to perform TOF calibration
32// optimized for pp runs
33// expect a maximum of 100 entries per channel
34// different ways to calibrate are implemented
35// C. Zampolli
36
37#include <TChain.h>
38#include <TObject.h>
39#include <TCanvas.h>
40#include <TROOT.h>
41#include <TSystem.h>
42#include <TProfile.h>
43#include <TF1.h>
44#include <TTree.h>
45#include <TFile.h>
46#include <TH1F.h>
47#include <TH1I.h>
48#include <TH1D.h>
49#include <TH2F.h>
e26353e8 50#include <TDirectory.h>
51#include <TGrid.h>
1ed77580 52
53#include "AliTOFCalibTask.h"
70ddd0ee 54#include "AliESDEvent.h"
1ed77580 55#include "AliESDtrack.h"
56#include "AliLog.h"
e26353e8 57#include "AliCDBManager.h"
58#include "AliCDBMetaData.h"
59#include "AliCDBId.h"
60#include "AliCDBStorage.h"
61
1ed77580 62//______________________________________________________________________________
63AliTOFCalibTask::AliTOFCalibTask(const char *name) :
64 AliAnalysisTask(name,""),
65 fdir(0),
66 fChain(0),
67 fESD(0),
1ed77580 68 fToT(0),
69 fTime(0),
70 fExpTimePi(0),
71 fExpTimeKa(0),
72 fExpTimePr(0),
1ed77580 73 fMinTime(0),
74 fbigarray(0x0),
75 findexarray(0x0),
76 fnESD(0),
77 fnESDselected(0),
78 fnESDkTOFout(0),
79 fnESDkTIME(0),
80 fnESDassTOFcl(0),
81 fnESDTIMEcut(0),
82 fnESDTRDcut(0),
83 fhToT(0),
84 fhTime(0),
85 fhExpTimePi(0),
86 fhExpTimeKa(0),
87 fhExpTimePr(0),
88 fhPID(0),
89 fhch(0),
e26353e8 90 fOutputContainer(0),
91 frun(0)
1ed77580 92{
93 // Constructor.
94 // Input slot #0 works with an Ntuple
95
96 DefineInput(0, TChain::Class());
97 // Output slot #0 writes into a TH1 container
98 DefineOutput(0, TObjArray::Class()) ;
99 fdir = gDirectory->GetPath();
70ddd0ee 100 fbigarray = new Float_t*[TOFCHANNELS];
101 findexarray = new Int_t[TOFCHANNELS];
102
e26353e8 103 for (Int_t i=0;i<11;i++){
104 fassparticle[i]=-1;
105 }
106
70ddd0ee 107 for (Int_t i =0;i<TOFCHANNELS;i++){
108 fbigarray[i] = new Float_t[CHENTRIES];
109 findexarray[i]=0;
110 for (Int_t j =0;j<CHENTRIES;j++){
111 fbigarray[i][j]=-1;
112 }
113 }
1ed77580 114}
115
116//______________________________________________________________________________
117AliTOFCalibTask::AliTOFCalibTask(const AliTOFCalibTask &calibtask) :
118 AliAnalysisTask("AliTOFCalibTask",""),
119 fdir(0),
120 fChain(0),
121 fESD(0),
1ed77580 122 fToT(0),
123 fTime(0),
124 fExpTimePi(0),
125 fExpTimeKa(0),
126 fExpTimePr(0),
1ed77580 127 fMinTime(0),
128 fbigarray(0x0),
129 findexarray(0x0),
130 fnESD(0),
131 fnESDselected(0),
132 fnESDkTOFout(0),
133 fnESDkTIME(0),
134 fnESDassTOFcl(0),
135 fnESDTIMEcut(0),
136 fnESDTRDcut(0),
137 fhToT(0),
138 fhTime(0),
139 fhExpTimePi(0),
140 fhExpTimeKa(0),
141 fhExpTimePr(0),
142 fhPID(0),
143 fhch(0),
e26353e8 144 fOutputContainer(0),
145 frun(0)
1ed77580 146{
147 // Copy Constructor.
148 fdir=calibtask.fdir;
149 fChain=calibtask.fChain;
150 fESD=calibtask.fESD;
1ed77580 151 fToT=calibtask.fToT;
152 fTime=calibtask.fTime;
153 fExpTimePi=calibtask.fExpTimePi;
154 fExpTimeKa=calibtask.fExpTimeKa;
155 fExpTimePr=calibtask.fExpTimePr;
1ed77580 156 fMinTime=calibtask.fMinTime;
1ed77580 157 fnESD=calibtask.fnESD;
158 fnESDselected=calibtask.fnESDselected;
159 fnESDkTOFout=calibtask.fnESDkTOFout;
160 fnESDkTIME=calibtask.fnESDkTIME;
161 fnESDassTOFcl=calibtask.fnESDassTOFcl;
162 fnESDTIMEcut=calibtask.fnESDTIMEcut;
163 fnESDTRDcut=calibtask.fnESDTRDcut;
164 fhToT=calibtask.fhToT;
165 fhTime=calibtask.fhTime;
166 fhExpTimePi=calibtask.fhExpTimePi;
167 fhExpTimeKa=calibtask.fhExpTimeKa;
168 fhExpTimePr=calibtask.fhExpTimePr;
169 fhPID=calibtask.fhPID;
170 fhch=calibtask.fhch;
e26353e8 171 frun=calibtask.frun;
1ed77580 172 fOutputContainer=calibtask.fOutputContainer;
70ddd0ee 173
174 fbigarray = new Float_t*[TOFCHANNELS];
175 findexarray = new Int_t[TOFCHANNELS];
176
e26353e8 177 for (Int_t i=0;i<11;i++){
178 fassparticle[i]=calibtask.fassparticle[i];
179 }
180
70ddd0ee 181 for (Int_t i =0;i<TOFCHANNELS;i++){
182 fbigarray[i] = new Float_t[CHENTRIES];
183 findexarray[i]=calibtask.findexarray[i];
184 for (Int_t j =0;j<CHENTRIES;j++){
185 fbigarray[i][j]=calibtask.fbigarray[i][j];
186 }
187 }
188
1ed77580 189}
190//______________________________________________________________________________
191AliTOFCalibTask:: ~AliTOFCalibTask()
192{
193 // destructor
194
195 AliInfo("TOF Calib Task: Deleting");
e88f3330 196 for (Int_t i =0;i<TOFCHANNELS;i++){
197 delete [] fbigarray[i];
198 }
1ed77580 199 delete[] fbigarray;
ea23eb9f 200 delete[] findexarray;
1ed77580 201 delete fOutputContainer;
202 delete fhToT;
203 delete fhTime;
204 delete fhExpTimePi;
205 delete fhExpTimeKa;
206 delete fhExpTimePr;
207 delete fhPID;
208 delete fhch;
209}
210//______________________________________________________________________________
211AliTOFCalibTask& AliTOFCalibTask::operator=(const AliTOFCalibTask &calibtask)
212{
213 //assignment operator
214 this->fdir=calibtask.fdir;
215 this->fChain=calibtask.fChain;
216 this->fESD=calibtask.fESD;
1ed77580 217 this->fToT=calibtask.fToT;
218 this->fTime=calibtask.fTime;
219 this->fExpTimePi=calibtask.fExpTimePi;
220 this->fExpTimeKa=calibtask.fExpTimeKa;
221 this->fExpTimePr=calibtask.fExpTimePr;
1ed77580 222 this->fMinTime=calibtask.fMinTime;
1ed77580 223 this->fnESD=calibtask.fnESD;
224 this->fnESDselected=calibtask.fnESDselected;
225 this->fnESDkTOFout=calibtask.fnESDkTOFout;
226 this->fnESDkTIME=calibtask.fnESDkTIME;
227 this->fnESDassTOFcl=calibtask.fnESDassTOFcl;
228 this->fnESDTIMEcut=calibtask.fnESDTIMEcut;
229 this->fnESDTRDcut=calibtask.fnESDTRDcut;
230 this->fhToT=calibtask.fhToT;
231 this->fhTime=calibtask.fhTime;
232 this->fhExpTimePi=calibtask.fhExpTimePi;
233 this->fhExpTimeKa=calibtask.fhExpTimeKa;
234 this->fhExpTimePr=calibtask.fhExpTimePr;
235 this->fOutputContainer=calibtask.fOutputContainer;
236 this->fhPID=calibtask.fhPID;
237 this->fhch=calibtask.fhch;
70ddd0ee 238
239 this->fbigarray = new Float_t*[TOFCHANNELS];
240 this->findexarray = new Int_t[TOFCHANNELS];
e26353e8 241
242 for (Int_t i=0;i<11;i++){
243 this->fassparticle[i]=calibtask.fassparticle[i];
244 }
70ddd0ee 245 for (Int_t i =0;i<TOFCHANNELS;i++){
246 this->fbigarray[i] = new Float_t[CHENTRIES];
247 this->findexarray[i]=calibtask.findexarray[i];
248 for (Int_t j =0;j<CHENTRIES;j++){
249 this->fbigarray[i][j]=calibtask.fbigarray[i][j];
250 }
251 }
1ed77580 252 return *this;
253}
254//--------------------------------------------------------------------------
255void AliTOFCalibTask::BookHistos(){
256
257 //booking histos
258
259 AliInfo(Form("*** Booking Histograms %s", GetName())) ;
260
261 fhToT=
262 new TH1F("hToT", " ToT distribution (ns) ", 400, 0, 40);
263 fhTime=
264 new TH1F("hTime", " Time distribution (ns) ", 400, 0, 40);
265 fhExpTimePi=
266 new TH1F("hExpTimePi", " Exp Time distribution, Pi (ns) ", 400, 0, 40);
267 fhExpTimeKa=
268 new TH1F("hExpTimeKa", " Exp Time distribution, Ka (ns) ", 400, 0, 40);
269 fhExpTimePr=
270 new TH1F("hExpTimePr", " Exp Time distribution, Pr (ns) ", 400, 0, 40);
271 fhPID=
272 new TH1I("hPID", " Combinatorial PID identity ", 3, 0, 3);
273 fhch=
274 new TH1D("hch", " TOF channel ", TOFCHANNELS, 0, TOFCHANNELS);
275
276 //create the putput container
277 fOutputContainer = new TObjArray(7) ;
278 fOutputContainer->SetName(GetName()) ;
279
280 fOutputContainer->AddAt(fhToT, 0) ;
281 fOutputContainer->AddAt(fhTime, 1) ;
282 fOutputContainer->AddAt(fhExpTimePi, 2) ;
283 fOutputContainer->AddAt(fhExpTimeKa, 3) ;
284 fOutputContainer->AddAt(fhExpTimePr, 4) ;
285 fOutputContainer->AddAt(fhPID, 5) ;
286 fOutputContainer->AddAt(fhch, 6) ;
287
288}
289//----------------------------------------------------------------------------
290void AliTOFCalibTask::DrawHistos(){
291
292 // drawing output histos
293
294 AliInfo(Form("*** Drawing Histograms %s", GetName())) ;
295
296 TCanvas * canvasToTTime = new TCanvas("canvasToTTime", " ToT and Time ",400, 30, 550, 630) ;
297 canvasToTTime->Divide(1,2);
298 canvasToTTime->cd(1);
299 fhToT->SetLineColor(4);
300 fhToT->GetXaxis()->SetTitle("ToT (ns)");
301 fhToT->Draw("hist");
302 canvasToTTime->cd(2);
303 fhTime->SetLineColor(4);
304 fhTime->GetXaxis()->SetTitle("Time (ns)");
305 fhTime->Draw("hist");
306 canvasToTTime->Update();
307 canvasToTTime->Print("ToTTime.gif");
308
309 TCanvas * canvasExpTime = new TCanvas("canvasExpTime", " Expected Times ",400, 30, 550, 630) ;
310 canvasExpTime->Divide(1,3);
311 canvasExpTime->cd(1);
312 fhExpTimePi->SetLineColor(4);
313 fhExpTimePi->GetXaxis()->SetTitle("Exp Time (ns), #pi");
314 fhExpTimePi->Draw("hist");
315 canvasExpTime->cd(2);
316 fhExpTimeKa->SetLineColor(4);
317 fhExpTimeKa->GetXaxis()->SetTitle("Exp Time (ns), K");
318 fhExpTimeKa->Draw("hist");
319 canvasExpTime->cd(3);
320 fhExpTimePr->SetLineColor(4);
321 fhExpTimePr->GetXaxis()->SetTitle("Exp Time (ns), p");
322 fhExpTimePr->Draw("hist");
323
324 canvasExpTime->Print("ExpTime.gif");
325
326 TCanvas * canvasPID = new TCanvas("canvasPID", " Combinatorial PID ",400, 30, 550, 400);
327 fhPID->GetXaxis()->SetTitle("Comb PID");
328 fhPID->GetXaxis()->SetBinLabel(1,"#pi");
329 fhPID->GetXaxis()->SetBinLabel(2,"K");
330 fhPID->GetXaxis()->SetBinLabel(3,"p");
331 fhPID->Draw("hist");
332
333 canvasPID->Print("PID.gif");
334
335 TCanvas * canvasrndch = new TCanvas("canvasrndch", " TOF channel ",400, 30, 550, 400);
336 fhch->GetXaxis()->SetTitle("TOF ch");
337 fhch->Draw("hist");
338 Float_t meanTOFch = 0;
339 for (Int_t ibin=0;ibin<TOFCHANNELS;ibin++){
340 meanTOFch+=(Float_t)fhch->GetBinContent(ibin+1);
341 }
342
343 meanTOFch/=TOFCHANNELS;
344 AliDebug(1,Form(" Mean number of tracks/channel = %f ",meanTOFch));
345
346 canvasrndch->Print("rndch.gif");
347
348 char line[1024] ;
349 sprintf(line, ".!tar -zcvf %s.tar.gz *.gif", GetName()) ;
350 gROOT->ProcessLine(line);
351 sprintf(line, ".!rm -fR *.gif");
352 gROOT->ProcessLine(line);
353 AliInfo(Form("*** TOF Calib Task: plots saved in %s.tar.gz...\n", GetName())) ;
354}
355
356//______________________________________________________________________________
357void AliTOFCalibTask::ConnectInputData(const Option_t*)
358{
359 // Initialization of branch container and histograms
360
361 // AliLog::SetClassDebugLevel("AliTOFCalibTask",1);
362 AliInfo(Form("*** Initialization of %s", GetName())) ;
363
364 // Get input data
365 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
366 if (!fChain) {
367 AliError(Form("Input 0 for %s not found\n", GetName()));
368 return ;
369 }
370
371 // One should first check if the branch address was taken by some other task
372 char ** address = (char **)GetBranchAddress(0, "ESD");
373 if (address) {
70ddd0ee 374 fESD = (AliESDEvent*)(*address);
1ed77580 375 } else {
70ddd0ee 376 fESD = new AliESDEvent();
1ed77580 377 }
1e76f1ec 378 fESD->ReadFromTree(fChain) ;
1ed77580 379
380 BookHistos();
381
382}
70ddd0ee 383//-----------------------------------------------------------------------
384Bool_t AliTOFCalibTask::Notify()
385{
386 // Initialisation of branch container and histograms
387
388 AliInfo(Form("*** We are in %s::Notify()", GetName())) ;
389
390 // Get input data
391 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
392 if (!fChain) {
393 AliError(Form("Input 0 for %s not found\n", GetName()));
394 return kFALSE;
395 }
396
1e76f1ec 397 char ** address = (char **)GetBranchAddress(0, "ESD");
398 if (address) {
399 fESD = (AliESDEvent*)(*address);
400 } else {
401 fESD = new AliESDEvent();
402 }
e26353e8 403 fESD->ReadFromTree(fChain) ;
70ddd0ee 404
405 return kTRUE;
406}
407
1ed77580 408//________________________________________________________________________
409void AliTOFCalibTask::CreateOutputObjects()
410{
411// Create histograms
412}
413
414//______________________________________________________________________________
415void AliTOFCalibTask::Exec(Option_t * opt)
416{
417
418 // main
419
420 AliInfo(Form("*** Executing %s", GetName())) ;
421
422//******* The loop over events -----------------------------------------------
423
424 // Processing of one event
425 Long64_t entry = fChain->GetReadEntry() ;
426 if (!fESD) {
427 AliError("fESD is not connected to the input!") ;
428 return ;
429 }
430
431 if ( !((entry-1)%100) )
432 AliDebug(1,Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
433
434 // ************************ TOF *************************************
435
436 fMinTime=22E3; //ns; not used
437 Int_t ntrk = fESD->GetNumberOfTracks() ;
438 fnESD+=ntrk;
439 Int_t nselected = 0;
440 Int_t itr = -1;
441 while ( ntrk-- ) {
442 itr++;
443 AliESDtrack * t = fESD->GetTrack(ntrk) ;
444 //selecting only good quality tracks
445 if (!Select(t)) continue;
446 nselected++;
447 Int_t ich = Int_t(t->GetTOFCalChannel());
448 fhch->Fill(ich);
b8047c99 449 // ich=3; //only for debug purpose
1ed77580 450 AliDebug(2,Form(" ESD in channel %i, filling array %i",t->GetTOFCalChannel(),ich));
451 Float_t tot = t->GetTOFsignalToT();
452 Float_t time = t->GetTOFsignalRaw();
e26353e8 453 AliDebug(2,Form(" track # %i in channel %i, time = %f \n",ntrk,ich,time));
1ed77580 454 Double_t expTime[10];
455 t->GetIntegratedTimes(expTime);
456 Float_t expTimePi = expTime[2]*1.E-3;
457 Float_t expTimeKa = expTime[3]*1.E-3;
458 Float_t expTimePr = expTime[4]*1.E-3;
459 if (findexarray[ich]==(Int_t)(CHENTRIES/NIDX)) {
460 AliInfo(Form("too many tracks in channel %i, not storing any more...",ich));
461 continue;
462 }
463 findexarray[ich]++;
464 AliDebug(2,Form("tracks in channel %i = %i, storing... ",ich, findexarray[ich] ));
465 Int_t ientry=(findexarray[ich]-1)*NIDX;
466 fbigarray[ich][ientry+DELTAIDXTOT]=tot; //in ns
467 fbigarray[ich][ientry+DELTAIDXTIME]=time*1E-3; // in ns
468 fbigarray[ich][ientry+DELTAIDXEXTIMEPI]=expTimePi;
469 fbigarray[ich][ientry+DELTAIDXEXTIMEKA]=expTimeKa;
470 fbigarray[ich][ientry+DELTAIDXEXTIMEPR]=expTimePr;
471 fhToT->Fill(fbigarray[ich][ientry+DELTAIDXTOT]);
472 fhTime->Fill(fbigarray[ich][ientry+DELTAIDXTIME]);
473 fhExpTimePi->Fill(fbigarray[ich][ientry+DELTAIDXEXTIMEPI]);
474 fhExpTimeKa->Fill(fbigarray[ich][ientry+DELTAIDXEXTIMEKA]);
475 fhExpTimePr->Fill(fbigarray[ich][ientry+DELTAIDXEXTIMEPR]);
476 AliDebug(2, Form("track = %i, tot = %f, time = %f, and Exp time in TOF: pi = %f, K = %f, p = %f",itr, fbigarray[ich][ientry+DELTAIDXTOT], fbigarray[ich][ientry+DELTAIDXTIME], expTimePi,expTimeKa,expTimePr));
477 }
478 fnESDselected+=nselected;
479
480 PostData(0, fOutputContainer);
481}
482
483//_____________________________________________________________________________
484void AliTOFCalibTask::Terminate(Option_t *)
485{
486 // Processing when the event loop is ended
487
488 // some plots
489
490 TH1::AddDirectory(0);
e26353e8 491 TDirectory *dir = gDirectory;
1ed77580 492 AliInfo("TOF Calib Task: End of events loop");
e26353e8 493 for (Int_t ich = 0;ich<TOFCHANNELS;ich++){
494 if (findexarray[ich]>0){
495 Int_t ncuttime = SelectOnTime(&fbigarray[ich][0],findexarray[ich],ich);
496 fnESDselected-=ncuttime;
497 }
498 }
499
1ed77580 500 AliInfo(Form(" Number of analyzed ESD tracks: %i\n",fnESD));
501 AliInfo(Form(" Number of selected ESD tracks: %i\n",fnESDselected));
502 AliInfo(Form(" Number of ESD tracks with kTOFout: %i\n",fnESDkTOFout));
503 AliInfo(Form(" Number of ESD tracks with kTIME: %i\n",fnESDkTIME));
504 AliInfo(Form(" Number of ESD tracks with TRDcut: %i\n",fnESDTRDcut));
505 AliInfo(Form(" Number of ESD tracks with TIMEcut: %i\n",fnESDTIMEcut));
506 AliInfo(Form(" Number of ESD tracks with assTOFcl: %i\n",fnESDassTOFcl));
507
508 for (Int_t i = 0;i<TOFCHANNELS;i++){
509 Int_t size=findexarray[i]*NIDX;
510 AliDebug(2, Form(" entries %i in channel %i ",findexarray[i],i));
e26353e8 511 if (findexarray[i]<=2) {
1ed77580 512 AliDebug(1, Form(" not enough statistics for combined PID for channel %i, putting all the tracks as if they were pions",i));
513 continue;
514 }
515 if (!CombPID(&fbigarray[i][0], size)) AliError("ERROR!!!!ERROR!!!");
516 }
70ddd0ee 517
e26353e8 518
1ed77580 519 DrawHistos();
520
521 // saving data in a tree
522
e26353e8 523 AliInfo("Building tree for Calibration");
524 TTree * tree = new TTree("T", "Tree for TOF calibration");
525 Float_t p[CHENTRIESSMALL];
1ed77580 526 Int_t nentries;
e26353e8 527 tree->Branch("nentries",&nentries,"nentries/I");
528 tree->Branch("TOFentries",p,"TOFentries[nentries]/F");
529 for (Int_t i=0;i<TOFCHANNELS;i++){
530 nentries=findexarray[i]*(NIDXSMALL); // when filling small array,
531 // only first 3 floats taken
532 // into account
533 for (Int_t j=0; j<findexarray[i];j++){
534 for (Int_t k=0; k<NIDXSMALL;k++){
535 Int_t index1= j*NIDXSMALL+k; // index in small array
536 Int_t index2=j*NIDX+k; // index in big array
537 p[index1]=fbigarray[i][index2];
538 }
1ed77580 539 }
e26353e8 540 tree->Fill();
1ed77580 541 }
1ed77580 542
e26353e8 543 AliInfo("Putting tree for calibration in Reference data");
544
545 // grid file option
546 Char_t filename[100];
547 sprintf(filename,"alien:///alice/cern.ch/user/c/czampolli/TOFCalibReference_%i.root",frun);
548 TGrid::Connect("alien://");
549 TFile *filegrid = TFile::Open(filename,"CREATE");
550 tree->Write();
551 dir->cd();
1ed77580 552
1ed77580 553}
1ed77580 554//_____________________________________________________________________________
555
556Bool_t AliTOFCalibTask::Select(AliESDtrack *t){
557
558 // selecting good quality tracks
559 //track selection: reconstrution to TOF:
560 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) {
561 return 0;
562 }
563 fnESDkTOFout++;
564 //IsStartedTimeIntegral
565 if ((t->GetStatus()&AliESDtrack::kTIME)==0) {
566 return 0;
567 }
568 fnESDkTIME++;
569 if (t->GetStatus() & AliESDtrack::kTRDbackup) {
570 Float_t xout= t->GetOuterParam()->GetX();
571 if (xout<364.25 && xout > 300.) return 0;
572 }
573 fnESDTRDcut++;
574 Double_t time=t->GetTOFsignal();
575 time*=1.E-3; // tof given in nanoseconds
576 if(time >= fMinTime){
577 return 0;
578 }
579 fnESDTIMEcut++;
580
581 Double_t mom=t->GetP();
582 if (!(mom<=UPPERMOMBOUND && mom>=LOWERMOMBOUND)){
583 // return 0; // skipping momentum cut
584 }
585
21a8ed8d 586 Int_t assignedTOFcluster=t->GetTOFcluster();//index of the assigned TOF cluster, >0 ?
587 if(assignedTOFcluster==-1){ // not matched
1ed77580 588 return 0;
589 }
590 fnESDassTOFcl++;
591 return 1;
592}
593//_____________________________________________________________________________
594
e26353e8 595Int_t AliTOFCalibTask::SelectOnTime(Float_t *charray, Int_t ntracks, Int_t ich){
596
597 // discarding tracks with time-mintime < MINTIME
598
599 Int_t ndeleted=0;
600 Float_t mintime = 1E6;
601 for (Int_t itr = 0;itr<ntracks;itr++){
602 Int_t ientry=itr*NIDX;
603 Float_t time = charray[ientry+DELTAIDXTIME];// in ns
604 if (time<mintime) mintime = time;
605 }
606 AliDebug(1,Form("Mintime for channel %i = %f",ich, mintime));
607 for (Int_t itr = 0;itr<ntracks;itr++){
608 Int_t ientry=itr*NIDX;
609 Float_t time = charray[ientry+DELTAIDXTIME];// in ns
610 if ((time-mintime)>MINTIME) {
611 ndeleted++;
612 AliDebug(1,Form("Deleting %i track from channel %i, time = %f",ndeleted, ich, time));
613 findexarray[ich]--;
614 for (Int_t j=itr+1;j<ntracks;j++){
615 Int_t ientry=j*NIDX;
616 Int_t idxtot = ientry+DELTAIDXTOT;
617 Int_t idxtime = ientry+DELTAIDXTIME;
618 Int_t idxextimePi = ientry+DELTAIDXEXTIMEPI;
619 Int_t idxextimeKa = ientry+DELTAIDXEXTIMEKA;
620 Int_t idxextimePr = ientry+DELTAIDXEXTIMEPR;
621 Int_t ientrydel=(j-1)*NIDX;
622 Int_t idxtotdel = ientrydel+DELTAIDXTOT;
623 Int_t idxtimedel = ientrydel+DELTAIDXTIME;
624 Int_t idxextimePidel = ientrydel+DELTAIDXEXTIMEPI;
625 Int_t idxextimeKadel = ientrydel+DELTAIDXEXTIMEKA;
626 Int_t idxextimePrdel = ientrydel+DELTAIDXEXTIMEPR;
627 charray[idxtotdel]=charray[idxtot];
628 charray[idxtimedel]=charray[idxtime];
629 charray[idxextimePidel]=charray[idxextimePi];
630 charray[idxextimeKadel]=charray[idxextimeKa];
631 charray[idxextimePrdel]=charray[idxextimePr];
632 }
633 }
634 }
635 return ndeleted;
636}
637//_____________________________________________________________________________
638
1ed77580 639Bool_t AliTOFCalibTask::CombPID(Float_t *smallarray, Int_t size){
640
e26353e8 641 // track Combinatorial PID for calibration, only when more than 2 particles
642 // fall in channel
1ed77580 643
1ed77580 644 Float_t t0offset=0;
e26353e8 645
646 if (size/NIDX<=2){
647 AliDebug(1,"Number of tracks in channel smaller than 2, identifying every particle as if it was a pion!");
1ed77580 648 return 0;
649 }
e26353e8 650
651 Int_t ntracksinchannel = size/NIDX;
652 Int_t ntrkinset=-1;
653 Int_t nset = ntracksinchannel/6;
654
655 if (nset ==0) {
656 nset=1;
657 if (ntracksinchannel < 6) {
658 AliDebug(2,"Number of tracks in set smaller than 6, Combinatorial PID still applied to this set.");
659 }
1ed77580 660 }
661
e26353e8 662 AliInfo(Form("number of sets = %i", nset));
663 // loop on sets
1ed77580 664 for (Int_t i=0; i< nset; i++) {
e26353e8 665 if (i<nset-1)ntrkinset=6;
666 else {
667 if (ntracksinchannel<6){
668 ntrkinset=ntracksinchannel;
669 }
670 else{
671 ntrkinset=6+Int_t((ntracksinchannel)%6);
672 }
673 }
674 AliInfo(Form("set = %i, number of tracks in set = %i", i,ntrkinset));
675
676 for (Int_t ii=0;ii<ntrkinset;ii++){
677 fassparticle[ii]=-1;
678 }
679 Float_t **exptof;
680 Float_t* timeofflight;
681 Float_t* texp;
682 Int_t* index;
683
684 exptof = new Float_t*[ntrkinset];
685 timeofflight = new Float_t[ntrkinset];
686 texp = new Float_t[ntrkinset];
687 index = new Int_t[ntrkinset];
688 for (Int_t ii=0;ii<ntrkinset;ii++){
689 exptof[ii] = new Float_t[3];
690 timeofflight[ii]=0;
691 texp[ii]=0;
692 index[ii]=-1;
693 }
694
695 for (Int_t j=0; j<ntrkinset; j++) {
696 Int_t idxtime = ((6*i+j)*NIDX)+DELTAIDXTIME;
697 Int_t idxextimePi = ((6*i+j)*NIDX)+DELTAIDXEXTIMEPI;
698 Int_t idxextimeKa = ((6*i+j)*NIDX)+DELTAIDXEXTIMEKA;
699 Int_t idxextimePr = ((6*i+j)*NIDX)+DELTAIDXEXTIMEPR;
700 AliDebug(2,Form("idxtime = %i, idxextimePi = %i, idxextimeKa = %i, idxextimePr = %i", idxtime, idxextimePi, idxextimeKa, idxextimePr));
1ed77580 701 Double_t time=smallarray[idxtime]; // TOF time in ns
702 timeofflight[j]=time+t0offset;
703 exptof[j][0]=smallarray[idxextimePi];
704 exptof[j][1]=smallarray[idxextimeKa];
705 exptof[j][2]=smallarray[idxextimePr];
706 AliDebug(2,Form("j = %i, Time = %f, and Exp time in PID: pi = %f, K = %f, p = %f",j,timeofflight[j],exptof[j][0],exptof[j][1],exptof[j][2]));
707 }
e26353e8 708
709
1ed77580 710 Float_t chisquarebest=999.;
e26353e8 711 AliInfo(Form(" Set = %i with %i tracks ", i,ntrkinset));
712 chisquarebest = LoopCombPID(ntrkinset, ntrkinset,exptof,&texp[0],&timeofflight[0], &index[0],chisquarebest);
713
1ed77580 714 Float_t confLevel=999;
715 if(chisquarebest<999.){
716 Double_t dblechisquare=(Double_t)chisquarebest;
e26353e8 717 confLevel=(Float_t)TMath::Prob(dblechisquare,ntrkinset-1);
1ed77580 718 }
e26353e8 719
1ed77580 720 // assume they are all pions for fake sets
721 if(confLevel<0.01 || confLevel==999. ){
e26353e8 722 for (Int_t itrk=0; itrk<ntrkinset; itrk++)fassparticle[itrk]=0;
1ed77580 723 }
e26353e8 724
725 AliDebug(2,Form(" Best Assignment, selection %i %i %i %i %i %i %i %i %i %i",fassparticle[0],fassparticle[1],fassparticle[2],fassparticle[3],fassparticle[4],fassparticle[5],fassparticle[6],fassparticle[7],fassparticle[8],fassparticle[9]));
1ed77580 726
e26353e8 727 for (Int_t kk=0;kk<ntrkinset;kk++){
728 Int_t idxextimePi = ((6*i+kk)*NIDX)+DELTAIDXEXTIMEPI;
729 Int_t idxextimeKa = ((6*i+kk)*NIDX)+DELTAIDXEXTIMEKA;
730 Int_t idxextimePr = ((6*i+kk)*NIDX)+DELTAIDXEXTIMEPR;
1ed77580 731 // storing in third slot of smallarray the assigned expected time
e26353e8 732 fhPID->Fill(fassparticle[kk]);
733 // fassparticle[kk]=0; //assuming all particles are pions
734 if (fassparticle[kk]==0){ //assigned to be a Pi
1ed77580 735 smallarray[idxextimeKa]=-1;
736 smallarray[idxextimePr]=-1;
737 }
e26353e8 738 else if (fassparticle[kk]==1){ //assigned to be a Ka
1ed77580 739 smallarray[idxextimePi]=smallarray[idxextimeKa];
740 smallarray[idxextimeKa]=-1;
741 smallarray[idxextimePr]=-1;
742 }
e26353e8 743 else if (fassparticle[kk]==2){ //assigned to be a Pr
1ed77580 744 smallarray[idxextimePi]=smallarray[idxextimePr];
745 smallarray[idxextimeKa]=-1;
746 smallarray[idxextimePr]=-1;
747 }
748 }
ea23eb9f 749 delete[] exptof;
750 delete[] timeofflight;
751 delete[] texp;
752 delete[] index;
e26353e8 753
1ed77580 754 }
755 return 1;
756}
e26353e8 757//---------------------------------------------------------------------
758Float_t AliTOFCalibTask::LoopCombPID(Int_t itrkinset, Int_t ntrkinset, Float_t **exptof, Float_t *texp, Float_t *timeofflight, Int_t *index, Float_t chisquarebest){
759
760 Int_t indextr = ntrkinset-itrkinset;
761
762 for (index[indextr]=0;index[indextr]<3;index[indextr]++){
763 Int_t ii = index[indextr];
764 texp[indextr]=exptof[indextr][ii];
765 if (indextr<ntrkinset-1){
766 chisquarebest = LoopCombPID(ntrkinset-indextr-1,ntrkinset,exptof,&texp[0],&timeofflight[0],&index[0], chisquarebest);
1ed77580 767 }
e26353e8 768
769 else {
770 Float_t *besttimezero = new Float_t[ntrkinset];
771 Float_t *bestchisquare = new Float_t[ntrkinset];
772 Float_t *bestweightedtimezero = new Float_t[ntrkinset];
773 Float_t sumAllweights=0.;
774 Float_t meantzero=0.;
775 Float_t *weightedtimezero = new Float_t[ntrkinset];
776 Float_t *timezero = new Float_t[ntrkinset];
777
778 AliDebug(2,Form("Configuration = %i, %i, %i, %i, %i, %i, %i, %i, %i, %i, so far chisquarebest = %f ",index[0],index[1],index[2],index[3],index[4],index[5],index[6],index[7],index[8],index[9],chisquarebest));
779 for (Int_t itz=0; itz<ntrkinset;itz++) {
780 timezero[itz]=texp[itz]-timeofflight[itz];
781 weightedtimezero[itz]=timezero[itz]/TRACKERROR;
782 sumAllweights+=1./TRACKERROR;
783 meantzero+=weightedtimezero[itz];
784 } // end loop for (Int_t itz=0; itz<ntrkinset;itz++)
785
786 meantzero=meantzero/sumAllweights; // it is given in [ns]
787
788 // calculate chisquare
789
790 Float_t chisquare=0.;
791 for (Int_t icsq=0; icsq<ntrkinset;icsq++) {
792 chisquare+=(timezero[icsq]-meantzero)*(timezero[icsq]-meantzero)/TRACKERROR;
793 } // end loop for (Int_t icsq=0; icsq<ntrkinset;icsq++)
794
795 Int_t npion=0;
796 for (Int_t j=0;j<ntrkinset;j++){
797 if(index[j]==0)npion++;
1ed77580 798 }
e26353e8 799
800 if(chisquare<=chisquarebest && ((Float_t) npion/ ((Float_t) ntrkinset)>0.3)){
801 for(Int_t iqsq = 0; iqsq<ntrkinset; iqsq++) {
802 besttimezero[iqsq]=timezero[iqsq];
803 bestweightedtimezero[iqsq]=weightedtimezero[iqsq];
804 bestchisquare[iqsq]=(timezero[iqsq]-meantzero)*(timezero[iqsq]-meantzero)/TRACKERROR;
805 }
806
807 for (Int_t j=0;j<ntrkinset;j++){
808 fassparticle[j]=index[j];
809 }
810
811 chisquarebest=chisquare;
1ed77580 812 }
ea23eb9f 813 delete[] besttimezero;
814 delete[] bestchisquare;
815 delete[] bestweightedtimezero;
816 delete[] weightedtimezero;
817 delete[] timezero;
1ed77580 818 }
819 }
e26353e8 820 return chisquarebest;
1ed77580 821}