AliPhysics  vAN-20150723 (baea2bf)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnalysisTaskEMCALTimeCalib.cxx
Go to the documentation of this file.
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 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 #include <TChain.h>
17 #include <TTree.h>
18 #include <TFile.h>
19 #include <TH1F.h>
20 #include <TH2D.h>
21 #include <TList.h>
22 #include <TCanvas.h>
23 #include <TGeoManager.h>
24 #include <TRefArray.h>
25 
26 #include "AliLog.h"
27 #include "AliAnalysisTask.h"
28 #include "AliAnalysisManager.h"
29 #include "AliESDEvent.h"
30 #include "AliAODEvent.h"
31 #include "AliVEvent.h"
32 #include "AliESDInputHandler.h"
33 #include "AliAODInputHandler.h"
34 #include "AliESDpid.h"
35 #include "AliTOFcalib.h"
36 #include "AliCDBManager.h"
37 #include "AliRunTag.h"
38 
39 #include "AliTOFT0maker.h"
40 #include "AliVCluster.h"
41 #include "AliESDCaloCluster.h"
42 #include "AliVCaloCells.h"
43 #include "AliESDCaloCells.h"
44 #include "AliAODCaloCluster.h"
45 #include "AliAODCaloCells.h"
46 #include "AliEMCALGeometry.h"
47 
49 
53 
54 using std::cout;
55 using std::endl;
56 
57 //________________________________________________________________________
60 : AliAnalysisTaskSE(name),
61  fRunNumber(-1),
62  fTOFmaker(0),
63  fOutputList(0x0),
64  fgeom(0),
65  fGeometryName(),
66  fMinClusterEnergy(0),
67  fMaxClusterEnergy(0),
68  fMinNcells(0),
69  fMaxNcells(0),
70  fMinLambda0(0),
71  fMaxLambda0(0),
72  fMinLambda0LG(0),
73  fMaxLambda0LG(0),
74  fMaxRtrack(0),
75  fMinCellEnergy(0),
76  fReferenceFileName(),
77  fPileupFromSPD(kFALSE),
78  fMinTime(0),
79  fMaxTime(0),
80  fhcalcEvtTime(0),
81  fhEvtTimeHeader(0),
82  fhEvtTimeDiff(0),
83  fhEventType(0),
84  fhTOFT0vsEventNumber(0),
85  fhTcellvsTOFT0(0),
86  fhTcellvsTOFT0HD(0),
87  fhTcellvsSM(0),
88  fhEneVsAbsIdHG(0),
89  fhEneVsAbsIdLG(0),
90  fhTimeVsBC(0),
91  fhTimeSumSq(),
92  fhTimeEnt(),
93  fhTimeSum(),
94  fhTimeLGSumSq(),
95  fhTimeLGEnt(),
96  fhTimeLGSum(),
97  fhAllAverageBC(),
98  fhAllAverageLGBC(),
99  fhTimeDsup(),
100  fhTimeDsupBC(),
101  fhRawTimeVsIdBC(),
102  fhRawTimeSumBC(),
103  fhRawTimeEntriesBC(),
104  fhRawTimeSumSqBC(),
105  fhRawTimeVsIdLGBC(),
106  fhRawTimeSumLGBC(),
107  fhRawTimeEntriesLGBC(),
108  fhRawTimeSumSqLGBC()
109 
110 {
111  for(Int_t i = 0; i < kNBCmask; i++)
112  {
113  fhAllAverageBC[i]=0;
114  fhAllAverageLGBC[i]=0;
115 
116  fhTimeSumSq[i]=0;
117  fhTimeEnt[i]=0;
118  fhTimeSum[i]=0;
119 
120  fhTimeLGSumSq[i]=0;
121  fhTimeLGEnt[i]=0;
122  fhTimeLGSum[i]=0;
123 
124  fhRawTimeVsIdBC[i]=0;
125  fhRawTimeSumBC[i]=0;
126  fhRawTimeEntriesBC[i]=0;
127  fhRawTimeSumSqBC[i]=0;
128 
129  fhRawTimeVsIdLGBC[i]=0;
130  fhRawTimeSumLGBC[i]=0;
132  fhRawTimeSumSqLGBC[i]=0;
133  }
134 
135  //set default cuts for calibration and geometry name
136  SetDefaultCuts();
137 
138  //T0 TOF time
140 
141  // Define input and output slots here
142  // Input slot #0 works with a TChain
143  DefineInput(0, TChain::Class());
144 
145  // Output slot #0 id reserved by the base class for AOD
146  // Output slot #1 writes into a TH1 container
147  DefineOutput(1, TList::Class());
148 
149 } // End ctor
150 
151 //_____________________________________________________________________
157 //void AliAnalysisTaskEMCALTimeCalib::LocalInit()
158 //{
159 // AliDebug(1,"AliAnalysisTaskEMCALTimeCalib::LocalInit()");
160 //}
161 
163 //_____________________________________________________________________
165 {
166  if(fReferenceFileName.Length()!=0){
167  TFile *myFile = TFile::Open(fReferenceFileName.Data());
168  AliInfo(Form("Reference file: %s, pointer %p",fReferenceFileName.Data(),myFile));
169  if(myFile==0x0)
170  AliDebug(1,"*** NO REFERENCE FILE");
171  else
172  {
173  AliDebug(1,"*** OK TFILE");
174  // connect ref run here
175  for(Int_t i = 0; i < kNBCmask; i++)
176  {
177  fhAllAverageBC[i]=(TH1D*) myFile->Get(Form("hAllTimeAvBC%d",i));
178  fhAllAverageLGBC[i]=(TH1D*) myFile->Get(Form("hAllTimeAvLGBC%d",i));
179  }
180 
181  AliDebug(1,Form("hAllAverage entries %d", (Int_t)fhAllAverageBC[0]->GetEntries() ));
182  AliDebug(1,Form("hAllAverage entries2 %d",(Int_t)fhAllAverageBC[2]->GetEntries() ));
183  AliDebug(1,Form("hAllAverageLG entries %d", (Int_t)fhAllAverageLGBC[0]->GetEntries() ));
184  AliDebug(1,Form("hAllAverageLG entries2 %d",(Int_t)fhAllAverageLGBC[2]->GetEntries() ));
185 
186  }
187  }//end of reference file is provided
188 } // End of AliAnalysisTaskEMCALTimeCalib::LoadReferenceHistos()
189 
190 //_____________________________________________________________________
194 {
195  AliDebug(1,"AnalysisTaskEMCalTimeCalib::Notify()");
196  AliDebug(2,Form("Notify(): EMCal geometry: fgeom = %p, fGeometryName=%s\n ",fgeom,fGeometryName.Data()));
197 
198  if (!InputEvent())
199  {
200  AliFatal("ERROR: InputEvent not set");
201  return kFALSE;
202  }
203  else AliDebug(1,"Good, InputEvent set");
204 
205  fRunNumber = InputEvent()->GetRunNumber();
206  AliDebug(1,Form("RunNumber %d", fRunNumber));
207 
208  // Init EMCAL geometry
209  if (!fgeom) SetEMCalGeometry();
210  //Init EMCAL geometry done
211 
212  return kTRUE;
213 }
214 
215 //_____________________________________________________________________
218 {
219  AliDebug(1,"AliAnalysisTaskEMCALTimeCalib::SetEMCalGeometry()");
220  if(fGeometryName.Length()==0){
221  fgeom=AliEMCALGeometry::GetInstanceFromRunNumber(fRunNumber);
222  AliInfo(Form("Get EMCAL geometry name <%s> for run %d",fgeom->GetName(),fRunNumber));
223  } else {
224  fgeom = AliEMCALGeometry::GetInstance(fGeometryName.Data());
225  AliInfo(Form("Set EMCAL geometry name to <%s>",fGeometryName.Data()));
226  }
227 
228  if (!fgeom){
229  AliWarning("Make sure the EMCal geometry is set properly !");
230  } else {
231  AliDebug(1,Form("EMCal geometry properly set: fGeom = %p, fGeometryName=%s",fgeom,fGeometryName.Data()));
232  }
233 
234  return kTRUE;
235 }
236 
237 //_____________________________________________________________________
240 {
241  //method under development
242  AliInfo(Form("<D> -- Run # = %d", fRunNumber));
243  AliInfo("prepare TOFT0maker!!");
244  //cout<<"Run "<< fRunNumber<<" in TOFT0maker"<<endl;
245 
246 
247  AliCDBManager * cdb = AliCDBManager::Instance();
248  cdb->SetDefaultStorage("raw://");
249  cdb->SetRun(fRunNumber);
250 
251  AliESDpid *extPID=new AliESDpid();
252 
253  // Wonder if some have to be declared as private variables??
254  // AliESDpid *extPID = new AliESDpid();
255  // AliTOFcalib * tofCalib = new AliTOFcalib();
256  // tofCalib->SetCalibrateTOFsignal(kTRUE);
257  // tofCalib->Init();
258 
259  fTOFmaker = new AliTOFT0maker(extPID);
260  fTOFmaker->SetTimeResolution(115.0); // if you want set the TOF res
261  // fTOFmaker = new AliTOFT0maker(extPID,tofCalib);
262  // fTOFmaker->SetTimeResolution(130.0);
263 
264  //cout<<"extPID "<<extPID<<" fTOFmaker "<<fTOFmaker<<endl;
265 
266 }// End PrepareTOFT0maker
267 
268 //________________________________________________________________________
272 {
273  AliDebug(1,"AliAnalysisTaskEMCALTimeCalib::UserCreateOutputObjects()");
274 
275  Double_t fineTmin = -500;
276  Double_t fineTmax = 400;
277  Double_t fineInterval = 0.20;
278  Int_t nfinebin = (fineTmax-fineTmin)/fineInterval;
279  //cout << "<D> nfinebin = " << nfinebin << endl;
280 
281  Int_t nChannels = 17664;
282  //book histograms
283  fhcalcEvtTime = new TH1F("fhcalcEvtTime","calculated event time from T0",nfinebin, fineTmin,fineTmax);
284  fhcalcEvtTime->GetXaxis()->SetTitle("T ");
285  fhcalcEvtTime->GetYaxis()->SetTitle("Counts (a.u.)");
286 
287  fhEvtTimeHeader = new TH1F("fhEvtTimeHeader","event time from header",nfinebin, fineTmin,fineTmax);
288  fhEvtTimeHeader->GetXaxis()->SetTitle("T ");
289  fhEvtTimeHeader->GetYaxis()->SetTitle("Counts (a.u.)");
290 
291  fhEvtTimeDiff = new TH1F("fhEvtTimeDiff","event time difference",nfinebin, fineTmin,fineTmax);
292  fhEvtTimeDiff->GetXaxis()->SetTitle("#Delta T ");
293  fhEvtTimeDiff->GetYaxis()->SetTitle("Counts (a.u.)");
294 
295  fhEventType = new TH1F("fhEventType","event type",10, 0.,10.);
296  fhEventType ->GetXaxis()->SetTitle("Type ");
297  fhEventType ->GetYaxis()->SetTitle("Counts (a.u.)");
298  fhTcellvsTOFT0 = new TH2F("hTcellvsTOFT0", " T_cell vs TOFT0", 500,-600.0,+400.0,1200,300.0,900.0);
299  fhTcellvsTOFT0HD = new TH2F("hTcellvsTOFT0HD", " T_cell vs TOFT0,HighEnergy", 500,-600.0,+400.0,4000,500.0,700.0);
300  fhTcellvsSM = new TH2F("hTcellvsSM", " T_cell vs SM", 20,0,20,300,300,900);
301  fhEneVsAbsIdHG = new TH2F("fhEneVsAbsIdHG", "energy vs ID for HG",1000,0,18000,200,0,10);
302  fhEneVsAbsIdLG = new TH2F("fhEneVsAbsIdLG", "energy vs ID for LG",1000,0,18000,100,0,20);
303 
304  for (Int_t i = 0; i < kNBCmask ; i++)
305  {
306  //already after correction
307  //high gain
308  fhTimeSumSq[i] = new TH1F(Form("hTimeSumSq%d", i),
309  Form("cell Sum Square time HG, BC %d ", i),
310  nChannels,0.,(Double_t)nChannels);
311  fhTimeSumSq[i]->SetYTitle("Sum Sq Time ");
312  fhTimeSumSq[i]->SetXTitle("AbsId");
313 
314  fhTimeSum[i] = new TH1F(Form("hTimeSum%d", i),
315  Form("cell Sum time HG, BC %d ", i),
316  nChannels,0.,(Double_t)nChannels);
317  fhTimeSum[i]->SetYTitle("Sum Time ");
318  fhTimeSum[i]->SetXTitle("AbsId");
319 
320  fhTimeEnt[i] = new TH1F(Form("hTimeEnt%d", i),
321  Form("cell Entries HG, BC %d ", i),
322  nChannels,0.,(Double_t)nChannels);
323  fhTimeEnt[i]->SetYTitle("Entries for Time ");
324  fhTimeEnt[i]->SetXTitle("AbsId");
325 
326  //low gain
327  fhTimeLGSumSq[i] = new TH1F(Form("hTimeLGSumSq%d", i),
328  Form("cell Sum Square time LG, BC %d ", i),
329  nChannels,0.,(Double_t)nChannels);
330  fhTimeLGSumSq[i]->SetYTitle("Sum Sq Time ");
331  fhTimeLGSumSq[i]->SetXTitle("AbsId");
332 
333  fhTimeLGSum[i] = new TH1F(Form("hTimeLGSum%d", i),
334  Form("cell Sum time LG, BC %d ", i),
335  nChannels,0.,(Double_t)nChannels);
336  fhTimeLGSum[i]->SetYTitle("Sum Time ");
337  fhTimeLGSum[i]->SetXTitle("AbsId");
338 
339  fhTimeLGEnt[i] = new TH1F(Form("hTimeLGEnt%d", i),
340  Form("cell Entries LG, BC %d ", i),
341  nChannels,0.,(Double_t)nChannels);
342  fhTimeLGEnt[i]->SetYTitle("Entries for Time ");
343  fhTimeLGEnt[i]->SetXTitle("AbsId");
344 
345  //raw time histograms
346  //high gain
347  fhRawTimeVsIdBC[i] = new TH2F(Form("RawTimeVsIdBC%d", i),
348  Form("cell raw time vs ID for high gain BC %d ", i),
349  nChannels,0.,(Double_t)nChannels,600,300,900);
350  fhRawTimeVsIdBC[i]->SetXTitle("AbsId");
351  fhRawTimeVsIdBC[i]->SetYTitle("Time");
352 
353  fhRawTimeSumBC[i] = new TH1F(Form("RawTimeSumBC%d", i),
354  Form("sum of cell raw time for high gain BC %d ", i),
355  nChannels,0.,(Double_t)nChannels);
356  fhRawTimeSumBC[i]->SetXTitle("AbsId");
357  fhRawTimeSumBC[i]->SetYTitle("Sum Time");
358 
359  fhRawTimeEntriesBC[i] = new TH1F(Form("RawTimeEntriesBC%d", i),
360  Form("No. entries of cells raw time for high gain BC %d ", i),
361  nChannels,0.,(Double_t)nChannels);
362  fhRawTimeEntriesBC[i]->SetXTitle("AbsId");
363  fhRawTimeEntriesBC[i]->SetYTitle("Entries for Time ");
364 
365  fhRawTimeSumSqBC[i] = new TH1F(Form("RawTimeSumSqBC%d", i),
366  Form("sum of (cell raw time)^2 for high gain BC %d ", i),
367  nChannels,0.,(Double_t)nChannels);
368  fhRawTimeSumSqBC[i]->SetXTitle("AbsId");
369  fhRawTimeSumSqBC[i]->SetYTitle("Sum Sq Time");
370 
371  //low gain
372  fhRawTimeVsIdLGBC[i] = new TH2F(Form("RawTimeVsIdLGBC%d", i),
373  Form("cell raw time vs ID for low gain BC %d ", i),
374  nChannels,0.,(Double_t)nChannels,600,300,900);
375  fhRawTimeVsIdLGBC[i]->SetXTitle("AbsId");
376  fhRawTimeVsIdLGBC[i]->SetYTitle("Time");
377 
378  fhRawTimeSumLGBC[i] = new TH1F(Form("RawTimeSumLGBC%d", i),
379  Form("sum of cell raw time for low gain BC %d ", i),
380  nChannels,0.,(Double_t)nChannels);
381  fhRawTimeSumLGBC[i]->SetXTitle("AbsId");
382  fhRawTimeSumLGBC[i]->SetYTitle("Sum Time");
383 
384  fhRawTimeEntriesLGBC[i] = new TH1F(Form("RawTimeEntriesLGBC%d", i),
385  Form("No. entries of cells raw time for low gain BC %d ", i),
386  nChannels,0.,(Double_t)nChannels);
387  fhRawTimeEntriesLGBC[i]->SetXTitle("AbsId");
388  fhRawTimeEntriesLGBC[i]->SetYTitle("Entries for Time ");
389 
390  fhRawTimeSumSqLGBC[i] = new TH1F(Form("RawTimeSumSqLGBC%d", i),
391  Form("sum of (cell raw time)^2 for low gain BC %d ", i),
392  nChannels,0.,(Double_t)nChannels);
393  fhRawTimeSumSqLGBC[i]->SetXTitle("AbsId");
394  fhRawTimeSumSqLGBC[i]->SetYTitle("Sum Sq Time");
395 
396 
397  for (Int_t j = 0; j < kNSM ; j++)
398  {
399  fhTimeDsupBC[j][i]= new TH2F(Form("SupMod%dBC%d",j,i), Form("SupMod %d time_vs_E BC %d",j,i),500,0.0,20.0,1400,-350.0,350.0);
400  fhTimeDsupBC[j][i]->SetYTitle(" Time (ns) ");
401  fhTimeDsupBC[j][i]->SetXTitle(" E (GeV) ");
402  }
403  }
404 
405  for (Int_t jj = 0; jj < kNSM ; jj++)
406  {
407  fhTimeDsup[jj] = new TH2F(Form("SupMod%d",jj), Form("SupMod %d time_vs_E ",jj),500,0.0,20.0,1400,-350.0,350.0);
408  fhTimeDsup[jj]->SetYTitle(" Time (ns) ");
409  fhTimeDsup[jj]->SetXTitle(" E (GeV) ");
410  }
411 
412  fhTimeVsBC = new TH2F("TimeVsBC"," SupMod time_vs_BC ", 4001,-0.5,4000.5,400,200.0,1000.0);
413 
414 
415  //add histos to list
416  fOutputList = new TList();
417 
421  fOutputList->Add(fhEventType);
424  fOutputList->Add(fhTcellvsSM);
427 
428  for (Int_t i = 0; i < kNBCmask ; i++)
429  {
430  fOutputList->Add(fhTimeSumSq[i]);
431  fOutputList->Add(fhTimeEnt[i]);
432  fOutputList->Add(fhTimeSum[i]);
433 
434  fOutputList->Add(fhTimeLGSumSq[i]);
435  fOutputList->Add(fhTimeLGEnt[i]);
436  fOutputList->Add(fhTimeLGSum[i]);
437 
438  fOutputList->Add(fhRawTimeVsIdBC[i]);
439  fOutputList->Add(fhRawTimeSumBC[i]);
441  fOutputList->Add(fhRawTimeSumSqBC[i]);
442 
444  fOutputList->Add(fhRawTimeSumLGBC[i]);
447 
448  for (Int_t j = 0; j < kNSM ; j++){
449  fOutputList->Add(fhTimeDsupBC[j][i]);
450  }
451  }
452 
453  for (Int_t j = 0; j < kNSM ; j++)
454  {
455  fOutputList->Add(fhTimeDsup[j]);
456  }
457 
458  fOutputList->Add(fhTimeVsBC);
459 
460  fOutputList->SetOwner(kTRUE);
461  PostData(1,fOutputList);
462 
463 
464 } // End of AliAnalysisTaskEMCALTimeCalib::UserCreateOuputObjects()
465 
466 //________________________________________________________________________
469 {
470  // Called for each event
471  AliDebug(2,Form("UserExec: EMCal geometry: fgeom = %p fGeometryName %s",fgeom,fGeometryName.Data()));
472  AliVEvent *event = InputEvent();
473  //cout<<"T0TOF "<<event->GetT0TOF()<<endl;//bad idea
474  //cout<< fEvent->GetTOFHeader()->GetDefaultEventTimeVal()<<endl;
475  AliDebug(2,Form("TOF time from header %f ps",event->GetTOFHeader()->GetDefaultEventTimeVal()));
476  fhEvtTimeHeader->Fill(event->GetTOFHeader()->GetDefaultEventTimeVal());
477 
478  //fEvent = dynamic_cast<AliESDEvent*>(event);
479  if (!event) {
480  AliError("ESD not available, exit");
481  fhEventType->Fill(0.5);
482  return;
483  }
484 
485  if(fPileupFromSPD==kTRUE){
486  if(event->IsPileupFromSPD(3,0.8,3.,2.,5.)){
487  AliDebug(1,"Event: PileUp skip.");
488  fhEventType->Fill(1.5);
489  return;
490  }
491  }
492 
493  TString triggerclasses = event->GetFiredTriggerClasses();
494  if(triggerclasses=="") {
495  fhEventType->Fill(2.5);
496  return;
497  }
498 
499  Int_t eventType = ((AliVHeader*)event->GetHeader())->GetEventType();
500  // physics events eventType=7, select only those
501  AliDebug(1,Form("Triggerclasses %s, eventType %d",triggerclasses.Data(),eventType));
502  if(eventType != 7) {
503  fhEventType->Fill(3.5);
504  return;
505  }
506 
507  // Check trigger
508  Bool_t bMB = kFALSE;
509  Bool_t bL0 = kFALSE;
510  Bool_t bL1G = kFALSE;
511  Bool_t bL1J = kFALSE;
512 
513  if(triggerclasses.Contains("CINT7-B-NOPF-ALLNOTRD") ||
514  triggerclasses.Contains("CINT7-I-NOPF-ALLNOTRD") ||
515  triggerclasses.Contains("CINT1-I-NOPF-ALLNOTRD") ||
516  triggerclasses.Contains("CINT1-B-NOPF-ALLNOTRD") ||
517  triggerclasses.Contains("CINT8") ||
518  triggerclasses.Contains("CINT7") ||
519  triggerclasses.Contains("CPBI2_B1-B-NOPF-ALLNOTRD") ) bMB = kTRUE;
520 
521  if(triggerclasses.Contains("CEMC7-B-NOPF-CENTNOTRD") ||
522  triggerclasses.Contains("CEMC1-B-NOPF-CENTNOTRD") ||
523  triggerclasses.Contains("CEMC7") ||
524  triggerclasses.Contains("CEMC8") ||
525  triggerclasses.Contains("CEMC8-B-NOPF-CENTNOTRD") ) bL0 = kTRUE;
526 
527  if(triggerclasses.Contains("CEMC7EG1-B-NOPF-CENTNOTRD") ||
528  triggerclasses.Contains("CEMC7EG2-B-NOPF-CENTNOTRD") ||
529  triggerclasses.Contains("CEMC8EG1-B-NOPF-CENTNOTRD") ||
530  triggerclasses.Contains("CEMC8EGA") ||
531  triggerclasses.Contains("CEMC7EGA") ||
532  triggerclasses.Contains("CPBI2EGA") ) bL1G = kTRUE;
533 
534 
535  if(triggerclasses.Contains("CEMC7EJ1-B-NOPF-CENTNOTRD") ||
536  triggerclasses.Contains("CEMC7EJ2-B-NOPF-CENTNOTRD") ||
537  triggerclasses.Contains("CEMC8EJ1-B-NOPF-CENTNOTRD") ||
538  triggerclasses.Contains("CEMC7EJE") ||
539  triggerclasses.Contains("CEMC8EJE") ||
540  triggerclasses.Contains("CPBI2EJE") ) bL1J = kTRUE;
541 
542  if( bL1G || bL1J || bL0 ){ fhEventType->Fill(4.5);}
543  if( bMB ){ fhEventType->Fill(5.5);}
544 
545 
546  // if(bL1G || bL1J || bL0){
547 
548  // Prepare TOFT0 maker at the beginning of a run
549 // if (event->GetRunNumber() != fRunNumber)
550 // {
551 // fRunNumber = event->GetRunNumber();
552 // // PrepareTOFT0maker();
553 // // cout<<"tofT0maker per run"<<fRunNumber<<endl;
554 // }// fi Check if run number has changed
555 
556  // --- Use of AliTOFT0maker
557  Double_t calcolot0=0.0;
558  if(!AODEvent()){
559  Double_t* timeTOFtable;
560  timeTOFtable=fTOFmaker->ComputeT0TOF(dynamic_cast<AliESDEvent*>(event));
561  AliDebug(2,Form("TOF time %f ps, resolution %f ps, tracks at TOF %f/used %f",timeTOFtable[0],timeTOFtable[1],timeTOFtable[3],timeTOFtable[7]));
562  //cout<<"event time "<<timeTOFtable[0]<<" resolution "<<timeTOFtable[1]<<"ps av. ev. time "<<timeTOFtable[2]<<" trks at TOF "<<timeTOFtable[3]<<" calc evnt time "<<timeTOFtable[4]<<" resolution "<<timeTOFtable[5]<<" tracks used "<<timeTOFtable[7]<<endl;
563  calcolot0=timeTOFtable[0];
564  }
565 
566  if (!fhcalcEvtTime) {
567  AliWarning("<E> fhcalcEvtTime not available");
568  return;
569  }// fi no simple histo present
570 
571  fhcalcEvtTime->Fill(calcolot0);
572  if(calcolot0 != 0 && event->GetTOFHeader()->GetDefaultEventTimeVal() != 0 )
573  fhEvtTimeDiff->Fill(calcolot0-event->GetTOFHeader()->GetDefaultEventTimeVal());
574 
575  TRefArray* caloClusters = new TRefArray();
576  event->GetEMCALClusters(caloClusters);
577  // cout << " ###########Bunch Cross nb = " << event->GetBunchCrossNumber() << endl;
578 
579  Int_t BunchCrossNumber =event->GetBunchCrossNumber();
580 
581  Float_t offset=0;
582 
583  Int_t nBC = 0;
584  nBC = BunchCrossNumber%4;
585  //Int_t nTriggerMask =event->GetTriggerMask();
586  // cout << " nBC " << nBC << " nTriggerMask " << nTriggerMask<< endl;
587  Float_t timeBCoffset = 0.; //manual offest
588  // if( nBC%4 ==0 || nBC%4==1) timeBCoffset = 100.; // correction was for LHC11 when BC was not corrected
589 
590  Int_t nclus = caloClusters->GetEntries();
591  AliDebug(1,Form("###########Bunch Cross nb = %d nclus = %d",nBC,nclus ));
592  //cout << " ###########Bunch Cross nb = " << nBC <<" nclus= "<< nclus<< endl;
593  //Int_t ntracks = event-> GetNumberOfTracks() ;
594 
595  AliVCaloCells &cells= *(event->GetEMCALCells());//it is cluster independent
596  //Variables used plenty in loops
597  Int_t nSupMod=-1, nModule=-1;
598  Int_t iphi=-1, ieta=-1, nIphi=-1, nIeta=-1;
599  Int_t absId=-1;
600  Float_t hkdtime=0.0;
601  Float_t amp=0.0;
602  Bool_t isHighGain=kTRUE;
603 
604  for (Int_t icl = 0; icl < nclus; icl++) {
605  //ESD and AOD CaloCells carries the same information
606  AliVCluster* clus = (AliVCluster*)caloClusters->At(icl);
607  if(!AcceptCluster(clus)) continue;
608 
609  //cout<<"nCells="<< clus->GetNCells();<<endl;
610 
611  UShort_t * index = clus->GetCellsAbsId() ;
612 
613  for(Int_t i = 0; i < clus->GetNCells() ; i++) {
614  absId = index[i]; // or clus->GetCellNumber(i) ;
615  hkdtime = cells.GetCellTime(absId) * 1.0e09; // to get ns
616  amp = cells.GetCellAmplitude(absId) ;
617  isHighGain = cells.GetCellHighGain(absId);
618  //cout<<"cell absID: "<<absId<<" cellTime: "<<hkdtime<<" cellaplit: "<< amp<<endl;
619 
620  //main histograms with raw time information
621  if(amp>fMinCellEnergy){
622  if(isHighGain){
623  fhRawTimeVsIdBC[nBC]->Fill(absId,hkdtime);
624  fhRawTimeSumBC[nBC]->Fill(absId,hkdtime);
625  fhRawTimeEntriesBC[nBC]->Fill(absId,1.);
626  fhRawTimeSumSqBC[nBC]->Fill(absId,hkdtime*hkdtime);
627  }else{
628  fhRawTimeVsIdLGBC[nBC]->Fill(absId,hkdtime);
629  fhRawTimeSumLGBC[nBC]->Fill(absId,hkdtime);
630  fhRawTimeEntriesLGBC[nBC]->Fill(absId,1.);
631  fhRawTimeSumSqLGBC[nBC]->Fill(absId,hkdtime*hkdtime);
632  }
633  }
634  //fgeom->PrintCellIndexes(absId);
635  //fgeom->PrintCellIndexes(absId,1);
636 
637  // GEOMETRY tranformations
638  fgeom->GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
639  fgeom->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
640 
641  // other histograms for cross-check
642  CheckCellRCU(nSupMod,ieta,iphi);//SM, column, row
643 
644  fhTcellvsSM->Fill(nSupMod,hkdtime);
645  if(isHighGain==kTRUE) {fhEneVsAbsIdHG->Fill(absId,amp);}
646  else {fhEneVsAbsIdLG->Fill(absId,amp);}
647 
648  fhTimeVsBC->Fill(1.*BunchCrossNumber,hkdtime-timeBCoffset);
649  if(isHighGain==kTRUE){
650  if(fhAllAverageBC[nBC]!=0) {//comming from file after the first iteration
651  offset = fhAllAverageBC[nBC]->GetBinContent(absId);
652  }
653  } else {
654  if(fhAllAverageLGBC[nBC]!=0) {//comming from file after the first iteration
655  offset = fhAllAverageLGBC[nBC]->GetBinContent(absId);
656  }
657  }
658  //if(offset==0)cout<<"offset 0 in SM "<<nSupMod<<endl;
659 
660 
661  if(amp>0.5) {
662  fhTimeDsup[nSupMod]->Fill(amp,hkdtime-offset);
663  fhTimeDsupBC[nSupMod][nBC]->Fill(amp,hkdtime-offset);
664  }
665 
666  if(amp>0.9) {
667  fhTcellvsTOFT0HD->Fill(calcolot0, hkdtime);
668  }
669 
670  fhTcellvsTOFT0->Fill(calcolot0, hkdtime-offset);
671 
672  hkdtime = hkdtime-timeBCoffset;//time corrected by manual offset (default=0)
673  Float_t hkdtimecorr;
674  hkdtimecorr= hkdtime-offset;//time after first iteration
675 
676  //main histograms after the first itereation for calibration constants
677  //if(hkdtimecorr>=-20. && hkdtimecorr<=20. && amp>0.9 ) {
678  if(hkdtimecorr>=fMinTime && hkdtimecorr<=fMaxTime && amp>fMinCellEnergy ) {
679  // per cell
680 // Float_t entriesTime=fhTimeEnt[nBC]->GetBinContent(absId)+1;
681 // Float_t sumTimeSq=(fhTimeSumSq[nBC]->GetBinContent(absId)+(hkdtime*hkdtime));
682 // Float_t sumTime=(fhTimeSum[nBC]->GetBinContent(absId)+hkdtime);
683 //
684 // fhTimeEnt[nBC]->SetBinContent(absId,entriesTime);
685 // fhTimeSumSq[nBC]->SetBinContent(absId,sumTimeSq);
686 // fhTimeSum[nBC]->SetBinContent(absId,sumTime);
687 
688  if(isHighGain){
689  fhTimeEnt[nBC]->Fill(absId,1.);
690  fhTimeSumSq[nBC]->Fill(absId,hkdtime*hkdtime);
691  fhTimeSum[nBC]->Fill(absId,hkdtime);
692  }else{
693  fhTimeLGEnt[nBC]->Fill(absId,1.);
694  fhTimeLGSumSq[nBC]->Fill(absId,hkdtime*hkdtime);
695  fhTimeLGSum[nBC]->Fill(absId,hkdtime);
696  }
697 
698 
699  } // hkdtime:[-20;20]
700  } // end icell
701  } //end cluster
702 
703 
704  // Post output data.
705  //cout<<"Post data and delete caloClusters"<<endl;
706  caloClusters->Delete();
707  delete caloClusters;
708 // } // end if trigger type
709 
710  PostData(1, fOutputList);
711 } // End of AliAnalysisTaskEMCALTimeCalib::UserExec()
712 
713 //________________________________________________________________________
717 {
718  fOutputList = dynamic_cast<TList*> (GetOutputData(1));
719 
720  if(fTOFmaker) delete fTOFmaker;
721 
722 
723  if (!fOutputList)
724  {
725  AliDebug(1,"ERROR: Output list not available");
726  return;
727  }
728 } // End of AliAnalysisTaskEMCALTimeCalib::Terminate
729 
730 //________________________________________________________________________
733 {
734  //fix with noisy EMCAL fee card
735  Int_t nCells = clus->GetNCells();
736 
737  if(clus->IsEMCAL())
738  {
739  if ((clus->E() > fMaxClusterEnergy && nCells > fMaxNcells ) || nCells > fMaxNcells){
740  AliDebug(1,"very big cluster with enormous energy - cluster rejected");
741  return kFALSE;
742  }
743  }
744 
745  // remove other than photonlike
746  Double_t lambda0=clus->GetM02();
747  if (lambda0>fMaxLambda0LG || lambda0<fMinLambda0LG){
748  AliDebug(1,"lambda0 loose cut failed - cluster rejected");
749  return kFALSE;
750  }
751 
752  // remove matched clusters
753  Double_t Dx=clus->GetTrackDx();
754  Double_t Dz=clus->GetTrackDz();
755  Double_t Rtrack = TMath::Sqrt(Dx*Dx+Dz*Dz);
756  if (Rtrack <fMaxRtrack)
757  {
758  AliDebug(1,"track matched - cluster rejected");
759  return kFALSE;
760  }
761 
762  if (nCells<fMinNcells)
763  {
764  AliDebug(1,"single cell cluster - cluster rejected");
765  return kFALSE;
766  }
767 
768  if(clus->E()<fMinClusterEnergy)
769  {
770  AliDebug(1,"cluster energy < 1 GeV- cluster rejected");
771  return kFALSE;
772  }
773 
774 
775  if(!IsLowGainCellInCluster(clus)) {//no low gain cell in cluster
776  //apply more strict lambda0^2 cut
777  if (lambda0>fMaxLambda0 || lambda0<fMinLambda0){
778  AliDebug(1,"lambda0 strict cut failed - cluster rejected");
779  return kFALSE;
780  }
781  }
782 
783 
784 
785 
786  return kTRUE;
787 }//End AliAnalysisTaskEMCALTimeCalib::AcceptCluster
788 
789 //________________________________________________________________________
792  UShort_t * index = clus->GetCellsAbsId() ;
793  AliVCaloCells &cells= *(InputEvent()->GetEMCALCells());
794  for(Int_t i = 0; i < clus->GetNCells() ; i++) {
795  if(cells.GetCellHighGain(index[i])==kFALSE) return kTRUE;//low gain cell
796  }
797  return kFALSE;
798 
799 }
800 
801 //________________________________________________________________________
803 Bool_t AliAnalysisTaskEMCALTimeCalib::CheckCellRCU(Int_t nSupMod,Int_t icol,Int_t irow)
804 {
805  Int_t iRCU;
806  if(nSupMod < 10 || (nSupMod >= 12 && nSupMod <18) )
807  {
808  if (0<=irow&&irow<8) iRCU=0; // first cable row
809  else if (8<=irow&&irow<16 && 0<=icol&&icol<24) iRCU=0; // first half;
810  //second cable row
811  //RCU1
812  else if (8<=irow&&irow<16 && 24<=icol&&icol<48) iRCU=1; // second half;
813  //second cable row
814  else if (16<=irow&&irow<24) iRCU=1; // third cable row
815 
816  if (nSupMod%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
817  }
818  else
819  {
820  // Last 2 SM have one single SRU, just assign RCU 0
821  iRCU = 0 ;
822  }
823 
824  //cout<<"RCU:"<<iRCU<<endl;
825  if (iRCU<0)
826  AliFatal(Form("Wrong EMCAL/DCAL RCU number = %d\n", iRCU));
827 
828  return kTRUE;
829 }//End AliAnalysisTaskEMCALTimeCalib::CheckCellRCU
830 
831 //________________________________________________________________________
834 {
835  fMinClusterEnergy=1.0;//0.5//0.7
836  fMaxClusterEnergy=500;
837  fMinNcells=2;
838  fMaxNcells=200;
839  fMinLambda0=0.1;
840  fMaxLambda0=0.4;
841  fMinLambda0LG=0.1;
842  fMaxLambda0LG=4.0;
843  fMaxRtrack=0.025;
844  fMinCellEnergy=0.4;//0.1//0.4
845  fReferenceFileName="";//Reference.root
846  fGeometryName="";//EMCAL_COMPLETE12SMV1_DCAL_8SM
847  fPileupFromSPD=kFALSE;
848  fMinTime=-20.;
849  fMaxTime=20.;
850 }
851 
852 //________________________________________________________________________
857 void AliAnalysisTaskEMCALTimeCalib::ProduceCalibConsts(TString inputFile,TString outputFile,Bool_t isFinal)
858 {
859  TFile *file =new TFile(inputFile.Data());
860  if(file==0x0) {
861  //AliWarning("Input file does not exist!");
862  return;
863  }
864 
865  TList *list=(TList*)file->Get("chistolist");
866  if(list==0x0)
867  {
868  //AliWarning("List chistolist does not exist in file!");
869  return;
870  }
871 
872  //high gain
873  TH1F *h1[4];
874  TH1F *h2[4];
875  TH1F *h3[4];
876  TH1F *hAllTimeAvBC[4];
877  TH1F *hAllTimeRMSBC[4];
878 
879  //low gain
880  TH1F *h4[4];
881  TH1F *h5[4];
882  TH1F *h6[4];
883  TH1F *hAllTimeAvLGBC[4];
884  TH1F *hAllTimeRMSLGBC[4];
885 
886  if(isFinal==kFALSE){//first itereation
887  for(Int_t i=0;i<4;i++){
888  h1[i]=(TH1F *)list->FindObject(Form("RawTimeSumBC%d",i));
889  h2[i]=(TH1F *)list->FindObject(Form("RawTimeEntriesBC%d",i));
890  h3[i]=(TH1F *)list->FindObject(Form("RawTimeSumSqBC%d",i));
891 
892  h4[i]=(TH1F *)list->FindObject(Form("RawTimeSumLGBC%d",i));
893  h5[i]=(TH1F *)list->FindObject(Form("RawTimeEntriesLGBC%d",i));
894  h6[i]=(TH1F *)list->FindObject(Form("RawTimeSumSqLGBC%d",i));
895  }
896  } else {//final iteration
897  for(Int_t i=0;i<4;i++){
898  h1[i]=(TH1F *)list->FindObject(Form("hTimeSum%d",i));
899  h2[i]=(TH1F *)list->FindObject(Form("hTimeEnt%d",i));
900  h3[i]=(TH1F *)list->FindObject(Form("hTimeSumSq%d",i));
901 
902  h4[i]=(TH1F *)list->FindObject(Form("hTimeLGSum%d",i));
903  h5[i]=(TH1F *)list->FindObject(Form("hTimeLGEnt%d",i));
904  h6[i]=(TH1F *)list->FindObject(Form("hTimeLGSumSq%d",i));
905  }
906  }
907  //AliWarning("Input histograms read.");
908 
909  for(Int_t i=0;i<4;i++){
910  hAllTimeAvBC[i]=new TH1F(Form("hAllTimeAvBC%d",i),Form("hAllTimeAvBC%d",i),h1[i]->GetNbinsX(),h1[i]->GetXaxis()->GetXmin(),h1[i]->GetXaxis()->GetXmax());
911  hAllTimeRMSBC[i]=new TH1F(Form("hAllTimeRMSBC%d",i),Form("hAllTimeRMSBC%d",i),h3[i]->GetNbinsX(),h3[i]->GetXaxis()->GetXmin(),h3[i]->GetXaxis()->GetXmax());
912 
913  hAllTimeAvLGBC[i]=new TH1F(Form("hAllTimeAvLGBC%d",i),Form("hAllTimeAvLGBC%d",i),h4[i]->GetNbinsX(),h4[i]->GetXaxis()->GetXmin(),h4[i]->GetXaxis()->GetXmax());
914  hAllTimeRMSLGBC[i]=new TH1F(Form("hAllTimeRMSLGBC%d",i),Form("hAllTimeRMSLGBC%d",i),h6[i]->GetNbinsX(),h6[i]->GetXaxis()->GetXmin(),h6[i]->GetXaxis()->GetXmax());
915  }
916 
917  //AliWarning("New histograms booked.");
918 
919  for(Int_t i=0;i<4;i++){
920  for(Int_t j=1;j<=h1[i]->GetNbinsX();j++){
921  //high gain
922  if(h2[i]->GetBinContent(j)!=0){
923  hAllTimeAvBC[i]->SetBinContent(j,h1[i]->GetBinContent(j)/h2[i]->GetBinContent(j));
924  hAllTimeRMSBC[i]->SetBinContent(j,TMath::Sqrt(h3[i]->GetBinContent(j)/h2[i]->GetBinContent(j)) );
925  } else {
926  hAllTimeAvBC[i]->SetBinContent(j,0.);
927  hAllTimeRMSBC[i]->SetBinContent(j,0.);
928  }
929  //low gain
930  if(h5[i]->GetBinContent(j)!=0){
931  hAllTimeAvLGBC[i]->SetBinContent(j,h4[i]->GetBinContent(j)/h5[i]->GetBinContent(j));
932  hAllTimeRMSLGBC[i]->SetBinContent(j,TMath::Sqrt(h6[i]->GetBinContent(j)/h5[i]->GetBinContent(j)) );
933  } else {
934  hAllTimeAvLGBC[i]->SetBinContent(j,0.);
935  hAllTimeRMSLGBC[i]->SetBinContent(j,0.);
936  }
937 
938  }
939  }
940 
941  //AliWarning("Average and rms calculated.");
942  TFile *fileNew=new TFile(outputFile.Data(),"recreate");
943  for(Int_t i=0;i<4;i++){
944  hAllTimeAvBC[i]->Write();
945  hAllTimeRMSBC[i]->Write();
946  hAllTimeAvLGBC[i]->Write();
947  hAllTimeRMSLGBC[i]->Write();
948  }
949 
950  //AliWarning(Form("Histograms saved in %s file.",outputFile.Data()));
951 
952  fileNew->Close();
953  delete fileNew;
954  file->Close();
955  delete file;
956 
957  //AliWarning("Pointers deleted. Memory cleaned.");
958 
959 }
Bool_t SetEMCalGeometry()
Set the EMCal Geometry.
TH1F * fhRawTimeEntriesLGBC[kNBCmask]
! 4 BCmask LG
TH2F * fhEneVsAbsIdHG
! energy of each cell for high gain cells with strange time
TH1F * fhRawTimeSumLGBC[kNBCmask]
! 4 BCmask LG
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Bool_t AcceptCluster(AliVCluster *clus)
Selection criteria of good cluster are set here.
TH2F * fhRawTimeVsIdBC[kNBCmask]
! 4 BCmask HG
TH1F * fhRawTimeSumSqBC[kNBCmask]
! 4 BCmask HG
Double_t fMinLambda0
minimum cluster lambda0
TH1F * fhcalcEvtTime
! spectrum calcolot0[0]
Double_t fMinTime
minimum cluster time after correction
AliTOFT0maker * fTOFmaker
pointer to get T0 from TOF
TList * list
void SetDefaultCuts()
Set default cuts for calibration.
TH1F * fhRawTimeSumSqLGBC[kNBCmask]
! 4 BCmask LG
Double_t fMaxRtrack
maximum cluster track distance
Double_t fMaxLambda0
maximum cluster lambda0
AliEMCALGeometry * fgeom
pointer to EMCal geometry
TH2F * fhTcellvsTOFT0HD
! time of cell vs TOF T0 time for higher energy threshold
Bool_t IsLowGainCellInCluster(AliVCluster *clus)
Check if low gain cell is in a cluster.
Double_t fMaxTime
maximum cluster time after correction
Double_t fMinClusterEnergy
minimum cluster energy
TString fReferenceFileName
! name of reference file
Int_t fMaxNcells
maximum number of cells in cluster
virtual void UserExec(Option_t *option)
Main loop executed for each event.
virtual void PrepareTOFT0maker()
Get T0 time from TOF.
Bool_t fPileupFromSPD
flag to set PileupFromSPD
TH1F * fhRawTimeSumBC[kNBCmask]
! 4 BCmask HG
Double_t fMinCellEnergy
minimum cell energy
Int_t fMinNcells
minimum number of cells in cluster
Task to work on Time Calibration for EMCal/DCal.
Double_t fMinLambda0LG
minimum cluster lambda0 Low Gain
Double_t fMaxClusterEnergy
maximum cluster energy
TH1F * fhEvtTimeHeader
! spectrum time from header
TH1F * fhEvtTimeDiff
! spectrum time difference
TH2F * fhTcellvsTOFT0
! time of cell vs TOF T0 time
Double_t fMaxLambda0LG
maximum cluster lambda0 Low Gain
TH1D * fhAllAverageLGBC[kNBCmask]
! 4 BCmask Low gain
static void ProduceCalibConsts(TString inputFile="time186319testWOL0.root", TString outputFile="Reference.root", Bool_t isFinal=kFALSE)
TH2F * fhEneVsAbsIdLG
! energy of each cell for low gain cells with strange time
TList * fOutputList
pointer to output list
TH1F * fhRawTimeEntriesBC[kNBCmask]
! 4 BCmask HG
TFile * file
Bool_t CheckCellRCU(Int_t nSupMod, Int_t icol, Int_t irow)
Check RCU for cell given by Super Module, column index, row index.
TH2F * fhRawTimeVsIdLGBC[kNBCmask]
! 4 BCmask LG
void LoadReferenceHistos()
Load reference Histograms from file.
TH1D * fhAllAverageBC[kNBCmask]
! 4 BCmask High gain
TH2F * fhTimeDsupBC[kNSM][kNBCmask]
! 20 x 4