1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is 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 **************************************************************************/
16 //------------------------------------------------------------------------//
17 // Fill histograms with basic QA information for EMCAL offline trigger //
18 // Author: Nicolas Arbor (LPSC-Grenoble), Rachid Guernane (LPSC-Grenoble)//
19 // Gustavo Conesa Balbastre (LPSC-Grenoble) //
21 //------------------------------------------------------------------------//
30 #include "AliVCluster.h"
31 #include "AliVCaloCells.h"
32 #include "AliVEvent.h"
33 #include "AliESDEvent.h"
34 #include "AliESDVZERO.h"
35 #include "AliESDCaloTrigger.h"
36 #include "AliEMCALGeometry.h"
37 #include "AliEMCALRecoUtils.h"
39 #include "AliAnalysisTaskEMCALTriggerQA.h"
41 ClassImp(AliAnalysisTaskEMCALTriggerQA)
43 //______________________________________________________________
44 AliAnalysisTaskEMCALTriggerQA::AliAnalysisTaskEMCALTriggerQA() :
47 fGeometry(0), fGeoName("EMCAL_COMPLETEV1"), fRecoUtils(new AliEMCALRecoUtils),
80 fNBinsSTUSignal (2000), fMaxSTUSignal (200000),
81 fNBinsTRUSignal (2000), fMaxTRUSignal (200000),
82 fNBinsV0Signal (2000), fMaxV0Signal (20000),
83 fNBinsSTUFEERatio(2000), fMaxSTUFEERatio(20000),
84 fNBinsSTUTRURatio(2000), fMaxSTUTRURatio(200),
85 fNBinsClusterE (500), fMaxClusterE (100)
92 //______________________________________________________________________________
93 AliAnalysisTaskEMCALTriggerQA::AliAnalysisTaskEMCALTriggerQA(const char *name) :
94 AliAnalysisTaskSE(name),
96 fGeometry(0), fGeoName("EMCAL_COMPLETEV1"), fRecoUtils(new AliEMCALRecoUtils),
129 fNBinsSTUSignal (2000), fMaxSTUSignal (200000),
130 fNBinsTRUSignal (2000), fMaxTRUSignal (200000),
131 fNBinsV0Signal (5000), fMaxV0Signal (50000),
132 fNBinsSTUFEERatio(2000), fMaxSTUFEERatio(20000),
133 fNBinsSTUTRURatio(2000), fMaxSTUTRURatio(200),
134 fNBinsClusterE (500), fMaxClusterE (100)
138 DefineOutput(1, TList::Class());
143 //___________________________________________________________
144 void AliAnalysisTaskEMCALTriggerQA::UserCreateOutputObjects()
146 // Init histograms and geometry
148 fGeometry = AliEMCALGeometry::GetInstance(fGeoName);
150 fOutputList = new TList;
151 fOutputList ->SetOwner(kTRUE);
153 fhNEvents = new TH1F("hNEvents","Number of selected events",7,0,7);
154 fhNEvents ->SetYTitle("N events");
155 fhNEvents ->GetXaxis()->SetBinLabel(1 ,"All");
156 fhNEvents ->GetXaxis()->SetBinLabel(2 ,"INT");
157 fhNEvents ->GetXaxis()->SetBinLabel(3 ,"L0");
158 fhNEvents ->GetXaxis()->SetBinLabel(4 ,"L1-G");
159 fhNEvents ->GetXaxis()->SetBinLabel(5 ,"L1-J");
160 fhNEvents ->GetXaxis()->SetBinLabel(6 ,"L1-G - !L1-J");
161 fhNEvents ->GetXaxis()->SetBinLabel(7 ,"L1-J - !L1-G");
164 fhFORAmp = new TH2F("hFORAmp", "FEE cells deposited energy, grouped like FastOR 2x2 per Row and Column",
165 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
166 fhFORAmp ->SetXTitle("Index #eta (columnns)");
167 fhFORAmp ->SetYTitle("Index #phi (rows)");
168 fhFORAmp ->SetZTitle("Amplitude");
170 fhFORAmpL1G = new TH2F("hFORAmpL1G", "FEE cells deposited energy, grouped like FastOR 2x2 per Row and Column, with L1G trigger condition",
171 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
172 fhFORAmpL1G ->SetXTitle("Index #eta (columnns)");
173 fhFORAmpL1G ->SetYTitle("Index #phi (rows)");
174 fhFORAmpL1G ->SetZTitle("Amplitude");
176 fhFORAmpL1J = new TH2F("hFORAmpL1J", "FEE cells deposited energy, grouped like FastOR 2x2 per Row and Column, with L1J trigger condition",
177 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
178 fhFORAmpL1J ->SetXTitle("Index #eta (columnns)");
179 fhFORAmpL1J ->SetYTitle("Index #phi (rows)");
180 fhFORAmpL1J ->SetZTitle("Amplitude");
182 fhL0Amp = new TH2F("hL0Amp","FALTRO signal per Row and Column",
183 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
184 fhL0Amp ->SetXTitle("Index #eta (columnns)");
185 fhL0Amp ->SetYTitle("Index #phi (rows)");
186 fhL0Amp ->SetZTitle("Amplitude");
188 fhL0AmpL1G = new TH2F("hL0AmpL1G","FALTRO signal per Row and Column, with L1G trigger condition",
189 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
190 fhL0AmpL1G ->SetXTitle("Index #eta (columnns)");
191 fhL0AmpL1G ->SetYTitle("Index #phi (rows)");
192 fhL0AmpL1G ->SetZTitle("Amplitude");
194 fhL0AmpL1J = new TH2F("hL0AmpL1J","FALTRO signal per Row and Column, with L1j trigger condition",
195 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
196 fhL0AmpL1J ->SetXTitle("Index #eta (columnns)");
197 fhL0AmpL1J ->SetYTitle("Index #phi (rows)");
198 fhL0AmpL1J ->SetZTitle("Amplitude");
200 fhL1Amp = new TH2F("hL1Amp","STU signal per Row and Column",
201 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
202 fhL1Amp ->SetXTitle("Index #eta (columnns)");
203 fhL1Amp ->SetYTitle("Index #phi (rows)");
204 fhL1Amp ->SetZTitle("Amplitude");
206 fhL1GAmp = new TH2F("hL1GAmp","STU signal per Row and Column for L1 Gamma",
207 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
208 fhL1GAmp ->SetXTitle("Index #eta (columnns)");
209 fhL1GAmp ->SetYTitle("Index #phi (rows)");
210 fhL1GAmp ->SetZTitle("Amplitude");
212 fhL1JAmp = new TH2F("hL1JAmp","STU signal per Row and Column for L1 Jet",
213 fgkFALTROCols/4,0,fgkFALTROCols,fgkFALTRORows/4,0,fgkFALTRORows);
214 fhL1JAmp ->SetXTitle("Index #eta (columnns)");
215 fhL1JAmp ->SetYTitle("Index #phi (rows)");
216 fhL1JAmp ->SetZTitle("Amplitude");
218 fhL0Patch = new TH2F("hL0Patch","FOR with associated L0 Patch",
219 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
220 fhL0Patch ->SetXTitle("Index #eta (columnns)");
221 fhL0Patch ->SetYTitle("Index #phi (rows)");
222 fhL0Patch ->SetZTitle("counts");
224 fhL1GPatch = new TH2F("hL1GPatch","FOR with associated L1 Gamma Patch",
225 fgkFALTROCols,0,fgkFALTROCols,fgkFALTRORows,0,fgkFALTRORows);
226 fhL1GPatch ->SetXTitle("Index #eta (columnns)");
227 fhL1GPatch ->SetYTitle("Index #phi (rows)");
228 fhL1GPatch ->SetZTitle("counts");
230 fhL1JPatch = new TH2F("hL1JPatch","FOR with associated L1 Jet Patch",
231 fgkFALTROCols/4,0,fgkFALTROCols,fgkFALTRORows/4,0,fgkFALTRORows);
232 fhL1JPatch ->SetXTitle("Index #eta (columnns)");
233 fhL1JPatch ->SetYTitle("Index #phi (rows)");
234 fhL1JPatch ->SetZTitle("counts");
236 fhFullTRUSTU = new TH2I("hFullTRUSTU","Total signal STU vs TRU",
237 fNBinsTRUSignal,0,fMaxTRUSignal,fNBinsSTUSignal,0,fMaxSTUSignal);
238 fhFullTRUSTU->SetXTitle("Total signal TRU");
239 fhFullTRUSTU->SetYTitle("Total signal STU");
240 fhFullTRUSTU->SetZTitle("counts");
242 fhV0STU = new TH2I("hV0STU","Total signal STU vs V0C+V0S",
243 fNBinsV0Signal,0,fMaxV0Signal,fNBinsSTUSignal,0,fMaxSTUSignal);
244 fhV0STU ->SetXTitle("Signal V0C+V0A");
245 fhV0STU ->SetYTitle("Total signal STU");
246 fhV0STU ->SetZTitle("counts");
248 fhSTUChecks = new TH2I("hSTUChecks","Check FEE/STU link",2,0,2,15,0,15);
249 fhSTUChecks ->SetXTitle("Index #eta");
250 fhSTUChecks ->SetYTitle("Index #phi");
252 fhClusMB = new TH1F("hClusMB","clusters distribution for MB trigger",fNBinsClusterE,0,fMaxClusterE);
253 fhClusMB ->SetXTitle("Energy (GeV)");
255 fhClusL0 = new TH1F("hClusL0","clusters distribution for L0 trigger",fNBinsClusterE,0,fMaxClusterE);
256 fhClusMB ->SetXTitle("Energy (GeV)");
258 fhClusL1G = new TH1F("hClusL1G","clusters distribution for L1G trigger",fNBinsClusterE,0,fMaxClusterE);
259 fhClusL1G ->SetXTitle("Energy (GeV)");
261 fhClusL1J = new TH1F("hClusL1J","clusters distribution for L1J trigger",fNBinsClusterE,0,fMaxClusterE);
262 fhClusL1J ->SetXTitle("Energy (GeV)");
264 fhClusL1GOnly = new TH1F("hClusL1GOnly","clusters distribution for L1G trigger and not L1J",fNBinsClusterE,0,fMaxClusterE);
265 fhClusL1GOnly->SetXTitle("Energy (GeV)");
267 fhClusL1JOnly = new TH1F("hClusL1JOnly","clusters distribution for L1J trigger and not L1G",fNBinsClusterE,0,fMaxClusterE);
268 fhClusL1JOnly->SetXTitle("Energy (GeV)");
270 fhClusMaxMB = new TH1F("hClusMaxMB","maximum energy cluster per event for MB trigger",fNBinsClusterE,0,fMaxClusterE);
271 fhClusMaxMB ->SetXTitle("Energy (GeV)");
273 fhClusMaxL0 = new TH1F("hClusMaxL0","maximum energy cluster per event for L0 trigger",fNBinsClusterE,0,fMaxClusterE);
274 fhClusMaxMB ->SetXTitle("Energy (GeV)");
276 fhClusMaxL1G = new TH1F("hClusMaxL1G","maximum energy cluster per event for L1G trigger",fNBinsClusterE,0,fMaxClusterE);
277 fhClusMaxL1G->SetXTitle("Energy (GeV)");
279 fhClusMaxL1J = new TH1F("hClusMaxL1J","maximum energy cluster per event for L1J trigger",fNBinsClusterE,0,fMaxClusterE);
280 fhClusMaxL1J->SetXTitle("Energy (GeV)");
282 fhClusMaxL1GOnly = new TH1F("hClusMaxL1GOnly","maximum energy cluster per event for L1G trigger and not L1J",fNBinsClusterE,0,fMaxClusterE);
283 fhClusMaxL1GOnly->SetXTitle("Energy (GeV)");
285 fhClusMaxL1JOnly = new TH1F("hClusMaxL1JOnly","maximum energy cluster per event for L1J trigger and not L1G",fNBinsClusterE,0,fMaxClusterE);
286 fhClusMaxL1JOnly->SetXTitle("Energy (GeV)");
288 fhFEESTU = new TH2F("hFEESTU","STU / FEE vs channel", fNBinsSTUFEERatio,0,fMaxSTUFEERatio,30,0,30);
289 fhFEESTU ->SetXTitle("STU/FEE signal");
290 fhFEESTU ->SetYTitle("channel");
291 fhFEESTU ->SetZTitle("counts");
293 fhTRUSTU = new TH2F("hTRUSTU","STU / TRU vs channel", fNBinsSTUTRURatio,0,fMaxSTUTRURatio,30,0,30);
294 fhTRUSTU ->SetXTitle("STU/TRU signal");
295 fhTRUSTU ->SetYTitle("channel");
296 fhTRUSTU ->SetZTitle("counts");
298 fhGPMaxVV0TT = new TH2F("hGPMaxVV0TT","Maximum patch of L1-Gamma vs V0 signal",fNBinsV0Signal,0,fMaxV0Signal,1000,0,1000);
299 fhGPMaxVV0TT ->SetXTitle("V0 TT");
300 fhGPMaxVV0TT ->SetYTitle("Patch Max");
302 fhJPMaxVV0TT = new TH2F("hJPMaxVV0TT","Maximum patch of L1-Jet vs V0 signal",fNBinsV0Signal,0,fMaxV0Signal,1000,0,1000);
303 fhJPMaxVV0TT ->SetXTitle("V0 TT");
304 fhJPMaxVV0TT ->SetYTitle("Patch Max");
306 fOutputList->Add(fhNEvents);
307 fOutputList->Add(fhV0STU);
308 fOutputList->Add(fhFORAmp);
309 fOutputList->Add(fhFORAmpL1G);
310 fOutputList->Add(fhFORAmpL1J);
311 fOutputList->Add(fhL0Amp);
312 fOutputList->Add(fhL0AmpL1G);
313 fOutputList->Add(fhL0AmpL1J);
314 fOutputList->Add(fhL1Amp);
315 fOutputList->Add(fhL1GAmp);
316 fOutputList->Add(fhL1JAmp);
317 fOutputList->Add(fhL0Patch);
318 fOutputList->Add(fhL1GPatch);
319 fOutputList->Add(fhL1JPatch);
320 fOutputList->Add(fhFullTRUSTU);
321 fOutputList->Add(fhSTUChecks);
323 fOutputList->Add(fhClusMB);
324 fOutputList->Add(fhClusL0);
325 fOutputList->Add(fhClusL1G);
326 fOutputList->Add(fhClusL1J);
327 fOutputList->Add(fhClusL1GOnly);
328 fOutputList->Add(fhClusL1JOnly);
330 fOutputList->Add(fhClusMaxMB);
331 fOutputList->Add(fhClusMaxL0);
332 fOutputList->Add(fhClusMaxL1G);
333 fOutputList->Add(fhClusMaxL1J);
334 fOutputList->Add(fhClusMaxL1GOnly);
335 fOutputList->Add(fhClusMaxL1JOnly);
337 fOutputList->Add(fhFEESTU);
338 fOutputList->Add(fhTRUSTU);
340 fOutputList->Add(fhGPMaxVV0TT);
341 fOutputList->Add(fhJPMaxVV0TT);
343 PostData(1, fOutputList);
346 //______________________________________________________
347 void AliAnalysisTaskEMCALTriggerQA::UserExec(Option_t *)
351 AliVEvent* event = InputEvent();
353 //Remove next lines when AODs ready
354 AliESDEvent *esdEvent = dynamic_cast<AliESDEvent*>(event);
358 AliError("Work only with ESDs, not available, exit");
362 //trigger configuration
363 TString triggerclasses = esdEvent->GetFiredTriggerClasses();
365 Int_t eventType = ((AliVHeader*)InputEvent()->GetHeader())->GetEventType();
366 // physics events eventType=7, select only those
368 if(triggerclasses=="" || eventType != 7) return;
370 fhNEvents->Fill(0.5);
371 if(triggerclasses.Contains("CINT7-B-NOPF-ALLNOTRD") ||
372 triggerclasses.Contains("CINT7-I-NOPF-ALLNOTRD") ||
373 triggerclasses.Contains("CINT1-I-NOPF-ALLNOTRD") ||
374 triggerclasses.Contains("CINT1-B-NOPF-ALLNOTRD") ||
375 triggerclasses.Contains("CPBI2_B1-B-NOPF-ALLNOTRD") ) fhNEvents->Fill(1.5);
376 if(triggerclasses.Contains("CEMC7-B-NOPF-ALLNOTRD") ||
377 triggerclasses.Contains("CEMC1-B-NOPF-ALLNOTRD") ) fhNEvents->Fill(2.5);
378 if(triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") ||
379 triggerclasses.Contains("CPBI2EGA") ) { fhNEvents->Fill(3.5);
380 if(!triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") &&
381 !triggerclasses.Contains("CPBI2EJE") ) fhNEvents->Fill(5.5); }
382 if(triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") ||
383 triggerclasses.Contains("CPBI2EJE") ) { fhNEvents->Fill(4.5);
384 if(!triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") &&
385 !triggerclasses.Contains("CPBI2EGA") ) fhNEvents->Fill(6.5); }
387 //std::cout << "trigger = " << triggerclasses << std::endl;
389 //map for cells and patches
391 Double_t emcalCell [fgkFALTRORows][fgkFALTROCols], emcalCellL1G [fgkFALTRORows][fgkFALTROCols];
392 Double_t emcalCellL1J [fgkFALTRORows][fgkFALTROCols], emcalTrigL0 [fgkFALTRORows][fgkFALTROCols];
393 Double_t emcalTrigL0L1G[fgkFALTRORows][fgkFALTROCols], emcalTrigL0L1J[fgkFALTRORows][fgkFALTROCols];
394 Double_t emcalTrigL1G [fgkFALTRORows][fgkFALTROCols], emcalTrigL1J [fgkFALTRORows][fgkFALTROCols], emcalTrigL1 [fgkFALTRORows][fgkFALTROCols];
395 Double_t emcalPatchL0 [fgkFALTRORows][fgkFALTROCols], emcalPatchL1G [fgkFALTRORows][fgkFALTROCols], emcalPatchL1J[fgkFALTRORows][fgkFALTROCols];
397 for (Int_t i = 0; i < fgkFALTRORows; i++)
399 for (Int_t j = 0; j < fgkFALTROCols; j++)
401 emcalTrigL0[i][j] = 0.;
402 emcalTrigL0L1G[i][j]= 0.;
403 emcalTrigL0L1J[i][j]= 0.;
404 emcalTrigL1G[i][j] = 0.;
405 emcalTrigL1J[i][j] = 0.;
406 emcalTrigL1[i][j] = 0.;
407 emcalCell[i][j] = 0.;
408 emcalCellL1G[i][j] = 0.;
409 emcalCellL1J[i][j] = 0.;
410 emcalPatchL0[i][j] = 0.;
411 emcalPatchL1G[i][j] = 0.;
412 emcalPatchL1J[i][j] = 0.;
416 // ---------------------------------
418 // Fill FEE energy per channel array
419 // ---------------------------------
421 Int_t posX = -1, posY = -1;
422 Int_t nSupMod = -1, ieta = -1, iphi = -1, nModule = -1, nIphi = -1, nIeta = -1;
426 AliVCaloCells& cells= *(event->GetEMCALCells());
430 for (Int_t icell = 0; icell < cells.GetNumberOfCells(); icell++)
434 Double_t amp =0., time = 0.;
436 cells.GetCell(icell, absId, amp, time);
438 fGeometry->GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
439 fGeometry->GetCellPhiEtaIndexInSModule(nSupMod, nModule, nIphi, nIeta, iphi, ieta);
441 posX = (nSupMod % 2) ? ieta + AliEMCALGeoParams::fgkEMCALCols : ieta;
442 posY = iphi + AliEMCALGeoParams::fgkEMCALRows * int(nSupMod / 2);
444 if(int(posX/2) > fgkFALTROCols || int(posY/2) > fgkFALTRORows ) {
445 if(DebugLevel() > 0) printf("AliAnalysisTaskEMCALTriggerQA::UserExec() - Wrong Position (x,y) = (%d,%d)\n",posX,posY);
449 emcalCell[int(posY/2)][int(posX/2)] += amp;
451 if(triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") || triggerclasses.Contains("CPBI2EGA")) emcalCellL1G[int(posY/2)][int(posX/2)] += amp;
452 if(triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") || triggerclasses.Contains("CPBI2EJE")) emcalCellL1J[int(posY/2)][int(posX/2)] += amp;
457 //-------------------------------------
458 // Trigger analysis, fill L0, L1 arrays
459 //-------------------------------------
461 AliESDCaloTrigger& trg= * (esdEvent->GetCaloTrigger("EMCAL"));
465 Double_t totSTU = 0.;
466 Double_t totTRU = 0.;
471 trg.GetPosition(posX,posY);
474 if (posX > -1 && posY > -1)
478 trg.GetNL0Times(nTimes);
481 trg.GetAmplitude(ampL0);
482 emcalTrigL0[posY][posX] += ampL0;
483 if(triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") || triggerclasses.Contains("CPBI2EGA")) emcalTrigL0L1G[posY][posX] += ampL0;
484 if(triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") || triggerclasses.Contains("CPBI2EJE")) emcalTrigL0L1J[posY][posX] += ampL0;
490 emcalPatchL0[posY][posX] = 1.;
491 fhL0Patch->Fill(posX,posY);
496 trg.GetTriggerBits(bit);
499 trg.GetL1TimeSum(ts);
500 emcalTrigL1 [posY][posX] = ts;
507 emcalPatchL1G[posY][posX] = 1.;
508 fhL1GPatch->Fill(posX,posY);
510 emcalTrigL1G[posY][posX] += ts;
512 //printf("Gamma STU patch %d, time sum %d, posX %d , posY %d\n",nL1Patch,ts,posX, posY);
519 emcalPatchL1J[posY][posX] = 1.;
520 fhL1JPatch->Fill(posX,posY);
522 emcalTrigL1J[posY][posX] += ts;
524 //printf("Jet STU patch %d, time sum %d, posX %d , posY %d\n",nL1Patch,ts,posX, posY);
531 if(totTRU > fMaxTRUSignal && DebugLevel() > 0) printf("AliAnalysisTaskEMCALTriggerQA::UserExec() - Large totTRU %f\n",totTRU);
532 if(totSTU > fMaxSTUSignal && DebugLevel() > 0) printf("AliAnalysisTaskEMCALTriggerQA::UserExec() - Large totSTU %f\n",totSTU);
534 if (totTRU != 0) fhFullTRUSTU->Fill(totTRU,totSTU);
537 AliESDVZERO* eventV0 = esdEvent->GetVZEROData();
539 Float_t v0C = 0, v0A = 0, v0TT = trg.GetL1V0(0)+trg.GetL1V0(1);
543 for (Int_t i = 0; i < 32; i++)
545 v0C += eventV0->GetAdcV0C(i);
546 v0A += eventV0->GetAdcV0A(i);
551 fhV0STU->Fill(v0A+v0C,totSTU);
552 if( v0A+v0C > fMaxV0Signal && DebugLevel() > 0) printf("AliAnalysisTaskEMCALTriggerQA::UserExec() - Large v0A+v0C %f\n",v0A+v0C);
555 //if(nL0Patch!=0 || nL1Patch!=0) printf("total TRU %f, total STU %f, V0C+V0A %f; nL0 %d, nL1 %d \n",
556 // totTRU,totSTU,v0A+v0C,nL0Patch,nL1Patch);
560 for (Int_t i = 0; i < 47; i++)
562 for (Int_t j = 0; j < 59; j++)
566 for (Int_t k = 0; k < 2; k++)
568 for (Int_t l = 0; l < 2; l++)
570 patchG += int(emcalTrigL1[i + k][j + l]);
574 if (patchG > patchMax) patchMax = patchG;
578 fhGPMaxVV0TT->Fill(v0TT, patchMax);
582 for (Int_t i = 0; i < 9; i++)
584 for (Int_t j = 0; j < 12; j++)
588 for (Int_t k = 0; k < 16; k++)
590 for (Int_t l = 0; l < 16; l++)
592 patchJ += int(emcalTrigL1[4*i + k][4*j + l]);
596 if (patchJ > patchMax) patchMax = patchJ;
600 fhJPMaxVV0TT->Fill(v0TT, patchMax);
602 //Matrix with signal per channel
603 for (Int_t i = 0; i < fgkFALTRORows; i++)
605 for (Int_t j = 0; j < fgkFALTROCols; j++) //check x,y direction for reading FOR ((0,0) = top left);
607 fhFORAmp ->Fill( j, i, emcalCell [i][j]);
608 fhFORAmpL1G->Fill( j, i, emcalCellL1G [i][j]);
609 fhFORAmpL1J->Fill( j, i, emcalCellL1J [i][j]);
610 fhL0Amp ->Fill( j, i, emcalTrigL0 [i][j]);
611 fhL0AmpL1G ->Fill( j, i, emcalTrigL0L1G[i][j]);
612 fhL0AmpL1J ->Fill( j, i, emcalTrigL0L1J[i][j]);
613 fhL1Amp ->Fill( j, i, emcalTrigL1 [i][j]);
614 fhL1GAmp ->Fill( j, i, emcalTrigL1G [i][j]);
615 fhL1JAmp ->Fill( j, i, emcalTrigL1J [i][j]);
619 //FEE-TRU-STU correlation checks
620 Double_t ampFOR[30] = {0.}, ampL0[30] = {0.}, ampL1[30] = {0.};
621 for (Int_t i = 0; i < fgkFALTRORows; i++)
623 for (Int_t j = 0; j < fgkFALTROCols; j++)
626 //method to get TRU number
628 fGeometry->GetAbsFastORIndexFromPositionInEMCAL(j,i,idFOR);
631 fGeometry->GetTRUFromAbsFastORIndex(idFOR,iTRU,iADC);
633 //printf("i %d, j %d, iTRU %d, iADC %d, idFOR %d; cell %f, L0 %f, L1 %f\n",
634 // i,j,iTRU,iADC,idFOR, emcalCell [i][j],emcalTrigL0[i][j],emcalTrigL1[i][j]);
638 ampFOR[iTRU] += emcalCell [i][j];
639 ampL0[iTRU] += emcalTrigL0[i][j];
640 ampL1[iTRU] += emcalTrigL1[i][j];
645 // FEE vs STU and TRU vs STU ratios
646 for (Int_t i = 0; i < 30; i++)
649 if (ampFOR[i] != 0 && ampL1[i] != 0) {
650 fhFEESTU->Fill(ampL1[i]/ampFOR[i],i);
651 if(ampL1[i]/ampFOR[i] > fMaxSTUFEERatio && DebugLevel() > 0 ) printf("AliAnalysisTaskEMCALTriggerQA::UserExec() - Large STU/FEE ratio %f\n",ampL1[i]/ampFOR[i]);
654 if (ampL0[i] != 0 && ampL1[i] != 0) {
655 fhTRUSTU->Fill(ampL1[i]/ampL0[i] ,i);
656 if(ampL1[i]/ampL0[i] > fMaxSTUTRURatio && DebugLevel() > 0 ) printf("AliAnalysisTaskEMCALTriggerQA::UserExec() - Large STU/TRU ratio %f\n",ampL1[i]/ampL0[i]);
662 Double_t v[3] = {0,0,0};
663 esdEvent->GetVertex()->GetXYZ(v);
665 //clusters distribution
666 TRefArray* caloClus = new TRefArray();
667 esdEvent->GetEMCALClusters(caloClus);
669 Int_t nCaloClusters = caloClus->GetEntriesFast();
671 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++)
673 AliESDCaloCluster *clus = (AliESDCaloCluster*) (caloClus->At(icalo));
675 if(!clus->IsEMCAL()) continue;
677 if(!fRecoUtils->IsGoodCluster(clus,fGeometry,InputEvent()->GetEMCALCells(),InputEvent()->GetBunchCrossNumber())){
681 if(clus->GetNCells() < 2) continue ; // Avoid 1 cell clusters, noisy.
683 if(clus->E() > emax) emax = clus->E();
685 if(triggerclasses.Contains("CINT7-B-NOPF-ALLNOTRD") ||
686 triggerclasses.Contains("CINT7-I-NOPF-ALLNOTRD") ||
687 triggerclasses.Contains("CINT1-I-NOPF-ALLNOTRD") ||
688 triggerclasses.Contains("CINT1-B-NOPF-ALLNOTRD") ||
689 triggerclasses.Contains("CPBI2_B1-B-NOPF-ALLNOTRD") ) fhClusMB ->Fill(clus->E());
690 if(triggerclasses.Contains("CEMC7-B-NOPF-ALLNOTRD") ||
691 triggerclasses.Contains("CEMC1-B-NOPF-ALLNOTRD") ) fhClusL0 ->Fill(clus->E());
692 if(triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") ||
693 triggerclasses.Contains("CPBI2EGA") ) { fhClusL1G ->Fill(clus->E());
694 if(!triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") &&
695 !triggerclasses.Contains("CPBI2EJE") ) fhClusL1GOnly->Fill(clus->E()); }
696 if(triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") ||
697 triggerclasses.Contains("CPBI2EJE") ) { fhClusL1J ->Fill(clus->E());
698 if(!triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") &&
699 !triggerclasses.Contains("CPBI2EGA") ) fhClusL1JOnly->Fill(clus->E()); }
703 if(triggerclasses.Contains("CINT7-B-NOPF-ALLNOTRD") ||
704 triggerclasses.Contains("CINT7-I-NOPF-ALLNOTRD") ||
705 triggerclasses.Contains("CINT1-I-NOPF-ALLNOTRD") ||
706 triggerclasses.Contains("CINT1-B-NOPF-ALLNOTRD") ||
707 triggerclasses.Contains("CPBI2_B1-B-NOPF-ALLNOTRD") ) fhClusMaxMB ->Fill(emax);
708 if(triggerclasses.Contains("CEMC7-B-NOPF-ALLNOTRD") ||
709 triggerclasses.Contains("CEMC1-B-NOPF-ALLNOTRD") ) fhClusMaxL0 ->Fill(emax);
710 if(triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") ||
711 triggerclasses.Contains("CPBI2EGA") ) { fhClusMaxL1G ->Fill(emax);
712 if(!triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") &&
713 !triggerclasses.Contains("CPBI2EJE") ) fhClusMaxL1GOnly->Fill(emax); }
714 if(triggerclasses.Contains("CEMC7EJE-B-NOPF-CENTNOTRD") ||
715 triggerclasses.Contains("CPBI2EJE") ) { fhClusMaxL1J ->Fill(emax);
716 if(!triggerclasses.Contains("CEMC7EGA-B-NOPF-CENTNOTRD") &&
717 !triggerclasses.Contains("CPBI2EGA") ) fhClusMaxL1JOnly->Fill(emax); }
719 PostData(1, fOutputList);
723 //_______________________________________________________
724 void AliAnalysisTaskEMCALTriggerQA::Terminate(Option_t *)
726 // Terminate analysis
730 for (Int_t i = 0; i < 30; i++)
732 checkSTU[i] = 1 ; // Init array
733 TH1F* hTRUSTUChannel = (TH1F*) fhTRUSTU->ProjectionY(Form("hTRUSTUChannel%d",i),i,i);
734 //printf("AliAnalysisTaskEMCALTriggerQA::Terminate() - Channel %d TRUSTU Entries %d, Integral(10,20) %f ",
735 // i, hTRUSTUChannel->GetEntries(), hTRUSTUChannel->Integral(10,20));
736 Int_t binMin = hTRUSTUChannel->FindBin(10);
737 Int_t binMax = hTRUSTUChannel->FindBin(20);
738 if (hTRUSTUChannel->GetEntries() > 0 &&
739 hTRUSTUChannel->Integral(binMin,binMax)/hTRUSTUChannel->GetEntries() < 0.9)
741 else if(hTRUSTUChannel->GetEntries() <= 0)
744 //printf(" - %d\n",checkSTU[i]);
745 delete hTRUSTUChannel;
748 for (Int_t i = 0; i < 30; i++)
750 if (i<15) fhSTUChecks->Fill(0.,i,checkSTU[i]);
751 else fhSTUChecks->Fill(1.,i,checkSTU[i]);