AliPhysics  vAN-20150822 (d56cf94)
 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  AliFatal("*** NO REFERENCE FILE");
171  } else {
172  AliDebug(1,"*** OK TFILE");
173  // connect ref run here
174  for(Int_t i = 0; i < kNBCmask; i++)
175  {
176  fhAllAverageBC[i]=(TH1D*) myFile->Get(Form("hAllTimeAvBC%d",i));
177  if(fhAllAverageBC[i]==0x0) AliFatal(Form("Reference histogram for BC%d does not exist",i));
178  if(fhAllAverageBC[i]->GetEntries()==0)AliWarning(Form("fhAllAverageLGBC[%d]->GetEntries() = 0",i));
179  fhAllAverageLGBC[i]=(TH1D*) myFile->Get(Form("hAllTimeAvLGBC%d",i));
180  if(fhAllAverageLGBC[i]==0x0) AliFatal(Form("Reference LG histogram for BC%d does not exist",i));
181  if(fhAllAverageLGBC[i]->GetEntries()==0)AliFatal(Form("fhAllAverageLGBC[%d]->GetEntries() = 0",i));
182  }
183 
184  AliDebug(1,Form("hAllAverage entries BC0 %d", (Int_t)fhAllAverageBC[0]->GetEntries() ));
185  AliDebug(1,Form("hAllAverage entries BC2 %d",(Int_t)fhAllAverageBC[2]->GetEntries() ));
186  AliDebug(1,Form("hAllAverageLG entries BC0 %d", (Int_t)fhAllAverageLGBC[0]->GetEntries() ));
187  AliDebug(1,Form("hAllAverageLG entries BC2 %d",(Int_t)fhAllAverageLGBC[2]->GetEntries() ));
188 
189  }
190  } else {
191 //end of reference file is provided
192  AliFatal("You require to load reference histos from file but FILENAME is not provided");
193  }
194 } // End of AliAnalysisTaskEMCALTimeCalib::LoadReferenceHistos()
195 
196 //_____________________________________________________________________
200 {
201  AliDebug(1,"AnalysisTaskEMCalTimeCalib::Notify()");
202  AliDebug(2,Form("Notify(): EMCal geometry: fgeom = %p, fGeometryName=%s\n ",fgeom,fGeometryName.Data()));
203 
204  if (!InputEvent())
205  {
206  AliFatal("ERROR: InputEvent not set");
207  return kFALSE;
208  }
209  else AliDebug(1,"Good, InputEvent set");
210 
211  fRunNumber = InputEvent()->GetRunNumber();
212  AliDebug(1,Form("RunNumber %d", fRunNumber));
213 
214  // Init EMCAL geometry
215  if (!fgeom) SetEMCalGeometry();
216  //Init EMCAL geometry done
217 
218  return kTRUE;
219 }
220 
221 //_____________________________________________________________________
224 {
225  AliDebug(1,"AliAnalysisTaskEMCALTimeCalib::SetEMCalGeometry()");
226  if(fGeometryName.Length()==0){
227  fgeom=AliEMCALGeometry::GetInstanceFromRunNumber(fRunNumber);
228  AliInfo(Form("Get EMCAL geometry name <%s> for run %d",fgeom->GetName(),fRunNumber));
229  } else {
230  fgeom = AliEMCALGeometry::GetInstance(fGeometryName.Data());
231  AliInfo(Form("Set EMCAL geometry name to <%s>",fGeometryName.Data()));
232  }
233 
234  if (!fgeom){
235  AliWarning("Make sure the EMCal geometry is set properly !");
236  } else {
237  AliDebug(1,Form("EMCal geometry properly set: fGeom = %p, fGeometryName=%s",fgeom,fGeometryName.Data()));
238  }
239 
240  return kTRUE;
241 }
242 
243 //_____________________________________________________________________
246 {
247  //method under development
248  AliInfo(Form("<D> -- Run # = %d", fRunNumber));
249  AliInfo("prepare TOFT0maker!!");
250  //cout<<"Run "<< fRunNumber<<" in TOFT0maker"<<endl;
251 
252 
253  AliCDBManager * cdb = AliCDBManager::Instance();
254  cdb->SetDefaultStorage("raw://");
255  cdb->SetRun(fRunNumber);
256 
257  AliESDpid *extPID=new AliESDpid();
258 
259  // Wonder if some have to be declared as private variables??
260  // AliESDpid *extPID = new AliESDpid();
261  // AliTOFcalib * tofCalib = new AliTOFcalib();
262  // tofCalib->SetCalibrateTOFsignal(kTRUE);
263  // tofCalib->Init();
264 
265  fTOFmaker = new AliTOFT0maker(extPID);
266  fTOFmaker->SetTimeResolution(115.0); // if you want set the TOF res
267  // fTOFmaker = new AliTOFT0maker(extPID,tofCalib);
268  // fTOFmaker->SetTimeResolution(130.0);
269 
270  //cout<<"extPID "<<extPID<<" fTOFmaker "<<fTOFmaker<<endl;
271 
272 }// End PrepareTOFT0maker
273 
274 //________________________________________________________________________
278 {
279  AliDebug(1,"AliAnalysisTaskEMCALTimeCalib::UserCreateOutputObjects()");
280 
281  Double_t fineTmin = -500;
282  Double_t fineTmax = 400;
283  Double_t fineInterval = 0.20;
284  Int_t nfinebin = (fineTmax-fineTmin)/fineInterval;
285  //cout << "<D> nfinebin = " << nfinebin << endl;
286 
287  Int_t nChannels = 17664;
288  //book histograms
289  fhcalcEvtTime = new TH1F("fhcalcEvtTime","calculated event time from T0",nfinebin, fineTmin,fineTmax);
290  fhcalcEvtTime->GetXaxis()->SetTitle("T ");
291  fhcalcEvtTime->GetYaxis()->SetTitle("Counts (a.u.)");
292 
293  fhEvtTimeHeader = new TH1F("fhEvtTimeHeader","event time from header",nfinebin, fineTmin,fineTmax);
294  fhEvtTimeHeader->GetXaxis()->SetTitle("T ");
295  fhEvtTimeHeader->GetYaxis()->SetTitle("Counts (a.u.)");
296 
297  fhEvtTimeDiff = new TH1F("fhEvtTimeDiff","event time difference",nfinebin, fineTmin,fineTmax);
298  fhEvtTimeDiff->GetXaxis()->SetTitle("#Delta T ");
299  fhEvtTimeDiff->GetYaxis()->SetTitle("Counts (a.u.)");
300 
301  fhEventType = new TH1F("fhEventType","event type",10, 0.,10.);
302  fhEventType ->GetXaxis()->SetTitle("Type ");
303  fhEventType ->GetYaxis()->SetTitle("Counts (a.u.)");
304  fhTcellvsTOFT0 = new TH2F("hTcellvsTOFT0", " T_cell vs TOFT0", 500,-600.0,+400.0,1200,300.0,900.0);
305  fhTcellvsTOFT0HD = new TH2F("hTcellvsTOFT0HD", " T_cell vs TOFT0,HighEnergy", 500,-600.0,+400.0,4000,500.0,700.0);
306  fhTcellvsSM = new TH2F("hTcellvsSM", " T_cell vs SM", 20,0,20,300,300,900);
307  fhEneVsAbsIdHG = new TH2F("fhEneVsAbsIdHG", "energy vs ID for HG",1000,0,18000,200,0,10);
308  fhEneVsAbsIdLG = new TH2F("fhEneVsAbsIdLG", "energy vs ID for LG",1000,0,18000,100,0,20);
309 
310  for (Int_t i = 0; i < kNBCmask ; i++)
311  {
312  //already after correction
313  //high gain
314  fhTimeSumSq[i] = new TH1F(Form("hTimeSumSq%d", i),
315  Form("cell Sum Square time HG, BC %d ", i),
316  nChannels,0.,(Double_t)nChannels);
317  fhTimeSumSq[i]->SetYTitle("Sum Sq Time ");
318  fhTimeSumSq[i]->SetXTitle("AbsId");
319 
320  fhTimeSum[i] = new TH1F(Form("hTimeSum%d", i),
321  Form("cell Sum time HG, BC %d ", i),
322  nChannels,0.,(Double_t)nChannels);
323  fhTimeSum[i]->SetYTitle("Sum Time ");
324  fhTimeSum[i]->SetXTitle("AbsId");
325 
326  fhTimeEnt[i] = new TH1F(Form("hTimeEnt%d", i),
327  Form("cell Entries HG, BC %d ", i),
328  nChannels,0.,(Double_t)nChannels);
329  fhTimeEnt[i]->SetYTitle("Entries for Time ");
330  fhTimeEnt[i]->SetXTitle("AbsId");
331 
332  //low gain
333  fhTimeLGSumSq[i] = new TH1F(Form("hTimeLGSumSq%d", i),
334  Form("cell Sum Square time LG, BC %d ", i),
335  nChannels,0.,(Double_t)nChannels);
336  fhTimeLGSumSq[i]->SetYTitle("Sum Sq Time ");
337  fhTimeLGSumSq[i]->SetXTitle("AbsId");
338 
339  fhTimeLGSum[i] = new TH1F(Form("hTimeLGSum%d", i),
340  Form("cell Sum time LG, BC %d ", i),
341  nChannels,0.,(Double_t)nChannels);
342  fhTimeLGSum[i]->SetYTitle("Sum Time ");
343  fhTimeLGSum[i]->SetXTitle("AbsId");
344 
345  fhTimeLGEnt[i] = new TH1F(Form("hTimeLGEnt%d", i),
346  Form("cell Entries LG, BC %d ", i),
347  nChannels,0.,(Double_t)nChannels);
348  fhTimeLGEnt[i]->SetYTitle("Entries for Time ");
349  fhTimeLGEnt[i]->SetXTitle("AbsId");
350 
351  //raw time histograms
352  //high gain
353  fhRawTimeVsIdBC[i] = new TH2F(Form("RawTimeVsIdBC%d", i),
354  Form("cell raw time vs ID for high gain BC %d ", i),
355  nChannels,0.,(Double_t)nChannels,600,300,900);
356  fhRawTimeVsIdBC[i]->SetXTitle("AbsId");
357  fhRawTimeVsIdBC[i]->SetYTitle("Time");
358 
359  fhRawTimeSumBC[i] = new TH1F(Form("RawTimeSumBC%d", i),
360  Form("sum of cell raw time for high gain BC %d ", i),
361  nChannels,0.,(Double_t)nChannels);
362  fhRawTimeSumBC[i]->SetXTitle("AbsId");
363  fhRawTimeSumBC[i]->SetYTitle("Sum Time");
364 
365  fhRawTimeEntriesBC[i] = new TH1F(Form("RawTimeEntriesBC%d", i),
366  Form("No. entries of cells raw time for high gain BC %d ", i),
367  nChannels,0.,(Double_t)nChannels);
368  fhRawTimeEntriesBC[i]->SetXTitle("AbsId");
369  fhRawTimeEntriesBC[i]->SetYTitle("Entries for Time ");
370 
371  fhRawTimeSumSqBC[i] = new TH1F(Form("RawTimeSumSqBC%d", i),
372  Form("sum of (cell raw time)^2 for high gain BC %d ", i),
373  nChannels,0.,(Double_t)nChannels);
374  fhRawTimeSumSqBC[i]->SetXTitle("AbsId");
375  fhRawTimeSumSqBC[i]->SetYTitle("Sum Sq Time");
376 
377  //low gain
378  fhRawTimeVsIdLGBC[i] = new TH2F(Form("RawTimeVsIdLGBC%d", i),
379  Form("cell raw time vs ID for low gain BC %d ", i),
380  nChannels,0.,(Double_t)nChannels,600,300,900);
381  fhRawTimeVsIdLGBC[i]->SetXTitle("AbsId");
382  fhRawTimeVsIdLGBC[i]->SetYTitle("Time");
383 
384  fhRawTimeSumLGBC[i] = new TH1F(Form("RawTimeSumLGBC%d", i),
385  Form("sum of cell raw time for low gain BC %d ", i),
386  nChannels,0.,(Double_t)nChannels);
387  fhRawTimeSumLGBC[i]->SetXTitle("AbsId");
388  fhRawTimeSumLGBC[i]->SetYTitle("Sum Time");
389 
390  fhRawTimeEntriesLGBC[i] = new TH1F(Form("RawTimeEntriesLGBC%d", i),
391  Form("No. entries of cells raw time for low gain BC %d ", i),
392  nChannels,0.,(Double_t)nChannels);
393  fhRawTimeEntriesLGBC[i]->SetXTitle("AbsId");
394  fhRawTimeEntriesLGBC[i]->SetYTitle("Entries for Time ");
395 
396  fhRawTimeSumSqLGBC[i] = new TH1F(Form("RawTimeSumSqLGBC%d", i),
397  Form("sum of (cell raw time)^2 for low gain BC %d ", i),
398  nChannels,0.,(Double_t)nChannels);
399  fhRawTimeSumSqLGBC[i]->SetXTitle("AbsId");
400  fhRawTimeSumSqLGBC[i]->SetYTitle("Sum Sq Time");
401 
402 
403  for (Int_t j = 0; j < kNSM ; j++)
404  {
405  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);
406  fhTimeDsupBC[j][i]->SetYTitle(" Time (ns) ");
407  fhTimeDsupBC[j][i]->SetXTitle(" E (GeV) ");
408  }
409  }
410 
411  for (Int_t jj = 0; jj < kNSM ; jj++)
412  {
413  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);
414  fhTimeDsup[jj]->SetYTitle(" Time (ns) ");
415  fhTimeDsup[jj]->SetXTitle(" E (GeV) ");
416  }
417 
418  fhTimeVsBC = new TH2F("TimeVsBC"," SupMod time_vs_BC ", 4001,-0.5,4000.5,400,200.0,1000.0);
419 
420 
421  //add histos to list
422  fOutputList = new TList();
423 
427  fOutputList->Add(fhEventType);
430  fOutputList->Add(fhTcellvsSM);
433 
434  for (Int_t i = 0; i < kNBCmask ; i++)
435  {
436  fOutputList->Add(fhTimeSumSq[i]);
437  fOutputList->Add(fhTimeEnt[i]);
438  fOutputList->Add(fhTimeSum[i]);
439 
440  fOutputList->Add(fhTimeLGSumSq[i]);
441  fOutputList->Add(fhTimeLGEnt[i]);
442  fOutputList->Add(fhTimeLGSum[i]);
443 
444  fOutputList->Add(fhRawTimeVsIdBC[i]);
445  fOutputList->Add(fhRawTimeSumBC[i]);
447  fOutputList->Add(fhRawTimeSumSqBC[i]);
448 
450  fOutputList->Add(fhRawTimeSumLGBC[i]);
453 
454  for (Int_t j = 0; j < kNSM ; j++){
455  fOutputList->Add(fhTimeDsupBC[j][i]);
456  }
457  }
458 
459  for (Int_t j = 0; j < kNSM ; j++)
460  {
461  fOutputList->Add(fhTimeDsup[j]);
462  }
463 
464  fOutputList->Add(fhTimeVsBC);
465 
466  fOutputList->SetOwner(kTRUE);
467  PostData(1,fOutputList);
468 
469 
470 } // End of AliAnalysisTaskEMCALTimeCalib::UserCreateOuputObjects()
471 
472 //________________________________________________________________________
475 {
476  // Called for each event
477  AliDebug(2,Form("UserExec: EMCal geometry: fgeom = %p fGeometryName %s",fgeom,fGeometryName.Data()));
478  AliVEvent *event = InputEvent();
479  //cout<<"T0TOF "<<event->GetT0TOF()<<endl;//bad idea
480  //cout<< fEvent->GetTOFHeader()->GetDefaultEventTimeVal()<<endl;
481  AliDebug(2,Form("TOF time from header %f ps",event->GetTOFHeader()->GetDefaultEventTimeVal()));
482  fhEvtTimeHeader->Fill(event->GetTOFHeader()->GetDefaultEventTimeVal());
483 
484  //fEvent = dynamic_cast<AliESDEvent*>(event);
485  if (!event) {
486  AliError("ESD not available, exit");
487  fhEventType->Fill(0.5);
488  return;
489  }
490 
491  if(fPileupFromSPD==kTRUE){
492  if(event->IsPileupFromSPD(3,0.8,3.,2.,5.)){
493  AliDebug(1,"Event: PileUp skip.");
494  fhEventType->Fill(1.5);
495  return;
496  }
497  }
498 
499  TString triggerclasses = event->GetFiredTriggerClasses();
500  if(triggerclasses=="") {
501  fhEventType->Fill(2.5);
502  return;
503  }
504 
505  Int_t eventType = ((AliVHeader*)event->GetHeader())->GetEventType();
506  // physics events eventType=7, select only those
507  AliDebug(1,Form("Triggerclasses %s, eventType %d",triggerclasses.Data(),eventType));
508  if(eventType != 7) {
509  fhEventType->Fill(3.5);
510  return;
511  }
512 
513  // Check trigger
514  Bool_t bMB = kFALSE;
515  Bool_t bL0 = kFALSE;
516  Bool_t bL1G = kFALSE;
517  Bool_t bL1J = kFALSE;
518 
519  if(triggerclasses.Contains("CINT7-B-NOPF-ALLNOTRD") ||
520  triggerclasses.Contains("CINT7-I-NOPF-ALLNOTRD") ||
521  triggerclasses.Contains("CINT1-I-NOPF-ALLNOTRD") ||
522  triggerclasses.Contains("CINT1-B-NOPF-ALLNOTRD") ||
523  triggerclasses.Contains("CINT8") ||
524  triggerclasses.Contains("CINT7") ||
525  triggerclasses.Contains("CPBI2_B1-B-NOPF-ALLNOTRD") ) bMB = kTRUE;
526 
527  if(triggerclasses.Contains("CEMC7-B-NOPF-CENTNOTRD") ||
528  triggerclasses.Contains("CEMC1-B-NOPF-CENTNOTRD") ||
529  triggerclasses.Contains("CEMC7") ||
530  triggerclasses.Contains("CEMC8") ||
531  triggerclasses.Contains("CEMC8-B-NOPF-CENTNOTRD") ) bL0 = kTRUE;
532 
533  if(triggerclasses.Contains("CEMC7EG1-B-NOPF-CENTNOTRD") ||
534  triggerclasses.Contains("CEMC7EG2-B-NOPF-CENTNOTRD") ||
535  triggerclasses.Contains("CEMC8EG1-B-NOPF-CENTNOTRD") ||
536  triggerclasses.Contains("CEMC8EGA") ||
537  triggerclasses.Contains("CEMC7EGA") ||
538  triggerclasses.Contains("CPBI2EGA") ) bL1G = kTRUE;
539 
540 
541  if(triggerclasses.Contains("CEMC7EJ1-B-NOPF-CENTNOTRD") ||
542  triggerclasses.Contains("CEMC7EJ2-B-NOPF-CENTNOTRD") ||
543  triggerclasses.Contains("CEMC8EJ1-B-NOPF-CENTNOTRD") ||
544  triggerclasses.Contains("CEMC7EJE") ||
545  triggerclasses.Contains("CEMC8EJE") ||
546  triggerclasses.Contains("CPBI2EJE") ) bL1J = kTRUE;
547 
548  if( bL1G || bL1J || bL0 ){ fhEventType->Fill(4.5);}
549  if( bMB ){ fhEventType->Fill(5.5);}
550 
551 
552  // if(bL1G || bL1J || bL0){
553 
554  // Prepare TOFT0 maker at the beginning of a run
555 // if (event->GetRunNumber() != fRunNumber)
556 // {
557 // fRunNumber = event->GetRunNumber();
558 // // PrepareTOFT0maker();
559 // // cout<<"tofT0maker per run"<<fRunNumber<<endl;
560 // }// fi Check if run number has changed
561 
562  // --- Use of AliTOFT0maker
563  Double_t calcolot0=0.0;
564  if(!AODEvent()){
565  Double_t* timeTOFtable;
566  timeTOFtable=fTOFmaker->ComputeT0TOF(dynamic_cast<AliESDEvent*>(event));
567  AliDebug(2,Form("TOF time %f ps, resolution %f ps, tracks at TOF %f/used %f",timeTOFtable[0],timeTOFtable[1],timeTOFtable[3],timeTOFtable[7]));
568  //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;
569  calcolot0=timeTOFtable[0];
570  }
571 
572  if (!fhcalcEvtTime) {
573  AliWarning("<E> fhcalcEvtTime not available");
574  return;
575  }// fi no simple histo present
576 
577  fhcalcEvtTime->Fill(calcolot0);
578  if(calcolot0 != 0 && event->GetTOFHeader()->GetDefaultEventTimeVal() != 0 )
579  fhEvtTimeDiff->Fill(calcolot0-event->GetTOFHeader()->GetDefaultEventTimeVal());
580 
581  TRefArray* caloClusters = new TRefArray();
582  event->GetEMCALClusters(caloClusters);
583  // cout << " ###########Bunch Cross nb = " << event->GetBunchCrossNumber() << endl;
584 
585  Int_t BunchCrossNumber =event->GetBunchCrossNumber();
586 
587  Float_t offset=0.;
588 
589  Int_t nBC = 0;
590  nBC = BunchCrossNumber%4;
591  //Int_t nTriggerMask =event->GetTriggerMask();
592  // cout << " nBC " << nBC << " nTriggerMask " << nTriggerMask<< endl;
593  Float_t timeBCoffset = 0.; //manual offest
594  // if( nBC%4 ==0 || nBC%4==1) timeBCoffset = 100.; // correction was for LHC11 when BC was not corrected
595 
596  Int_t nclus = caloClusters->GetEntries();
597  AliDebug(1,Form("###########Bunch Cross nb = %d nclus = %d",nBC,nclus ));
598  //cout << " ###########Bunch Cross nb = " << nBC <<" nclus= "<< nclus<< endl;
599  //Int_t ntracks = event-> GetNumberOfTracks() ;
600 
601  AliVCaloCells &cells= *(event->GetEMCALCells());//it is cluster independent
602  //Variables used plenty in loops
603  Int_t nSupMod=-1, nModule=-1;
604  Int_t iphi=-1, ieta=-1, nIphi=-1, nIeta=-1;
605  Int_t absId=-1;
606  Float_t hkdtime=0.0;
607  Float_t amp=0.0;
608  Bool_t isHighGain=kTRUE;
609 
610  for (Int_t icl = 0; icl < nclus; icl++) {
611  //ESD and AOD CaloCells carries the same information
612  AliVCluster* clus = (AliVCluster*)caloClusters->At(icl);
613  if(!AcceptCluster(clus)) continue;
614 
615  //cout<<"nCells="<< clus->GetNCells();<<endl;
616 
617  UShort_t * index = clus->GetCellsAbsId() ;
618 
619  for(Int_t i = 0; i < clus->GetNCells() ; i++) {
620  absId = index[i]; // or clus->GetCellNumber(i) ;
621  hkdtime = cells.GetCellTime(absId) * 1.0e09; // to get ns
622  amp = cells.GetCellAmplitude(absId) ;
623  isHighGain = cells.GetCellHighGain(absId);
624  //cout<<"cell absID: "<<absId<<" cellTime: "<<hkdtime<<" cellaplit: "<< amp<<endl;
625 
626  //main histograms with raw time information
627  if(amp>fMinCellEnergy){
628  if(isHighGain){
629  fhRawTimeVsIdBC[nBC]->Fill(absId,hkdtime);
630  fhRawTimeSumBC[nBC]->Fill(absId,hkdtime);
631  fhRawTimeEntriesBC[nBC]->Fill(absId,1.);
632  fhRawTimeSumSqBC[nBC]->Fill(absId,hkdtime*hkdtime);
633  }else{
634  fhRawTimeVsIdLGBC[nBC]->Fill(absId,hkdtime);
635  fhRawTimeSumLGBC[nBC]->Fill(absId,hkdtime);
636  fhRawTimeEntriesLGBC[nBC]->Fill(absId,1.);
637  fhRawTimeSumSqLGBC[nBC]->Fill(absId,hkdtime*hkdtime);
638  }
639  }
640  //fgeom->PrintCellIndexes(absId);
641  //fgeom->PrintCellIndexes(absId,1);
642 
643  // GEOMETRY tranformations
644  fgeom->GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
645  fgeom->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
646 
647  // other histograms for cross-check
648  CheckCellRCU(nSupMod,ieta,iphi);//SM, column, row
649 
650  fhTcellvsSM->Fill(nSupMod,hkdtime);
651  if(isHighGain==kTRUE) {fhEneVsAbsIdHG->Fill(absId,amp);}
652  else {fhEneVsAbsIdLG->Fill(absId,amp);}
653 
654  fhTimeVsBC->Fill(1.*BunchCrossNumber,hkdtime-timeBCoffset);
655  if(isHighGain==kTRUE){
656  if(fhAllAverageBC[nBC]!=0) {//comming from file after the first iteration
657  offset = (Float_t)(fhAllAverageBC[nBC]->GetBinContent(absId+1));//channel absId=0 has histogram bin=1
658  } else if(fReferenceFileName.Length()!=0){//protection against missing reference histogram
659  AliFatal(Form("Reference histogram for BC%d not properly loaded",nBC));
660  }
661  } else {
662  if(fhAllAverageLGBC[nBC]!=0) {//comming from file after the first iteration
663  offset = (Float_t)(fhAllAverageLGBC[nBC]->GetBinContent(absId+1));//channel absId=0 has histogram bin=1
664  } else if(fReferenceFileName.Length()!=0){//protection against missing reference histogram
665  AliFatal(Form("Reference LG histogram for BC%d not properly loaded",nBC));
666  }
667  }
668  //if(offset==0)cout<<"offset 0 in SM "<<nSupMod<<endl;
669 
670 
671  if(amp>0.5) {
672  fhTimeDsup[nSupMod]->Fill(amp,hkdtime-offset);
673  fhTimeDsupBC[nSupMod][nBC]->Fill(amp,hkdtime-offset);
674  }
675 
676  if(amp>0.9) {
677  fhTcellvsTOFT0HD->Fill(calcolot0, hkdtime);
678  }
679 
680  fhTcellvsTOFT0->Fill(calcolot0, hkdtime-offset);
681 
682  hkdtime = hkdtime-timeBCoffset;//time corrected by manual offset (default=0)
683  Float_t hkdtimecorr;
684  hkdtimecorr= hkdtime-offset;//time after first iteration
685 
686  //main histograms after the first itereation for calibration constants
687  //if(hkdtimecorr>=-20. && hkdtimecorr<=20. && amp>0.9 ) {
688  if(hkdtimecorr>=fMinTime && hkdtimecorr<=fMaxTime && amp>fMinCellEnergy ) {
689  // per cell
690 // Float_t entriesTime=fhTimeEnt[nBC]->GetBinContent(absId)+1;
691 // Float_t sumTimeSq=(fhTimeSumSq[nBC]->GetBinContent(absId)+(hkdtime*hkdtime));
692 // Float_t sumTime=(fhTimeSum[nBC]->GetBinContent(absId)+hkdtime);
693 //
694 // fhTimeEnt[nBC]->SetBinContent(absId,entriesTime);
695 // fhTimeSumSq[nBC]->SetBinContent(absId,sumTimeSq);
696 // fhTimeSum[nBC]->SetBinContent(absId,sumTime);
697 
698  if(isHighGain){
699  fhTimeEnt[nBC]->Fill(absId,1.);
700  fhTimeSumSq[nBC]->Fill(absId,hkdtime*hkdtime);
701  fhTimeSum[nBC]->Fill(absId,hkdtime);
702  }else{
703  fhTimeLGEnt[nBC]->Fill(absId,1.);
704  fhTimeLGSumSq[nBC]->Fill(absId,hkdtime*hkdtime);
705  fhTimeLGSum[nBC]->Fill(absId,hkdtime);
706  }
707 
708 
709  } // hkdtime:[-20;20]
710  } // end icell
711  } //end cluster
712 
713 
714  // Post output data.
715  //cout<<"Post data and delete caloClusters"<<endl;
716  caloClusters->Delete();
717  delete caloClusters;
718 // } // end if trigger type
719 
720  PostData(1, fOutputList);
721 } // End of AliAnalysisTaskEMCALTimeCalib::UserExec()
722 
723 //________________________________________________________________________
727 {
728  fOutputList = dynamic_cast<TList*> (GetOutputData(1));
729 
730  if(fTOFmaker) delete fTOFmaker;
731 
732 
733  if (!fOutputList)
734  {
735  AliDebug(1,"ERROR: Output list not available");
736  return;
737  }
738 } // End of AliAnalysisTaskEMCALTimeCalib::Terminate
739 
740 //________________________________________________________________________
743 {
744  //fix with noisy EMCAL fee card
745  Int_t nCells = clus->GetNCells();
746 
747  if(clus->IsEMCAL())
748  {
749  if ((clus->E() > fMaxClusterEnergy && nCells > fMaxNcells ) || nCells > fMaxNcells){
750  AliDebug(1,"very big cluster with enormous energy - cluster rejected");
751  return kFALSE;
752  }
753  }
754 
755  // remove other than photonlike
756  Double_t lambda0=clus->GetM02();
757  if (lambda0>fMaxLambda0LG || lambda0<fMinLambda0LG){
758  AliDebug(1,"lambda0 loose cut failed - cluster rejected");
759  return kFALSE;
760  }
761 
762  // remove matched clusters
763  Double_t Dx=clus->GetTrackDx();
764  Double_t Dz=clus->GetTrackDz();
765  Double_t Rtrack = TMath::Sqrt(Dx*Dx+Dz*Dz);
766  if (Rtrack <fMaxRtrack)
767  {
768  AliDebug(1,"track matched - cluster rejected");
769  return kFALSE;
770  }
771 
772  if (nCells<fMinNcells)
773  {
774  AliDebug(1,"single cell cluster - cluster rejected");
775  return kFALSE;
776  }
777 
778  if(clus->E()<fMinClusterEnergy)
779  {
780  AliDebug(1,"cluster energy < 1 GeV- cluster rejected");
781  return kFALSE;
782  }
783 
784 
785  if(!IsLowGainCellInCluster(clus)) {//no low gain cell in cluster
786  //apply more strict lambda0^2 cut
787  if (lambda0>fMaxLambda0 || lambda0<fMinLambda0){
788  AliDebug(1,"lambda0 strict cut failed - cluster rejected");
789  return kFALSE;
790  }
791  }
792 
793 
794 
795 
796  return kTRUE;
797 }//End AliAnalysisTaskEMCALTimeCalib::AcceptCluster
798 
799 //________________________________________________________________________
802  UShort_t * index = clus->GetCellsAbsId() ;
803  AliVCaloCells &cells= *(InputEvent()->GetEMCALCells());
804  for(Int_t i = 0; i < clus->GetNCells() ; i++) {
805  if(cells.GetCellHighGain(index[i])==kFALSE) return kTRUE;//low gain cell
806  }
807  return kFALSE;
808 
809 }
810 
811 //________________________________________________________________________
813 Bool_t AliAnalysisTaskEMCALTimeCalib::CheckCellRCU(Int_t nSupMod,Int_t icol,Int_t irow)
814 {
815  Int_t iRCU;
816  if(nSupMod < 10 || (nSupMod >= 12 && nSupMod <18) )
817  {
818  if (0<=irow&&irow<8) iRCU=0; // first cable row
819  else if (8<=irow&&irow<16 && 0<=icol&&icol<24) iRCU=0; // first half;
820  //second cable row
821  //RCU1
822  else if (8<=irow&&irow<16 && 24<=icol&&icol<48) iRCU=1; // second half;
823  //second cable row
824  else if (16<=irow&&irow<24) iRCU=1; // third cable row
825 
826  if (nSupMod%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
827  }
828  else
829  {
830  // Last 2 SM have one single SRU, just assign RCU 0
831  iRCU = 0 ;
832  }
833 
834  //cout<<"RCU:"<<iRCU<<endl;
835  if (iRCU<0)
836  AliFatal(Form("Wrong EMCAL/DCAL RCU number = %d\n", iRCU));
837 
838  return kTRUE;
839 }//End AliAnalysisTaskEMCALTimeCalib::CheckCellRCU
840 
841 //________________________________________________________________________
844 {
845  fMinClusterEnergy=1.0;//0.5//0.7
846  fMaxClusterEnergy=500;
847  fMinNcells=2;
848  fMaxNcells=200;
849  fMinLambda0=0.1;
850  fMaxLambda0=0.4;
851  fMinLambda0LG=0.1;
852  fMaxLambda0LG=4.0;
853  fMaxRtrack=0.025;
854  fMinCellEnergy=0.4;//0.1//0.4
855  fReferenceFileName="";//Reference.root
856  fGeometryName="";//EMCAL_COMPLETE12SMV1_DCAL_8SM
857  fPileupFromSPD=kFALSE;
858  fMinTime=-20.;
859  fMaxTime=20.;
860 }
861 
862 //________________________________________________________________________
867 void AliAnalysisTaskEMCALTimeCalib::ProduceCalibConsts(TString inputFile,TString outputFile,Bool_t isFinal)
868 {
869  TFile *file =new TFile(inputFile.Data());
870  if(file==0x0) {
871  //AliWarning("Input file does not exist!");
872  return;
873  }
874 
875  TList *list=(TList*)file->Get("chistolist");
876  if(list==0x0)
877  {
878  //AliWarning("List chistolist does not exist in file!");
879  return;
880  }
881 
882  //high gain
883  TH1F *h1[4];
884  TH1F *h2[4];
885  TH1F *h3[4];
886  TH1F *hAllTimeAvBC[4];
887  TH1F *hAllTimeRMSBC[4];
888 
889  //low gain
890  TH1F *h4[4];
891  TH1F *h5[4];
892  TH1F *h6[4];
893  TH1F *hAllTimeAvLGBC[4];
894  TH1F *hAllTimeRMSLGBC[4];
895 
896  if(isFinal==kFALSE){//first itereation
897  for(Int_t i=0;i<4;i++){
898  h1[i]=(TH1F *)list->FindObject(Form("RawTimeSumBC%d",i));
899  h2[i]=(TH1F *)list->FindObject(Form("RawTimeEntriesBC%d",i));
900  h3[i]=(TH1F *)list->FindObject(Form("RawTimeSumSqBC%d",i));
901 
902  h4[i]=(TH1F *)list->FindObject(Form("RawTimeSumLGBC%d",i));
903  h5[i]=(TH1F *)list->FindObject(Form("RawTimeEntriesLGBC%d",i));
904  h6[i]=(TH1F *)list->FindObject(Form("RawTimeSumSqLGBC%d",i));
905  }
906  } else {//final iteration
907  for(Int_t i=0;i<4;i++){
908  h1[i]=(TH1F *)list->FindObject(Form("hTimeSum%d",i));
909  h2[i]=(TH1F *)list->FindObject(Form("hTimeEnt%d",i));
910  h3[i]=(TH1F *)list->FindObject(Form("hTimeSumSq%d",i));
911 
912  h4[i]=(TH1F *)list->FindObject(Form("hTimeLGSum%d",i));
913  h5[i]=(TH1F *)list->FindObject(Form("hTimeLGEnt%d",i));
914  h6[i]=(TH1F *)list->FindObject(Form("hTimeLGSumSq%d",i));
915  }
916  }
917  //AliWarning("Input histograms read.");
918 
919  for(Int_t i=0;i<4;i++){
920  hAllTimeAvBC[i]=new TH1F(Form("hAllTimeAvBC%d",i),Form("hAllTimeAvBC%d",i),h1[i]->GetNbinsX(),h1[i]->GetXaxis()->GetXmin(),h1[i]->GetXaxis()->GetXmax());
921  hAllTimeRMSBC[i]=new TH1F(Form("hAllTimeRMSBC%d",i),Form("hAllTimeRMSBC%d",i),h3[i]->GetNbinsX(),h3[i]->GetXaxis()->GetXmin(),h3[i]->GetXaxis()->GetXmax());
922 
923  hAllTimeAvLGBC[i]=new TH1F(Form("hAllTimeAvLGBC%d",i),Form("hAllTimeAvLGBC%d",i),h4[i]->GetNbinsX(),h4[i]->GetXaxis()->GetXmin(),h4[i]->GetXaxis()->GetXmax());
924  hAllTimeRMSLGBC[i]=new TH1F(Form("hAllTimeRMSLGBC%d",i),Form("hAllTimeRMSLGBC%d",i),h6[i]->GetNbinsX(),h6[i]->GetXaxis()->GetXmin(),h6[i]->GetXaxis()->GetXmax());
925  }
926 
927  //AliWarning("New histograms booked.");
928 
929  for(Int_t i=0;i<4;i++){
930  for(Int_t j=1;j<=h1[i]->GetNbinsX();j++){
931  //high gain
932  if(h2[i]->GetBinContent(j)!=0){
933  hAllTimeAvBC[i]->SetBinContent(j,h1[i]->GetBinContent(j)/h2[i]->GetBinContent(j));
934  hAllTimeRMSBC[i]->SetBinContent(j,TMath::Sqrt(h3[i]->GetBinContent(j)/h2[i]->GetBinContent(j)) );
935  } else {
936  hAllTimeAvBC[i]->SetBinContent(j,0.);
937  hAllTimeRMSBC[i]->SetBinContent(j,0.);
938  }
939  //low gain
940  if(h5[i]->GetBinContent(j)!=0){
941  hAllTimeAvLGBC[i]->SetBinContent(j,h4[i]->GetBinContent(j)/h5[i]->GetBinContent(j));
942  hAllTimeRMSLGBC[i]->SetBinContent(j,TMath::Sqrt(h6[i]->GetBinContent(j)/h5[i]->GetBinContent(j)) );
943  } else {
944  hAllTimeAvLGBC[i]->SetBinContent(j,0.);
945  hAllTimeRMSLGBC[i]->SetBinContent(j,0.);
946  }
947 
948  }
949  }
950 
951  //AliWarning("Average and rms calculated.");
952  TFile *fileNew=new TFile(outputFile.Data(),"recreate");
953  for(Int_t i=0;i<4;i++){
954  hAllTimeAvBC[i]->Write();
955  hAllTimeRMSBC[i]->Write();
956  hAllTimeAvLGBC[i]->Write();
957  hAllTimeRMSLGBC[i]->Write();
958  }
959 
960  //AliWarning(Form("Histograms saved in %s file.",outputFile.Data()));
961 
962  fileNew->Close();
963  delete fileNew;
964  file->Close();
965  delete file;
966 
967  //AliWarning("Pointers deleted. Memory cleaned.");
968 
969 }
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