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