TOF/AliTOFCalibTask.cxx

   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 $Log$
  17 Revision 1.6  2007/10/05 10:44:11  zampolli
  18 Oversight for debugging purpose fixed
  19
  20 Revision 1.5  2007/10/05 10:38:10  zampolli
  21 Oversight fixed
  22
  23 Revision 1.4  2007/10/04 15:36:44  zampolli
  24 Updates to new TOF offline calibration schema
  25
  26 Revision 1.3  2007/07/31 07:26:16  zampolli
  27 Bug fixed in the $ field
  28
  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>
  50 #include <TDirectory.h>
  51 #include <TGrid.h>
  52
  53 #include "AliTOFCalibTask.h"
  54 #include "AliESDEvent.h"
  55 #include "AliESDtrack.h"
  56 #include "AliLog.h"
  57 #include "AliCDBManager.h"
  58 #include "AliCDBMetaData.h"
  59 #include "AliCDBId.h"
  60 #include "AliCDBStorage.h"
  61
  62 //______________________________________________________________________________
  63 AliTOFCalibTask::AliTOFCalibTask(const char *name) :
  64   AliAnalysisTask(name,""),
  65   fdir(0),
  66   fChain(0),
  67   fESD(0),
  68   fToT(0),
  69   fTime(0),
  70   fExpTimePi(0),
  71   fExpTimeKa(0),
  72   fExpTimePr(0),
  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),
  90   fOutputContainer(0),
  91   frun(0)
  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();
 100   fbigarray = new Float_t*[TOFCHANNELS];
 101   findexarray = new Int_t[TOFCHANNELS];
 102
 103   for (Int_t i=0;i<11;i++){
 104     fassparticle[i]=-1;
 105   }
 106
 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   }
 114 }
 115
 116 //______________________________________________________________________________
 117 AliTOFCalibTask::AliTOFCalibTask(const AliTOFCalibTask &calibtask) :
 118   AliAnalysisTask("AliTOFCalibTask",""),
 119   fdir(0),
 120   fChain(0),
 121   fESD(0),
 122   fToT(0),
 123   fTime(0),
 124   fExpTimePi(0),
 125   fExpTimeKa(0),
 126   fExpTimePr(0),
 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),
 144   fOutputContainer(0),
 145   frun(0)
 146 {
 147   // Copy Constructor.
 148   fdir=calibtask.fdir;
 149   fChain=calibtask.fChain;
 150   fESD=calibtask.fESD;
 151   fToT=calibtask.fToT;
 152   fTime=calibtask.fTime;
 153   fExpTimePi=calibtask.fExpTimePi;
 154   fExpTimeKa=calibtask.fExpTimeKa;
 155   fExpTimePr=calibtask.fExpTimePr;
 156   fMinTime=calibtask.fMinTime;
 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;
 171   frun=calibtask.frun;
 172   fOutputContainer=calibtask.fOutputContainer;
 173
 174   fbigarray = new Float_t*[TOFCHANNELS];
 175   findexarray = new Int_t[TOFCHANNELS];
 176
 177   for (Int_t i=0;i<11;i++){
 178     fassparticle[i]=calibtask.fassparticle[i];
 179   }
 180
 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
 189 }
 190 //______________________________________________________________________________
 191 AliTOFCalibTask:: ~AliTOFCalibTask()
 192 {
 193   // destructor
 194
 195   AliInfo("TOF Calib Task: Deleting");
 196   delete[] fbigarray;
 197   delete[] findexarray;
 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 //______________________________________________________________________________
 208 AliTOFCalibTask& AliTOFCalibTask::operator=(const AliTOFCalibTask &calibtask)
 209 {
 210    //assignment operator
 211   this->fdir=calibtask.fdir;
 212   this->fChain=calibtask.fChain;
 213   this->fESD=calibtask.fESD;
 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;
 219   this->fMinTime=calibtask.fMinTime;
 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;
 235
 236   this->fbigarray = new Float_t*[TOFCHANNELS];
 237   this->findexarray = new Int_t[TOFCHANNELS];
 238
 239   for (Int_t i=0;i<11;i++){
 240     this->fassparticle[i]=calibtask.fassparticle[i];
 241   }
 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   }
 249   return *this;
 250 }
 251 //--------------------------------------------------------------------------
 252 void 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 //----------------------------------------------------------------------------
 287 void 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 //______________________________________________________________________________
 354 void 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) {
 371     fESD = (AliESDEvent*)(*address);
 372   } else {
 373     fESD = new AliESDEvent();
 374   }
 375   fESD->ReadFromTree(fChain) ;
 376
 377   BookHistos();
 378
 379 }
 380 //-----------------------------------------------------------------------
 381 Bool_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
 394   char ** address = (char **)GetBranchAddress(0, "ESD");
 395   if (address) {
 396     fESD = (AliESDEvent*)(*address);
 397   } else {
 398     fESD = new AliESDEvent();
 399   }
 400   fESD->ReadFromTree(fChain) ;
 401
 402   return kTRUE;
 403 }
 404
 405 //________________________________________________________________________
 406 void AliTOFCalibTask::CreateOutputObjects()
 407 {
 408 // Create histograms
 409 }
 410
 411 //______________________________________________________________________________
 412 void 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);
 446     //    ich=3; //only for debug purpose
 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();
 450     AliDebug(2,Form(" track # %i in channel %i, time = %f \n",ntrk,ich,time));
 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 //_____________________________________________________________________________
 481 void AliTOFCalibTask::Terminate(Option_t *)
 482 {
 483   // Processing when the event loop is ended
 484
 485   // some plots
 486
 487   TH1::AddDirectory(0);
 488   TDirectory *dir = gDirectory;
 489   AliInfo("TOF Calib Task: End of events loop");
 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
 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));
 508     if (findexarray[i]<=2) {
 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   }
 514
 515
 516   DrawHistos();
 517
 518   // saving data in a tree
 519
 520   AliInfo("Building tree for Calibration");
 521   TTree * tree = new TTree("T", "Tree for TOF calibration");
 522   Float_t p[CHENTRIESSMALL];
 523   Int_t nentries;
 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       }
 536     }
 537     tree->Fill();
 538   }
 539
 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();
 549
 550 }
 551 //_____________________________________________________________________________
 552
 553 Bool_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
 592 Int_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
 636 Bool_t AliTOFCalibTask::CombPID(Float_t *smallarray, Int_t size){
 637
 638   // track Combinatorial PID for calibration, only when more than 2 particles
 639   // fall in channel
 640
 641   Float_t t0offset=0;
 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!");
 645     return 0;
 646   }
 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     }
 657   }
 658
 659   AliInfo(Form("number of sets = %i", nset));
 660   // loop on sets
 661   for (Int_t i=0; i< nset; i++) {
 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));
 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     }
 705
 706
 707     Float_t chisquarebest=999.;
 708     AliInfo(Form(" Set = %i with %i tracks ", i,ntrkinset));
 709     chisquarebest = LoopCombPID(ntrkinset, ntrkinset,exptof,&texp[0],&timeofflight[0], &index[0],chisquarebest);
 710
 711     Float_t confLevel=999;
 712     if(chisquarebest<999.){
 713       Double_t dblechisquare=(Double_t)chisquarebest;
 714       confLevel=(Float_t)TMath::Prob(dblechisquare,ntrkinset-1);
 715     }
 716
 717     // assume they are all pions for fake sets
 718     if(confLevel<0.01 || confLevel==999. ){
 719       for (Int_t itrk=0; itrk<ntrkinset; itrk++)fassparticle[itrk]=0;
 720     }
 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]));
 723
 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;
 728       // storing in third slot of smallarray the assigned expected time
 729       fhPID->Fill(fassparticle[kk]);
 730       //      fassparticle[kk]=0;  //assuming all particles are pions
 731       if (fassparticle[kk]==0){       //assigned to be a Pi
 732         smallarray[idxextimeKa]=-1;
 733         smallarray[idxextimePr]=-1;
 734       }
 735       else if (fassparticle[kk]==1){  //assigned to be a Ka
 736         smallarray[idxextimePi]=smallarray[idxextimeKa];
 737         smallarray[idxextimeKa]=-1;
 738         smallarray[idxextimePr]=-1;
 739       }
 740       else if (fassparticle[kk]==2){  //assigned to be a Pr
 741         smallarray[idxextimePi]=smallarray[idxextimePr];
 742         smallarray[idxextimeKa]=-1;
 743         smallarray[idxextimePr]=-1;
 744       }
 745     }
 746     delete[] exptof;
 747     delete[] timeofflight;
 748     delete[] texp;
 749     delete[] index;
 750
 751   }
 752   return 1;
 753 }
 754 //---------------------------------------------------------------------
 755 Float_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);
 764     }
 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++;
 795       }
 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;
 809       }
 810       delete[] besttimezero;
 811       delete[] bestchisquare;
 812       delete[] bestweightedtimezero;
 813       delete[] weightedtimezero;
 814       delete[] timezero;
 815     }
 816   }
 817   return chisquarebest;
 818 }