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