]> 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
AliMC::BeginEvent gets run number from CDB manager
[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");
196 delete[] fbigarray;
ea23eb9f 197 delete[] findexarray;
1ed77580 198 delete fOutputContainer;
199 delete fhToT;
200 delete fhTime;
201 delete fhExpTimePi;
202 delete fhExpTimeKa;
203 delete fhExpTimePr;
204 delete fhPID;
205 delete fhch;
206}
207//______________________________________________________________________________
208AliTOFCalibTask& AliTOFCalibTask::operator=(const AliTOFCalibTask &calibtask)
209{
210 //assignment operator
211 this->fdir=calibtask.fdir;
212 this->fChain=calibtask.fChain;
213 this->fESD=calibtask.fESD;
1ed77580 214 this->fToT=calibtask.fToT;
215 this->fTime=calibtask.fTime;
216 this->fExpTimePi=calibtask.fExpTimePi;
217 this->fExpTimeKa=calibtask.fExpTimeKa;
218 this->fExpTimePr=calibtask.fExpTimePr;
1ed77580 219 this->fMinTime=calibtask.fMinTime;
1ed77580 220 this->fnESD=calibtask.fnESD;
221 this->fnESDselected=calibtask.fnESDselected;
222 this->fnESDkTOFout=calibtask.fnESDkTOFout;
223 this->fnESDkTIME=calibtask.fnESDkTIME;
224 this->fnESDassTOFcl=calibtask.fnESDassTOFcl;
225 this->fnESDTIMEcut=calibtask.fnESDTIMEcut;
226 this->fnESDTRDcut=calibtask.fnESDTRDcut;
227 this->fhToT=calibtask.fhToT;
228 this->fhTime=calibtask.fhTime;
229 this->fhExpTimePi=calibtask.fhExpTimePi;
230 this->fhExpTimeKa=calibtask.fhExpTimeKa;
231 this->fhExpTimePr=calibtask.fhExpTimePr;
232 this->fOutputContainer=calibtask.fOutputContainer;
233 this->fhPID=calibtask.fhPID;
234 this->fhch=calibtask.fhch;
70ddd0ee 235
236 this->fbigarray = new Float_t*[TOFCHANNELS];
237 this->findexarray = new Int_t[TOFCHANNELS];
e26353e8 238
239 for (Int_t i=0;i<11;i++){
240 this->fassparticle[i]=calibtask.fassparticle[i];
241 }
70ddd0ee 242 for (Int_t i =0;i<TOFCHANNELS;i++){
243 this->fbigarray[i] = new Float_t[CHENTRIES];
244 this->findexarray[i]=calibtask.findexarray[i];
245 for (Int_t j =0;j<CHENTRIES;j++){
246 this->fbigarray[i][j]=calibtask.fbigarray[i][j];
247 }
248 }
1ed77580 249 return *this;
250}
251//--------------------------------------------------------------------------
252void AliTOFCalibTask::BookHistos(){
253
254 //booking histos
255
256 AliInfo(Form("*** Booking Histograms %s", GetName())) ;
257
258 fhToT=
259 new TH1F("hToT", " ToT distribution (ns) ", 400, 0, 40);
260 fhTime=
261 new TH1F("hTime", " Time distribution (ns) ", 400, 0, 40);
262 fhExpTimePi=
263 new TH1F("hExpTimePi", " Exp Time distribution, Pi (ns) ", 400, 0, 40);
264 fhExpTimeKa=
265 new TH1F("hExpTimeKa", " Exp Time distribution, Ka (ns) ", 400, 0, 40);
266 fhExpTimePr=
267 new TH1F("hExpTimePr", " Exp Time distribution, Pr (ns) ", 400, 0, 40);
268 fhPID=
269 new TH1I("hPID", " Combinatorial PID identity ", 3, 0, 3);
270 fhch=
271 new TH1D("hch", " TOF channel ", TOFCHANNELS, 0, TOFCHANNELS);
272
273 //create the putput container
274 fOutputContainer = new TObjArray(7) ;
275 fOutputContainer->SetName(GetName()) ;
276
277 fOutputContainer->AddAt(fhToT, 0) ;
278 fOutputContainer->AddAt(fhTime, 1) ;
279 fOutputContainer->AddAt(fhExpTimePi, 2) ;
280 fOutputContainer->AddAt(fhExpTimeKa, 3) ;
281 fOutputContainer->AddAt(fhExpTimePr, 4) ;
282 fOutputContainer->AddAt(fhPID, 5) ;
283 fOutputContainer->AddAt(fhch, 6) ;
284
285}
286//----------------------------------------------------------------------------
287void AliTOFCalibTask::DrawHistos(){
288
289 // drawing output histos
290
291 AliInfo(Form("*** Drawing Histograms %s", GetName())) ;
292
293 TCanvas * canvasToTTime = new TCanvas("canvasToTTime", " ToT and Time ",400, 30, 550, 630) ;
294 canvasToTTime->Divide(1,2);
295 canvasToTTime->cd(1);
296 fhToT->SetLineColor(4);
297 fhToT->GetXaxis()->SetTitle("ToT (ns)");
298 fhToT->Draw("hist");
299 canvasToTTime->cd(2);
300 fhTime->SetLineColor(4);
301 fhTime->GetXaxis()->SetTitle("Time (ns)");
302 fhTime->Draw("hist");
303 canvasToTTime->Update();
304 canvasToTTime->Print("ToTTime.gif");
305
306 TCanvas * canvasExpTime = new TCanvas("canvasExpTime", " Expected Times ",400, 30, 550, 630) ;
307 canvasExpTime->Divide(1,3);
308 canvasExpTime->cd(1);
309 fhExpTimePi->SetLineColor(4);
310 fhExpTimePi->GetXaxis()->SetTitle("Exp Time (ns), #pi");
311 fhExpTimePi->Draw("hist");
312 canvasExpTime->cd(2);
313 fhExpTimeKa->SetLineColor(4);
314 fhExpTimeKa->GetXaxis()->SetTitle("Exp Time (ns), K");
315 fhExpTimeKa->Draw("hist");
316 canvasExpTime->cd(3);
317 fhExpTimePr->SetLineColor(4);
318 fhExpTimePr->GetXaxis()->SetTitle("Exp Time (ns), p");
319 fhExpTimePr->Draw("hist");
320
321 canvasExpTime->Print("ExpTime.gif");
322
323 TCanvas * canvasPID = new TCanvas("canvasPID", " Combinatorial PID ",400, 30, 550, 400);
324 fhPID->GetXaxis()->SetTitle("Comb PID");
325 fhPID->GetXaxis()->SetBinLabel(1,"#pi");
326 fhPID->GetXaxis()->SetBinLabel(2,"K");
327 fhPID->GetXaxis()->SetBinLabel(3,"p");
328 fhPID->Draw("hist");
329
330 canvasPID->Print("PID.gif");
331
332 TCanvas * canvasrndch = new TCanvas("canvasrndch", " TOF channel ",400, 30, 550, 400);
333 fhch->GetXaxis()->SetTitle("TOF ch");
334 fhch->Draw("hist");
335 Float_t meanTOFch = 0;
336 for (Int_t ibin=0;ibin<TOFCHANNELS;ibin++){
337 meanTOFch+=(Float_t)fhch->GetBinContent(ibin+1);
338 }
339
340 meanTOFch/=TOFCHANNELS;
341 AliDebug(1,Form(" Mean number of tracks/channel = %f ",meanTOFch));
342
343 canvasrndch->Print("rndch.gif");
344
345 char line[1024] ;
346 sprintf(line, ".!tar -zcvf %s.tar.gz *.gif", GetName()) ;
347 gROOT->ProcessLine(line);
348 sprintf(line, ".!rm -fR *.gif");
349 gROOT->ProcessLine(line);
350 AliInfo(Form("*** TOF Calib Task: plots saved in %s.tar.gz...\n", GetName())) ;
351}
352
353//______________________________________________________________________________
354void AliTOFCalibTask::ConnectInputData(const Option_t*)
355{
356 // Initialization of branch container and histograms
357
358 // AliLog::SetClassDebugLevel("AliTOFCalibTask",1);
359 AliInfo(Form("*** Initialization of %s", GetName())) ;
360
361 // Get input data
362 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
363 if (!fChain) {
364 AliError(Form("Input 0 for %s not found\n", GetName()));
365 return ;
366 }
367
368 // One should first check if the branch address was taken by some other task
369 char ** address = (char **)GetBranchAddress(0, "ESD");
370 if (address) {
70ddd0ee 371 fESD = (AliESDEvent*)(*address);
1ed77580 372 } else {
70ddd0ee 373 fESD = new AliESDEvent();
1ed77580 374 }
1e76f1ec 375 fESD->ReadFromTree(fChain) ;
1ed77580 376
377 BookHistos();
378
379}
70ddd0ee 380//-----------------------------------------------------------------------
381Bool_t AliTOFCalibTask::Notify()
382{
383 // Initialisation of branch container and histograms
384
385 AliInfo(Form("*** We are in %s::Notify()", GetName())) ;
386
387 // Get input data
388 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
389 if (!fChain) {
390 AliError(Form("Input 0 for %s not found\n", GetName()));
391 return kFALSE;
392 }
393
1e76f1ec 394 char ** address = (char **)GetBranchAddress(0, "ESD");
395 if (address) {
396 fESD = (AliESDEvent*)(*address);
397 } else {
398 fESD = new AliESDEvent();
399 }
e26353e8 400 fESD->ReadFromTree(fChain) ;
70ddd0ee 401
402 return kTRUE;
403}
404
1ed77580 405//________________________________________________________________________
406void AliTOFCalibTask::CreateOutputObjects()
407{
408// Create histograms
409}
410
411//______________________________________________________________________________
412void AliTOFCalibTask::Exec(Option_t * opt)
413{
414
415 // main
416
417 AliInfo(Form("*** Executing %s", GetName())) ;
418
419//******* The loop over events -----------------------------------------------
420
421 // Processing of one event
422 Long64_t entry = fChain->GetReadEntry() ;
423 if (!fESD) {
424 AliError("fESD is not connected to the input!") ;
425 return ;
426 }
427
428 if ( !((entry-1)%100) )
429 AliDebug(1,Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
430
431 // ************************ TOF *************************************
432
433 fMinTime=22E3; //ns; not used
434 Int_t ntrk = fESD->GetNumberOfTracks() ;
435 fnESD+=ntrk;
436 Int_t nselected = 0;
437 Int_t itr = -1;
438 while ( ntrk-- ) {
439 itr++;
440 AliESDtrack * t = fESD->GetTrack(ntrk) ;
441 //selecting only good quality tracks
442 if (!Select(t)) continue;
443 nselected++;
444 Int_t ich = Int_t(t->GetTOFCalChannel());
445 fhch->Fill(ich);
b8047c99 446 // ich=3; //only for debug purpose
1ed77580 447 AliDebug(2,Form(" ESD in channel %i, filling array %i",t->GetTOFCalChannel(),ich));
448 Float_t tot = t->GetTOFsignalToT();
449 Float_t time = t->GetTOFsignalRaw();
e26353e8 450 AliDebug(2,Form(" track # %i in channel %i, time = %f \n",ntrk,ich,time));
1ed77580 451 Double_t expTime[10];
452 t->GetIntegratedTimes(expTime);
453 Float_t expTimePi = expTime[2]*1.E-3;
454 Float_t expTimeKa = expTime[3]*1.E-3;
455 Float_t expTimePr = expTime[4]*1.E-3;
456 if (findexarray[ich]==(Int_t)(CHENTRIES/NIDX)) {
457 AliInfo(Form("too many tracks in channel %i, not storing any more...",ich));
458 continue;
459 }
460 findexarray[ich]++;
461 AliDebug(2,Form("tracks in channel %i = %i, storing... ",ich, findexarray[ich] ));
462 Int_t ientry=(findexarray[ich]-1)*NIDX;
463 fbigarray[ich][ientry+DELTAIDXTOT]=tot; //in ns
464 fbigarray[ich][ientry+DELTAIDXTIME]=time*1E-3; // in ns
465 fbigarray[ich][ientry+DELTAIDXEXTIMEPI]=expTimePi;
466 fbigarray[ich][ientry+DELTAIDXEXTIMEKA]=expTimeKa;
467 fbigarray[ich][ientry+DELTAIDXEXTIMEPR]=expTimePr;
468 fhToT->Fill(fbigarray[ich][ientry+DELTAIDXTOT]);
469 fhTime->Fill(fbigarray[ich][ientry+DELTAIDXTIME]);
470 fhExpTimePi->Fill(fbigarray[ich][ientry+DELTAIDXEXTIMEPI]);
471 fhExpTimeKa->Fill(fbigarray[ich][ientry+DELTAIDXEXTIMEKA]);
472 fhExpTimePr->Fill(fbigarray[ich][ientry+DELTAIDXEXTIMEPR]);
473 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));
474 }
475 fnESDselected+=nselected;
476
477 PostData(0, fOutputContainer);
478}
479
480//_____________________________________________________________________________
481void AliTOFCalibTask::Terminate(Option_t *)
482{
483 // Processing when the event loop is ended
484
485 // some plots
486
487 TH1::AddDirectory(0);
e26353e8 488 TDirectory *dir = gDirectory;
1ed77580 489 AliInfo("TOF Calib Task: End of events loop");
e26353e8 490 for (Int_t ich = 0;ich<TOFCHANNELS;ich++){
491 if (findexarray[ich]>0){
492 Int_t ncuttime = SelectOnTime(&fbigarray[ich][0],findexarray[ich],ich);
493 fnESDselected-=ncuttime;
494 }
495 }
496
1ed77580 497 AliInfo(Form(" Number of analyzed ESD tracks: %i\n",fnESD));
498 AliInfo(Form(" Number of selected ESD tracks: %i\n",fnESDselected));
499 AliInfo(Form(" Number of ESD tracks with kTOFout: %i\n",fnESDkTOFout));
500 AliInfo(Form(" Number of ESD tracks with kTIME: %i\n",fnESDkTIME));
501 AliInfo(Form(" Number of ESD tracks with TRDcut: %i\n",fnESDTRDcut));
502 AliInfo(Form(" Number of ESD tracks with TIMEcut: %i\n",fnESDTIMEcut));
503 AliInfo(Form(" Number of ESD tracks with assTOFcl: %i\n",fnESDassTOFcl));
504
505 for (Int_t i = 0;i<TOFCHANNELS;i++){
506 Int_t size=findexarray[i]*NIDX;
507 AliDebug(2, Form(" entries %i in channel %i ",findexarray[i],i));
e26353e8 508 if (findexarray[i]<=2) {
1ed77580 509 AliDebug(1, Form(" not enough statistics for combined PID for channel %i, putting all the tracks as if they were pions",i));
510 continue;
511 }
512 if (!CombPID(&fbigarray[i][0], size)) AliError("ERROR!!!!ERROR!!!");
513 }
70ddd0ee 514
e26353e8 515
1ed77580 516 DrawHistos();
517
518 // saving data in a tree
519
e26353e8 520 AliInfo("Building tree for Calibration");
521 TTree * tree = new TTree("T", "Tree for TOF calibration");
522 Float_t p[CHENTRIESSMALL];
1ed77580 523 Int_t nentries;
e26353e8 524 tree->Branch("nentries",&nentries,"nentries/I");
525 tree->Branch("TOFentries",p,"TOFentries[nentries]/F");
526 for (Int_t i=0;i<TOFCHANNELS;i++){
527 nentries=findexarray[i]*(NIDXSMALL); // when filling small array,
528 // only first 3 floats taken
529 // into account
530 for (Int_t j=0; j<findexarray[i];j++){
531 for (Int_t k=0; k<NIDXSMALL;k++){
532 Int_t index1= j*NIDXSMALL+k; // index in small array
533 Int_t index2=j*NIDX+k; // index in big array
534 p[index1]=fbigarray[i][index2];
535 }
1ed77580 536 }
e26353e8 537 tree->Fill();
1ed77580 538 }
1ed77580 539
e26353e8 540 AliInfo("Putting tree for calibration in Reference data");
541
542 // grid file option
543 Char_t filename[100];
544 sprintf(filename,"alien:///alice/cern.ch/user/c/czampolli/TOFCalibReference_%i.root",frun);
545 TGrid::Connect("alien://");
546 TFile *filegrid = TFile::Open(filename,"CREATE");
547 tree->Write();
548 dir->cd();
1ed77580 549
1ed77580 550}
1ed77580 551//_____________________________________________________________________________
552
553Bool_t AliTOFCalibTask::Select(AliESDtrack *t){
554
555 // selecting good quality tracks
556 //track selection: reconstrution to TOF:
557 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) {
558 return 0;
559 }
560 fnESDkTOFout++;
561 //IsStartedTimeIntegral
562 if ((t->GetStatus()&AliESDtrack::kTIME)==0) {
563 return 0;
564 }
565 fnESDkTIME++;
566 if (t->GetStatus() & AliESDtrack::kTRDbackup) {
567 Float_t xout= t->GetOuterParam()->GetX();
568 if (xout<364.25 && xout > 300.) return 0;
569 }
570 fnESDTRDcut++;
571 Double_t time=t->GetTOFsignal();
572 time*=1.E-3; // tof given in nanoseconds
573 if(time >= fMinTime){
574 return 0;
575 }
576 fnESDTIMEcut++;
577
578 Double_t mom=t->GetP();
579 if (!(mom<=UPPERMOMBOUND && mom>=LOWERMOMBOUND)){
580 // return 0; // skipping momentum cut
581 }
582
583 UInt_t assignedTOFcluster=t->GetTOFcluster();//index of the assigned TOF cluster, >0 ?
584 if(assignedTOFcluster==0){ // not matched
585 return 0;
586 }
587 fnESDassTOFcl++;
588 return 1;
589}
590//_____________________________________________________________________________
591
e26353e8 592Int_t AliTOFCalibTask::SelectOnTime(Float_t *charray, Int_t ntracks, Int_t ich){
593
594 // discarding tracks with time-mintime < MINTIME
595
596 Int_t ndeleted=0;
597 Float_t mintime = 1E6;
598 for (Int_t itr = 0;itr<ntracks;itr++){
599 Int_t ientry=itr*NIDX;
600 Float_t time = charray[ientry+DELTAIDXTIME];// in ns
601 if (time<mintime) mintime = time;
602 }
603 AliDebug(1,Form("Mintime for channel %i = %f",ich, mintime));
604 for (Int_t itr = 0;itr<ntracks;itr++){
605 Int_t ientry=itr*NIDX;
606 Float_t time = charray[ientry+DELTAIDXTIME];// in ns
607 if ((time-mintime)>MINTIME) {
608 ndeleted++;
609 AliDebug(1,Form("Deleting %i track from channel %i, time = %f",ndeleted, ich, time));
610 findexarray[ich]--;
611 for (Int_t j=itr+1;j<ntracks;j++){
612 Int_t ientry=j*NIDX;
613 Int_t idxtot = ientry+DELTAIDXTOT;
614 Int_t idxtime = ientry+DELTAIDXTIME;
615 Int_t idxextimePi = ientry+DELTAIDXEXTIMEPI;
616 Int_t idxextimeKa = ientry+DELTAIDXEXTIMEKA;
617 Int_t idxextimePr = ientry+DELTAIDXEXTIMEPR;
618 Int_t ientrydel=(j-1)*NIDX;
619 Int_t idxtotdel = ientrydel+DELTAIDXTOT;
620 Int_t idxtimedel = ientrydel+DELTAIDXTIME;
621 Int_t idxextimePidel = ientrydel+DELTAIDXEXTIMEPI;
622 Int_t idxextimeKadel = ientrydel+DELTAIDXEXTIMEKA;
623 Int_t idxextimePrdel = ientrydel+DELTAIDXEXTIMEPR;
624 charray[idxtotdel]=charray[idxtot];
625 charray[idxtimedel]=charray[idxtime];
626 charray[idxextimePidel]=charray[idxextimePi];
627 charray[idxextimeKadel]=charray[idxextimeKa];
628 charray[idxextimePrdel]=charray[idxextimePr];
629 }
630 }
631 }
632 return ndeleted;
633}
634//_____________________________________________________________________________
635
1ed77580 636Bool_t AliTOFCalibTask::CombPID(Float_t *smallarray, Int_t size){
637
e26353e8 638 // track Combinatorial PID for calibration, only when more than 2 particles
639 // fall in channel
1ed77580 640
1ed77580 641 Float_t t0offset=0;
e26353e8 642
643 if (size/NIDX<=2){
644 AliDebug(1,"Number of tracks in channel smaller than 2, identifying every particle as if it was a pion!");
1ed77580 645 return 0;
646 }
e26353e8 647
648 Int_t ntracksinchannel = size/NIDX;
649 Int_t ntrkinset=-1;
650 Int_t nset = ntracksinchannel/6;
651
652 if (nset ==0) {
653 nset=1;
654 if (ntracksinchannel < 6) {
655 AliDebug(2,"Number of tracks in set smaller than 6, Combinatorial PID still applied to this set.");
656 }
1ed77580 657 }
658
e26353e8 659 AliInfo(Form("number of sets = %i", nset));
660 // loop on sets
1ed77580 661 for (Int_t i=0; i< nset; i++) {
e26353e8 662 if (i<nset-1)ntrkinset=6;
663 else {
664 if (ntracksinchannel<6){
665 ntrkinset=ntracksinchannel;
666 }
667 else{
668 ntrkinset=6+Int_t((ntracksinchannel)%6);
669 }
670 }
671 AliInfo(Form("set = %i, number of tracks in set = %i", i,ntrkinset));
672
673 for (Int_t ii=0;ii<ntrkinset;ii++){
674 fassparticle[ii]=-1;
675 }
676 Float_t **exptof;
677 Float_t* timeofflight;
678 Float_t* texp;
679 Int_t* index;
680
681 exptof = new Float_t*[ntrkinset];
682 timeofflight = new Float_t[ntrkinset];
683 texp = new Float_t[ntrkinset];
684 index = new Int_t[ntrkinset];
685 for (Int_t ii=0;ii<ntrkinset;ii++){
686 exptof[ii] = new Float_t[3];
687 timeofflight[ii]=0;
688 texp[ii]=0;
689 index[ii]=-1;
690 }
691
692 for (Int_t j=0; j<ntrkinset; j++) {
693 Int_t idxtime = ((6*i+j)*NIDX)+DELTAIDXTIME;
694 Int_t idxextimePi = ((6*i+j)*NIDX)+DELTAIDXEXTIMEPI;
695 Int_t idxextimeKa = ((6*i+j)*NIDX)+DELTAIDXEXTIMEKA;
696 Int_t idxextimePr = ((6*i+j)*NIDX)+DELTAIDXEXTIMEPR;
697 AliDebug(2,Form("idxtime = %i, idxextimePi = %i, idxextimeKa = %i, idxextimePr = %i", idxtime, idxextimePi, idxextimeKa, idxextimePr));
1ed77580 698 Double_t time=smallarray[idxtime]; // TOF time in ns
699 timeofflight[j]=time+t0offset;
700 exptof[j][0]=smallarray[idxextimePi];
701 exptof[j][1]=smallarray[idxextimeKa];
702 exptof[j][2]=smallarray[idxextimePr];
703 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]));
704 }
e26353e8 705
706
1ed77580 707 Float_t chisquarebest=999.;
e26353e8 708 AliInfo(Form(" Set = %i with %i tracks ", i,ntrkinset));
709 chisquarebest = LoopCombPID(ntrkinset, ntrkinset,exptof,&texp[0],&timeofflight[0], &index[0],chisquarebest);
710
1ed77580 711 Float_t confLevel=999;
712 if(chisquarebest<999.){
713 Double_t dblechisquare=(Double_t)chisquarebest;
e26353e8 714 confLevel=(Float_t)TMath::Prob(dblechisquare,ntrkinset-1);
1ed77580 715 }
e26353e8 716
1ed77580 717 // assume they are all pions for fake sets
718 if(confLevel<0.01 || confLevel==999. ){
e26353e8 719 for (Int_t itrk=0; itrk<ntrkinset; itrk++)fassparticle[itrk]=0;
1ed77580 720 }
e26353e8 721
722 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 723
e26353e8 724 for (Int_t kk=0;kk<ntrkinset;kk++){
725 Int_t idxextimePi = ((6*i+kk)*NIDX)+DELTAIDXEXTIMEPI;
726 Int_t idxextimeKa = ((6*i+kk)*NIDX)+DELTAIDXEXTIMEKA;
727 Int_t idxextimePr = ((6*i+kk)*NIDX)+DELTAIDXEXTIMEPR;
1ed77580 728 // storing in third slot of smallarray the assigned expected time
e26353e8 729 fhPID->Fill(fassparticle[kk]);
730 // fassparticle[kk]=0; //assuming all particles are pions
731 if (fassparticle[kk]==0){ //assigned to be a Pi
1ed77580 732 smallarray[idxextimeKa]=-1;
733 smallarray[idxextimePr]=-1;
734 }
e26353e8 735 else if (fassparticle[kk]==1){ //assigned to be a Ka
1ed77580 736 smallarray[idxextimePi]=smallarray[idxextimeKa];
737 smallarray[idxextimeKa]=-1;
738 smallarray[idxextimePr]=-1;
739 }
e26353e8 740 else if (fassparticle[kk]==2){ //assigned to be a Pr
1ed77580 741 smallarray[idxextimePi]=smallarray[idxextimePr];
742 smallarray[idxextimeKa]=-1;
743 smallarray[idxextimePr]=-1;
744 }
745 }
ea23eb9f 746 delete[] exptof;
747 delete[] timeofflight;
748 delete[] texp;
749 delete[] index;
e26353e8 750
1ed77580 751 }
752 return 1;
753}
e26353e8 754//---------------------------------------------------------------------
755Float_t AliTOFCalibTask::LoopCombPID(Int_t itrkinset, Int_t ntrkinset, Float_t **exptof, Float_t *texp, Float_t *timeofflight, Int_t *index, Float_t chisquarebest){
756
757 Int_t indextr = ntrkinset-itrkinset;
758
759 for (index[indextr]=0;index[indextr]<3;index[indextr]++){
760 Int_t ii = index[indextr];
761 texp[indextr]=exptof[indextr][ii];
762 if (indextr<ntrkinset-1){
763 chisquarebest = LoopCombPID(ntrkinset-indextr-1,ntrkinset,exptof,&texp[0],&timeofflight[0],&index[0], chisquarebest);
1ed77580 764 }
e26353e8 765
766 else {
767 Float_t *besttimezero = new Float_t[ntrkinset];
768 Float_t *bestchisquare = new Float_t[ntrkinset];
769 Float_t *bestweightedtimezero = new Float_t[ntrkinset];
770 Float_t sumAllweights=0.;
771 Float_t meantzero=0.;
772 Float_t *weightedtimezero = new Float_t[ntrkinset];
773 Float_t *timezero = new Float_t[ntrkinset];
774
775 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));
776 for (Int_t itz=0; itz<ntrkinset;itz++) {
777 timezero[itz]=texp[itz]-timeofflight[itz];
778 weightedtimezero[itz]=timezero[itz]/TRACKERROR;
779 sumAllweights+=1./TRACKERROR;
780 meantzero+=weightedtimezero[itz];
781 } // end loop for (Int_t itz=0; itz<ntrkinset;itz++)
782
783 meantzero=meantzero/sumAllweights; // it is given in [ns]
784
785 // calculate chisquare
786
787 Float_t chisquare=0.;
788 for (Int_t icsq=0; icsq<ntrkinset;icsq++) {
789 chisquare+=(timezero[icsq]-meantzero)*(timezero[icsq]-meantzero)/TRACKERROR;
790 } // end loop for (Int_t icsq=0; icsq<ntrkinset;icsq++)
791
792 Int_t npion=0;
793 for (Int_t j=0;j<ntrkinset;j++){
794 if(index[j]==0)npion++;
1ed77580 795 }
e26353e8 796
797 if(chisquare<=chisquarebest && ((Float_t) npion/ ((Float_t) ntrkinset)>0.3)){
798 for(Int_t iqsq = 0; iqsq<ntrkinset; iqsq++) {
799 besttimezero[iqsq]=timezero[iqsq];
800 bestweightedtimezero[iqsq]=weightedtimezero[iqsq];
801 bestchisquare[iqsq]=(timezero[iqsq]-meantzero)*(timezero[iqsq]-meantzero)/TRACKERROR;
802 }
803
804 for (Int_t j=0;j<ntrkinset;j++){
805 fassparticle[j]=index[j];
806 }
807
808 chisquarebest=chisquare;
1ed77580 809 }
ea23eb9f 810 delete[] besttimezero;
811 delete[] bestchisquare;
812 delete[] bestweightedtimezero;
813 delete[] weightedtimezero;
814 delete[] timezero;
1ed77580 815 }
816 }
e26353e8 817 return chisquarebest;
1ed77580 818}