AliPhysics  251aa1e (251aa1e)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliEmcalPatchFromCellMaker.cxx
Go to the documentation of this file.
1 // $Id$
2 //
3 // Class to put cells into trigger patches
4 //
5 // Author: M. Verweij
6 
7 #include <TClonesArray.h>
8 #include <TRandom3.h>
9 #include <TProfile.h>
10 #include <TH3F.h>
11 
12 #include "AliAnalysisManager.h"
13 #include "AliLog.h"
14 #include "AliEMCALGeometry.h"
15 #include "AliEMCALTriggerPatchInfo.h"
16 #include "AliEMCALTriggerBitConfig.h"
17 #include "AliEMCALTriggerConstants.h"
18 
20 
22 
23 //________________________________________________________________________
26  fCaloTriggersOutName("EmcalPatches32x32"),
27  fCaloTriggersOut(0),
28  fPatchDim(32),
29  fMinCellE(0.05),
30  fCellTimeMin(485e-9),
31  fCellTimeMax(685e-9),
32  fL1Slide(0),
33  fTriggerBitConfig(0x0),
34  fh3EEtaPhiCell(0),
35  fh2CellEnergyVsTime(0),
36  fh1CellEnergySum(0)
37 {
38  // Constructor.
39  for (Int_t i = 0; i < kPatchCols; i++) {
40  for (Int_t j = 0; j < kPatchRows; j++) {
41  fPatchADCSimple[i][j] = 0.;
42  fPatchESimple[i][j] = 0.;
43  }
44  }
45 
46  SetMakeGeneralHistograms(kTRUE);
47 }
48 
49 //________________________________________________________________________
51  AliAnalysisTaskEmcal(name,kTRUE),
52  fCaloTriggersOutName("EmcalPatches32x32"),
53  fCaloTriggersOut(0),
54  fPatchDim(32),
55  fMinCellE(0.05),
56  fCellTimeMin(485e-9),
57  fCellTimeMax(685e-9),
58  fL1Slide(0),
59  fTriggerBitConfig(0x0),
60  fh3EEtaPhiCell(0),
61  fh2CellEnergyVsTime(0),
62  fh1CellEnergySum(0)
63 {
64  // Constructor.
65  for (Int_t i = 0; i < kPatchCols; i++) {
66  for (Int_t j = 0; j < kPatchRows; j++) {
67  fPatchADCSimple[i][j] = 0.;
68  fPatchESimple[i][j] = 0.;
69  }
70  }
71 
73 }
74 
75 //________________________________________________________________________
77 {
78  // Destructor.
79 }
80 
81 //________________________________________________________________________
83 {
84  // Init the analysis.
85 
87 
88  if (!fLocalInitialized)
89  return;
90 
91  if (!fCaloTriggersOutName.IsNull()) {
92  fCaloTriggersOut = new TClonesArray("AliEMCALTriggerPatchInfo");
94 
95  if (!(InputEvent()->FindListObject(fCaloTriggersOutName))) {
96  InputEvent()->AddObject(fCaloTriggersOut);
97  }
98  else {
99  fLocalInitialized = kFALSE;
100  AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fCaloTriggersOutName.Data()));
101  return;
102  }
103  }
104 
105  if(!fTriggerBitConfig)
106  fTriggerBitConfig = new AliEMCALTriggerBitConfigNew();
107 
108 }
109 
110 //________________________________________________________________________
112 {
113  // Create user output.
114 
116 
117  Int_t fgkNPhiBins = 18*8;
118  Float_t kMinPhi = 0.;
119  Float_t kMaxPhi = 2.*TMath::Pi();
120  Double_t *binsPhi = new Double_t[fgkNPhiBins+1];
121  for(Int_t i=0; i<=fgkNPhiBins; i++) binsPhi[i]=(Double_t)kMinPhi + (kMaxPhi-kMinPhi)/fgkNPhiBins*(Double_t)i ;
122 
123  Int_t fgkNEtaBins = 100;
124  Float_t fgkEtaMin = -1.;
125  Float_t fgkEtaMax = 1.;
126  Double_t *binsEta=new Double_t[fgkNEtaBins+1];
127  for(Int_t i=0; i<=fgkNEtaBins; i++) binsEta[i]=(Double_t)fgkEtaMin + (fgkEtaMax-fgkEtaMin)/fgkNEtaBins*(Double_t)i ;
128 
129  Int_t fgkNTimeBins = 600;
130  Float_t kMinTime = -200.;
131  Float_t kMaxTime = 1000;
132  Double_t *binsTime = new Double_t[fgkNTimeBins+1];
133  for(Int_t i=0; i<=fgkNTimeBins; i++) binsTime[i]=(Double_t)kMinTime + (kMaxTime-kMinTime)/fgkNTimeBins*(Double_t)i ;
134 
135  Double_t enBinEdges[3][2];
136  enBinEdges[0][0] = 1.; //10 bins
137  enBinEdges[0][1] = 0.1;
138  enBinEdges[1][0] = 5.; //8 bins
139  enBinEdges[1][1] = 0.5;
140  enBinEdges[2][0] = 100.;//95 bins
141  enBinEdges[2][1] = 1.;
142 
143  const Float_t enmin1 = 0;
144  const Float_t enmax1 = enBinEdges[0][0];
145  const Float_t enmin2 = enmax1 ;
146  const Float_t enmax2 = enBinEdges[1][0];
147  const Float_t enmin3 = enmax2 ;
148  const Float_t enmax3 = enBinEdges[2][0];//fgkEnMax;
149  const Int_t nbin11 = (int)((enmax1-enmin1)/enBinEdges[0][1]);
150  const Int_t nbin12 = (int)((enmax2-enmin2)/enBinEdges[1][1])+nbin11;
151  const Int_t nbin13 = (int)((enmax3-enmin3)/enBinEdges[2][1])+nbin12;
152 
153  Int_t fgkNEnBins=nbin13;
154  Double_t *binsEn=new Double_t[fgkNEnBins+1];
155  for(Int_t i=0; i<=fgkNEnBins; i++) {
156  if(i<=nbin11) binsEn[i]=(Double_t)enmin1 + (enmax1-enmin1)/nbin11*(Double_t)i ;
157  if(i<=nbin12 && i>nbin11) binsEn[i]=(Double_t)enmin2 + (enmax2-enmin2)/(nbin12-nbin11)*((Double_t)i-(Double_t)nbin11) ;
158  if(i<=nbin13 && i>nbin12) binsEn[i]=(Double_t)enmin3 + (enmax3-enmin3)/(nbin13-nbin12)*((Double_t)i-(Double_t)nbin12) ;
159  }
160 
161  fh3EEtaPhiCell = new TH3F("fh3EEtaPhiCell","fh3EEtaPhiCell;E_{cell};#eta;#phi",fgkNEnBins,binsEn,fgkNEtaBins,binsEta,fgkNPhiBins,binsPhi);
162  fOutput->Add(fh3EEtaPhiCell);
163 
164  fh2CellEnergyVsTime = new TH2F("fh2CellEnergyVsTime","fh2CellEnergyVsTime;E_{cell};time",fgkNEnBins,binsEn,fgkNTimeBins,binsTime);
166 
167  fh1CellEnergySum = new TH1F("fh1CellEnergySum","fh1CellEnergySum;E_{cell};time",fgkNEnBins,binsEn);
169 
170  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
171 
172  if(binsEn) delete [] binsEn;
173  if(binsPhi) delete [] binsPhi;
174  if(binsEta) delete [] binsEta;
175  if(binsTime) delete [] binsTime;
176 }
177 
178 //________________________________________________________________________
180 {
181  // Main loop, called for each event.
182 
183  fCaloTriggersOut->Delete();
184 
185  if (!fCaloCells) {
186  AliError(Form("Calo cells container %s not available.", fCaloCellsName.Data()));
187  return kFALSE;
188  }
189 
190  for (Int_t i = 0; i < kPatchCols; i++) {
191  for (Int_t j = 0; j < kPatchRows; j++) {
192  fPatchADCSimple[i][j] = 0.;
193  fPatchESimple[i][j] = 0.;
194  }
195  }
196 
197  if(!FillPatchADCSimple()) {
198  AliError(Form("%s Could not create simple ADC patches",GetName()));
199  return kFALSE;
200  }
201 
203 
204  Double_t sum = 0.;
205  for (Int_t i = 0; i < kPatchCols; i++) {
206  for (Int_t j = 0; j < kPatchRows; j++) {
207  sum+=fPatchESimple[i][j];
208  }
209  }
210 
211  return kTRUE;
212 }
213 
214 //________________________________________________________________________
216 {
217 
218  // fill the array for offline trigger processing
219 
220  // fill the patch ADCs from cells
221  Double_t sum = 0.;
222  Int_t nCell = fCaloCells->GetNumberOfCells();
223  for(Int_t iCell = 0; iCell < nCell; ++iCell) {
224  // get the cell info, based in index in array
225  Short_t cellId = fCaloCells->GetCellNumber(iCell);
226 
227  Double_t cellT = fCaloCells->GetCellTime(cellId);
228  Double_t amp = fCaloCells->GetAmplitude(iCell);
229  fh2CellEnergyVsTime->Fill(amp,cellT*1e9);
230 
231  //timing cuts
232  if(cellT<fCellTimeMin || cellT>fCellTimeMax) continue;
233  //energy cut
234  if(amp<fMinCellE) continue;
235  sum+=amp;
236 
237  // get position
238  Int_t absId=-1;
239  fGeom->GetFastORIndexFromCellIndex(cellId, absId);
240  Int_t globCol=-1, globRow=-1;
241  fGeom->GetPositionInEMCALFromAbsFastORIndex(absId, globCol, globRow);
242  // add
243  fPatchADCSimple[globCol][globRow] += amp/EMCALTrigger::kEMCL1ADCtoGeV;
244  fPatchESimple[globCol][globRow] += amp;
245 
246  TVector3 pos;
247  fGeom->GetGlobal(cellId, pos);
248  TLorentzVector lv(pos,amp);
249  Double_t cellEta = lv.Eta();
250  Double_t cellPhi = lv.Phi();
251  if(cellPhi<0.) cellPhi+=TMath::TwoPi();
252  if(cellPhi>TMath::TwoPi()) cellPhi-=TMath::TwoPi();
253  fh3EEtaPhiCell->Fill(amp,cellEta,cellPhi);
254  }
255  fh1CellEnergySum->Fill(sum);
256 
257  return kTRUE;
258 }
259 
260 //________________________________________________________________________
262 {
263  // Runs a simple offline trigger algorithm.
264  // It creates separate patches with dimension fPatchDim
265 
266  // run the trigger algo, stepping by stepsize (in trigger tower units)
267  Int_t itrig = 0;
268  Int_t patchSize = GetDimFastor();
269  Int_t stepSize = GetSlidingStepSizeFastor();
270  Int_t maxCol = kPatchCols - patchSize;
271  // Int_t maxRow = kPatchRows - patchSize;
272  Int_t maxRow = fGeom->GetNTotalTRU()*2 - patchSize; //apparently this is the total TRU in phi ??
273  // Printf("fGeom->GetNTotalTRU(): %d %d = %d x %d ; %d x %d",fGeom->GetNTRU(),fGeom->GetNTotalTRU(),fGeom->GetNTRUEta(),fGeom->GetNTRUPhi(),fGeom->GetNModulesInTRUEta(),fGeom->GetNModulesInTRUPhi());
274 
275  for (Int_t i = 0; i <= maxCol; i += stepSize) {
276  for (Int_t j = 0; j <= maxRow; j += stepSize) {
277  // get the trigger towers composing the patch
278  Int_t adcAmp = 0;
279  Double_t enAmp = 0.;
280  // window
281  for (Int_t k = 0; k < patchSize; ++k) {
282  for (Int_t l = 0; l < patchSize; ++l) {
283  // add amplitudes
284  adcAmp += (ULong64_t)fPatchADCSimple[i+k][j+l];
285  enAmp += fPatchESimple[i+k][j+l];
286  }
287  }
288 
289  if (adcAmp == 0) {
290  AliDebug(2,"EMCal trigger patch with 0 ADC counts.");
291  continue;
292  }
293 
294  Int_t absId=-1;
295  Int_t cellAbsId[4]={-1,-1,-1,-1};
296 
297  // get low left edge (eta max, phi min)
298  fGeom->GetAbsFastORIndexFromPositionInEMCAL(i, j, absId);
299  // convert to the 4 absId of the cells composing the trigger channel
300  fGeom->GetCellIndexFromFastORIndex(absId, cellAbsId);
301  TVector3 edge1;
302  fGeom->GetGlobal(cellAbsId[0], edge1);
303 
304  // get up right edge (eta min, phi max)
305  fGeom->GetAbsFastORIndexFromPositionInEMCAL(i+patchSize-1, j+patchSize-1, absId);
306  fGeom->GetCellIndexFromFastORIndex(absId, cellAbsId);
307  TVector3 edge2;
308  fGeom->GetGlobal(cellAbsId[3], edge2);
309 
310  // get the center of the patch
311  Int_t offsetCenter = TMath::FloorNint(0.5*patchSize);
312  fGeom->GetAbsFastORIndexFromPositionInEMCAL(i+offsetCenter-1, j+offsetCenter-1, absId);
313  fGeom->GetCellIndexFromFastORIndex(absId, cellAbsId);
314  TVector3 center1;
315  fGeom->GetGlobal(cellAbsId[3], center1);
316 
317  fGeom->GetAbsFastORIndexFromPositionInEMCAL(i+offsetCenter, j+offsetCenter, absId);
318  fGeom->GetCellIndexFromFastORIndex(absId, cellAbsId);
319  TVector3 center2;
320  fGeom->GetGlobal(cellAbsId[0], center2);
321 
322  TVector3 centerGeo(center1);
323  centerGeo += center2;
324  centerGeo *= 0.5;
325 
326  // save the trigger object
327  AliEMCALTriggerPatchInfo *trigger =
328  new ((*fCaloTriggersOut)[itrig]) AliEMCALTriggerPatchInfo();
329  itrig++;
330  trigger->SetTriggerBitConfig(fTriggerBitConfig);
331  trigger->SetCenterGeo(centerGeo, enAmp);
332  trigger->SetEdge1(edge1, enAmp);
333  trigger->SetEdge2(edge2, enAmp);
334  trigger->SetADCAmp(adcAmp);
335  trigger->SetEdgeCell(i*2, j*2); // from triggers to cells
336  }
337  } // trigger algo
338  AliDebug(2,Form("Created %d trigger patches (%d) in this event",itrig,patchSize));
339 
340 }
341 
342 //________________________________________________________________________
344 
345  Int_t dim = TMath::FloorNint((Double_t)(fPatchDim/2.));
346  return dim;
347 }
348 
349 //________________________________________________________________________
351 
352  Int_t dim = GetDimFastor();
353  if(!fL1Slide) return dim;
354 
355  if(dim==2) return 2;
356  else if(dim==4) return 4;
357  else if(dim==8) return 4;
358  else if(dim==16) return 8;
359  else return -1;
360 }
361 
362 
363 
364 
double Double_t
Definition: External.C:58
Definition: External.C:260
Definition: External.C:236
Double_t fPatchESimple[kPatchCols][kPatchRows]
Base task in the EMCAL framework.
Bool_t fLocalInitialized
whether or not the task has been already initialized
AliEMCALTriggerBitConfig * fTriggerBitConfig
int Int_t
Definition: External.C:63
float Float_t
Definition: External.C:68
AliEMCALGeometry * fGeom
!emcal geometry
TString fCaloCellsName
name of calo cell collection
TH2F * fh2CellEnergyVsTime
cell E, eta, phi
short Short_t
Definition: External.C:23
AliVCaloCells * fCaloCells
!cells
Double_t fPatchADCSimple[kPatchCols][kPatchRows]
trigger array out
ClassImp(AliEmcalPatchFromCellMaker) AliEmcalPatchFromCellMaker
AliEmcalList * fOutput
!output list
void SetMakeGeneralHistograms(Bool_t g)
bool Bool_t
Definition: External.C:53
TH1F * fh1CellEnergySum
emcal cell energy vs time