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