-add aod debug level to avoid huge logfiles when testing lego trains
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliAnalysisTaskHFECal.cxx
CommitLineData
bfc7c23b 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 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 **************************************************************************/
15
16// Class for heavy-flavour electron with EMCal triggered events
17// Author: Shingo Sakai
18
19
20#include "TChain.h"
21#include "TTree.h"
22#include "TH2F.h"
23#include "TMath.h"
24#include "TCanvas.h"
25#include "THnSparse.h"
26#include "TLorentzVector.h"
27#include "TString.h"
28#include "TFile.h"
29
30#include "AliAnalysisTask.h"
31#include "AliAnalysisManager.h"
32
33#include "AliESDEvent.h"
34#include "AliESDHandler.h"
35#include "AliAODEvent.h"
36#include "AliAODHandler.h"
37
e7d87aef 38#include "AliMCEventHandler.h"
39#include "AliMCEvent.h"
40#include "AliMCParticle.h"
41
bfc7c23b 42#include "AliAnalysisTaskHFECal.h"
43#include "TGeoGlobalMagField.h"
44#include "AliLog.h"
45#include "AliAnalysisTaskSE.h"
46#include "TRefArray.h"
47#include "TVector.h"
48#include "AliESDInputHandler.h"
49#include "AliESDpid.h"
50#include "AliESDtrackCuts.h"
51#include "AliPhysicsSelection.h"
52#include "AliESDCaloCluster.h"
53#include "AliAODCaloCluster.h"
54#include "AliEMCALRecoUtils.h"
55#include "AliEMCALGeometry.h"
56#include "AliGeomManager.h"
57#include "stdio.h"
58#include "TGeoManager.h"
59#include "iostream"
60#include "fstream"
61
62#include "AliEMCALTrack.h"
63#include "AliMagF.h"
64
65#include "AliKFParticle.h"
66#include "AliKFVertex.h"
67
68#include "AliPID.h"
69#include "AliPIDResponse.h"
70#include "AliHFEcontainer.h"
71#include "AliHFEcuts.h"
72#include "AliHFEpid.h"
73#include "AliHFEpidBase.h"
74#include "AliHFEpidQAmanager.h"
75#include "AliHFEtools.h"
76#include "AliCFContainer.h"
77#include "AliCFManager.h"
78
e7d87aef 79#include "AliStack.h"
80
bfc7c23b 81#include "AliCentrality.h"
82
83ClassImp(AliAnalysisTaskHFECal)
84//________________________________________________________________________
85AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name)
86 : AliAnalysisTaskSE(name)
87 ,fESD(0)
e7d87aef 88 ,fMC(0)
fb87d707 89 ,fGeom(0)
bfc7c23b 90 ,fOutputList(0)
91 ,fTrackCuts(0)
92 ,fCuts(0)
93 ,fIdentifiedAsOutInz(kFALSE)
94 ,fPassTheEventCut(kFALSE)
95 ,fRejectKinkMother(kFALSE)
e7d87aef 96 ,fmcData(kFALSE)
bfc7c23b 97 ,fVz(0.0)
98 ,fCFM(0)
99 ,fPID(0)
100 ,fPIDqa(0)
101 ,fOpeningAngleCut(0.1)
102 ,fInvmassCut(0.01)
103 ,fNoEvents(0)
104 ,fEMCAccE(0)
105 ,fTrkpt(0)
106 ,fTrkEovPBef(0)
107 ,fTrkEovPAft(0)
108 ,fdEdxBef(0)
109 ,fdEdxAft(0)
110 ,fIncpT(0)
fb87d707 111 ,fIncpTM20(0)
bfc7c23b 112 ,fInvmassLS(0)
113 ,fInvmassULS(0)
114 ,fOpeningAngleLS(0)
115 ,fOpeningAngleULS(0)
116 ,fPhotoElecPt(0)
117 ,fPhoElecPt(0)
fb87d707 118 ,fPhoElecPtM20(0)
f09766dd 119 ,fSameElecPt(0)
fb87d707 120 ,fSameElecPtM20(0)
bfc7c23b 121 ,fTrackPtBefTrkCuts(0)
122 ,fTrackPtAftTrkCuts(0)
123 ,fTPCnsigma(0)
124 ,fCent(0)
125 ,fEleInfo(0)
fb87d707 126 ,fClsEBftTrigCut(0)
127 ,fClsEAftTrigCut(0)
128 ,fClsEAftTrigCut1(0)
129 ,fClsEAftTrigCut2(0)
130 ,fClsEAftTrigCut3(0)
131 ,fClsEAftTrigCut4(0)
132 ,fClsETime(0)
133 ,fClsETime1(0)
134 ,fTrigTimes(0)
42c75692 135 ,fCellCheck(0)
e7d87aef 136 ,fInputHFEMC(0)
137 ,fIncpTMChfe(0)
138 ,fIncpTMCM20hfe(0)
139 ,fIncpTMCpho(0)
140 ,fIncpTMCM20pho(0)
141 ,fPhoElecPtMC(0)
142 ,fPhoElecPtMCM20(0)
143 ,fSameElecPtMC(0)
144 ,fSameElecPtMCM20(0)
bfc7c23b 145{
146 //Named constructor
147
148 fPID = new AliHFEpid("hfePid");
149 fTrackCuts = new AliESDtrackCuts();
150
151 // Define input and output slots here
152 // Input slot #0 works with a TChain
153 DefineInput(0, TChain::Class());
154 // Output slot #0 id reserved by the base class for AOD
155 // Output slot #1 writes into a TH1 container
156 // DefineOutput(1, TH1I::Class());
157 DefineOutput(1, TList::Class());
158 // DefineOutput(3, TTree::Class());
159}
160
161//________________________________________________________________________
162AliAnalysisTaskHFECal::AliAnalysisTaskHFECal()
163 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskHFECal")
164 ,fESD(0)
e7d87aef 165 ,fMC(0)
fb87d707 166 ,fGeom(0)
bfc7c23b 167 ,fOutputList(0)
168 ,fTrackCuts(0)
169 ,fCuts(0)
170 ,fIdentifiedAsOutInz(kFALSE)
171 ,fPassTheEventCut(kFALSE)
172 ,fRejectKinkMother(kFALSE)
e7d87aef 173 ,fmcData(kFALSE)
bfc7c23b 174 ,fVz(0.0)
175 ,fCFM(0)
176 ,fPID(0)
177 ,fPIDqa(0)
178 ,fOpeningAngleCut(0.1)
179 ,fInvmassCut(0.01)
180 ,fNoEvents(0)
181 ,fEMCAccE(0)
182 ,fTrkpt(0)
183 ,fTrkEovPBef(0)
184 ,fTrkEovPAft(0)
185 ,fdEdxBef(0)
186 ,fdEdxAft(0)
187 ,fIncpT(0)
fb87d707 188 ,fIncpTM20(0)
bfc7c23b 189 ,fInvmassLS(0)
190 ,fInvmassULS(0)
191 ,fOpeningAngleLS(0)
192 ,fOpeningAngleULS(0)
193 ,fPhotoElecPt(0)
194 ,fPhoElecPt(0)
fb87d707 195 ,fPhoElecPtM20(0)
f09766dd 196 ,fSameElecPt(0)
fb87d707 197 ,fSameElecPtM20(0)
bfc7c23b 198 ,fTrackPtBefTrkCuts(0)
199 ,fTrackPtAftTrkCuts(0)
200 ,fTPCnsigma(0)
201 ,fCent(0)
202 ,fEleInfo(0)
fb87d707 203 ,fClsEBftTrigCut(0)
204 ,fClsEAftTrigCut(0)
205 ,fClsEAftTrigCut1(0)
206 ,fClsEAftTrigCut2(0)
207 ,fClsEAftTrigCut3(0)
208 ,fClsEAftTrigCut4(0)
209 ,fClsETime(0)
210 ,fClsETime1(0)
211 ,fTrigTimes(0)
42c75692 212 ,fCellCheck(0)
e7d87aef 213 ,fInputHFEMC(0)
214 ,fIncpTMChfe(0)
215 ,fIncpTMCM20hfe(0)
216 ,fIncpTMCpho(0)
217 ,fIncpTMCM20pho(0)
218 ,fPhoElecPtMC(0)
219 ,fPhoElecPtMCM20(0)
220 ,fSameElecPtMC(0)
221 ,fSameElecPtMCM20(0)
bfc7c23b 222{
223 //Default constructor
224 fPID = new AliHFEpid("hfePid");
225
226 fTrackCuts = new AliESDtrackCuts();
227
228 // Constructor
229 // Define input and output slots here
230 // Input slot #0 works with a TChain
231 DefineInput(0, TChain::Class());
232 // Output slot #0 id reserved by the base class for AOD
233 // Output slot #1 writes into a TH1 container
234 // DefineOutput(1, TH1I::Class());
235 DefineOutput(1, TList::Class());
236 //DefineOutput(3, TTree::Class());
237}
238//_________________________________________
239
240AliAnalysisTaskHFECal::~AliAnalysisTaskHFECal()
241{
242 //Destructor
243
244 delete fOutputList;
fb87d707 245 delete fGeom;
bfc7c23b 246 delete fPID;
247 delete fCFM;
248 delete fPIDqa;
249 delete fTrackCuts;
250}
251//_________________________________________
252
253void AliAnalysisTaskHFECal::UserExec(Option_t*)
254{
255 //Main loop
256 //Called for each event
257
258 // create pointer to event
259 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
260 if (!fESD) {
261 printf("ERROR: fESD not available\n");
262 return;
263 }
264
265 if(!fCuts){
266 AliError("HFE cuts not available");
267 return;
268 }
269
270 if(!fPID->IsInitialized()){
271 // Initialize PID with the given run number
272 AliWarning("PID not initialised, get from Run no");
273 fPID->InitializePID(fESD->GetRunNumber());
274 }
275
e7d87aef 276 if(fmcData)fMC = MCEvent();
277 AliStack* stack = NULL;
278 if(fmcData && fMC)stack = fMC->Stack();
279
280 Float_t cent = -1.;
281 AliCentrality *centrality = fESD->GetCentrality();
282 cent = centrality->GetCentralityPercentile("V0M");
283
284 //---- fill MC track info
285 if(fmcData && fMC)
286 {
287 Int_t nParticles = stack->GetNtrack();
288 for (Int_t iParticle = 0; iParticle < nParticles; iParticle++) {
289 TParticle* particle = stack->Particle(iParticle);
290 int fPDG = particle->GetPdgCode();
291 double mcZvertex = fMC->GetPrimaryVertex()->GetZ();
292 double pTMC = particle->Pt();
293 Bool_t mcInDtoE= kFALSE;
294 Bool_t mcInBtoE= kFALSE;
295
296 if(particle->GetFirstMother()>-1 && fabs(fPDG)==11)
297 {
298 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
299 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(fPDG)==11)mcInDtoE = kTRUE;
300 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(fPDG)==11)mcInBtoE = kTRUE;
301 if((mcInBtoE || mcInDtoE) && fabs(mcZvertex)<10.0)fInputHFEMC->Fill(cent,pTMC);
302 }
303
304 }
305 }
306
bfc7c23b 307 fNoEvents->Fill(0);
308
309 Int_t fNOtrks = fESD->GetNumberOfTracks();
310 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
311
312 Double_t pVtxZ = -999;
313 pVtxZ = pVtx->GetZ();
314
315 if(TMath::Abs(pVtxZ)>10) return;
316 fNoEvents->Fill(1);
317
318 if(fNOtrks<2) return;
319 fNoEvents->Fill(2);
320
321 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
322 if(!pidResponse){
323 AliDebug(1, "Using default PID Response");
324 pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
325 }
326
327 fPID->SetPIDResponse(pidResponse);
328
329 fCFM->SetRecEventInfo(fESD);
330
e7d87aef 331 //Float_t cent = -1.;
332 //AliCentrality *centrality = fESD->GetCentrality();
333 //cent = centrality->GetCentralityPercentile("V0M");
bfc7c23b 334 fCent->Fill(cent);
335
336 if(cent>90.) return;
337
338 // Calorimeter info.
339
42c75692 340 FindTriggerClusters();
341
bfc7c23b 342 // make EMCAL array
343 for(Int_t iCluster=0; iCluster<fESD->GetNumberOfCaloClusters(); iCluster++)
344 {
345 AliESDCaloCluster *clust = fESD->GetCaloCluster(iCluster);
346 if(clust->IsEMCAL())
347 {
348 double clustE = clust->E();
349 float emcx[3]; // cluster pos
350 clust->GetPosition(emcx);
351 TVector3 clustpos(emcx[0],emcx[1],emcx[2]);
352 double emcphi = clustpos.Phi();
353 double emceta = clustpos.Eta();
fb87d707 354 double calInfo[5];
355 calInfo[0] = emcphi; calInfo[1] = emceta; calInfo[2] = clustE; calInfo[3] = cent; calInfo[4] = clust->Chi2();
f09766dd 356 if(clustE>1.5)fEMCAccE->Fill(calInfo);
bfc7c23b 357 }
358 }
359
360 // Track loop
361 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
362 AliESDtrack* track = fESD->GetTrack(iTracks);
363 if (!track) {
364 printf("ERROR: Could not receive track %d\n", iTracks);
365 continue;
366 }
e7d87aef 367
368 Bool_t mcPho = kFALSE;
369 Bool_t mcDtoE= kFALSE;
370 Bool_t mcBtoE= kFALSE;
371 double mcele = -1.;
372 if(fmcData && fMC && stack)
373 {
374 Int_t label = track->GetLabel();
375 TParticle* particle = stack->Particle(label);
376 int mcpid = particle->GetPdgCode();
377 printf("MCpid = %d",mcpid);
378 if(particle->GetFirstMother()>-1)
379 {
380 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
381 if((parentPID==22 || parentPID==111 || parentPID==221)&& fabs(mcpid)==11)mcPho = kTRUE;
382 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(mcpid)==11)mcDtoE = kTRUE;
383 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(mcpid)==11)mcBtoE = kTRUE;
384 }
385 if(fabs(mcpid)==11)mcele= 1.;
386 }
387
bfc7c23b 388 if(TMath::Abs(track->Eta())>0.7) continue;
f09766dd 389 if(TMath::Abs(track->Pt()<2.0)) continue;
bfc7c23b 390
391 fTrackPtBefTrkCuts->Fill(track->Pt());
392 // RecKine: ITSTPC cuts
393 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
394
395 //RecKink
396 if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
397 if(track->GetKinkIndex(0) != 0) continue;
398 }
399
400 // RecPrim
401 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
402
403 // HFEcuts: ITS layers cuts
404 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
405
f09766dd 406 // HFE cuts: TPC PID cleanup
407 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
408
bfc7c23b 409
410 fTrackPtAftTrkCuts->Fill(track->Pt());
411
412 Double_t mom = -999., eop=-999., pt = -999., dEdx=-999., fTPCnSigma=-10, phi=-999., eta=-999.;
fb87d707 413 //pt = track->Pt();
414 //if(pt<2.0)continue;
bfc7c23b 415
416 // Track extrapolation
417
bfc7c23b 418 Int_t charge = track->Charge();
419 fTrkpt->Fill(pt);
420 mom = track->P();
421 phi = track->Phi();
422 eta = track->Eta();
423 dEdx = track->GetTPCsignal();
424 fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
425
426 double ncells = -1.0;
427 double m20 = -1.0;
428 double m02 = -1.0;
429 double disp = -1.0;
430 double rmatch = -1.0;
431 double nmatch = -1.0;
f09766dd 432 double oppstatus = 0.0;
433 double samestatus = 0.0;
bfc7c23b 434
f09766dd 435 Bool_t fFlagPhotonicElec = kFALSE;
436 Bool_t fFlagConvinatElec = kFALSE;
437 SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicElec,fFlagConvinatElec);
438 if(fFlagPhotonicElec)oppstatus = 1.0;
439 if(fFlagConvinatElec)samestatus = 1.0;
bfc7c23b 440
441 Int_t clsId = track->GetEMCALcluster();
442 if (clsId>0){
443 AliESDCaloCluster *clust = fESD->GetCaloCluster(clsId);
444 if(clust && clust->IsEMCAL()){
445
446 double clustE = clust->E();
447 eop = clustE/fabs(mom);
448 //double clustT = clust->GetTOF();
449 ncells = clust->GetNCells();
450 m02 = clust->GetM02();
451 m20 = clust->GetM20();
452 disp = clust->GetDispersion();
453 double delphi = clust->GetTrackDx();
454 double deleta = clust->GetTrackDz();
455 rmatch = sqrt(pow(delphi,2)+pow(deleta,2));
456 nmatch = clust->GetNTracksMatched();
457
e7d87aef 458 double valdedx[17];
fb87d707 459 valdedx[0] = pt; valdedx[1] = dEdx; valdedx[2] = phi; valdedx[3] = eta; valdedx[4] = fTPCnSigma;
460 valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = ncells, valdedx[8] = m02; valdedx[9] = m20; valdedx[10] = disp;
461 valdedx[11] = cent; valdedx[12] = charge; valdedx[13] = oppstatus; valdedx[14] = samestatus; valdedx[15] = clust->Chi2();
e7d87aef 462 valdedx[16] = mcele;
bfc7c23b 463 fEleInfo->Fill(valdedx);
464
465
466 }
467 }
468
42c75692 469 fdEdxBef->Fill(mom,fTPCnSigma);
bfc7c23b 470 fTPCnsigma->Fill(mom,fTPCnSigma);
f09766dd 471 if(fTPCnSigma >= -1.0 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,eop);
bfc7c23b 472
bfc7c23b 473 Int_t pidpassed = 1;
474
475
476 //--- track accepted
477 AliHFEpidObject hfetrack;
478 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
479 hfetrack.SetRecTrack(track);
fb87d707 480 int centf = (int)cent;
481 hfetrack.SetCentrality(centf); //added
bfc7c23b 482 hfetrack.SetPbPb();
483 if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
484
485 if(pidpassed==0) continue;
486
487 fTrkEovPAft->Fill(pt,eop);
42c75692 488 fdEdxAft->Fill(mom,fTPCnSigma);
bfc7c23b 489
e7d87aef 490 // Fill real data
fb87d707 491 fIncpT->Fill(cent,pt);
bfc7c23b 492 if(fFlagPhotonicElec) fPhoElecPt->Fill(cent,pt);
f09766dd 493 if(fFlagConvinatElec) fSameElecPt->Fill(cent,pt);
fb87d707 494
495 if(m20>0.0 && m20<0.3)
496 {
497 fIncpTM20->Fill(cent,pt);
498 if(fFlagPhotonicElec) fPhoElecPtM20->Fill(cent,pt);
499 if(fFlagConvinatElec) fSameElecPtM20->Fill(cent,pt);
500 }
e7d87aef 501
502 // MC
503 if(mcele==1)
504 {
505 if(mcBtoE || mcDtoE)
506 {
507 fIncpTMChfe->Fill(cent,pt);
508 if(m20>0.0 && m20<0.3)fIncpTMCM20hfe->Fill(cent,pt);
509 }
510 if(mcPho)
511 {
512 fIncpTMCpho->Fill(cent,pt);
513 if(fFlagPhotonicElec) fPhoElecPtMC->Fill(cent,pt);
514 if(fFlagConvinatElec) fSameElecPtMC->Fill(cent,pt);
515
516 if(m20>0.0 && m20<0.3)
517 {
518 fIncpTMCM20pho->Fill(cent,pt);
519 if(fFlagPhotonicElec) fPhoElecPtMCM20->Fill(cent,pt);
520 if(fFlagConvinatElec) fSameElecPtMCM20->Fill(cent,pt);
521 }
522 }
523 }
bfc7c23b 524 }
525 PostData(1, fOutputList);
526}
527//_________________________________________
528void AliAnalysisTaskHFECal::UserCreateOutputObjects()
529{
42c75692 530 //--- Check MC
531
e7d87aef 532 //Bool_t mcData = kFALSE;
42c75692 533 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
534 {
e7d87aef 535 fmcData = kTRUE;
42c75692 536 printf("+++++ MC Data available");
537 }
e7d87aef 538 if(fmcData)
42c75692 539 {
540 printf("++++++++= MC analysis \n");
541 }
542 else
543 {
544 printf("++++++++ real data analysis \n");
545 }
546
fb87d707 547 //---- Geometry
548 fGeom = AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
549
bfc7c23b 550 //--------Initialize PID
42c75692 551 //fPID->SetHasMCData(kFALSE);
e7d87aef 552 fPID->SetHasMCData(fmcData);
bfc7c23b 553 if(!fPID->GetNumberOfPIDdetectors())
554 {
555 fPID->AddDetector("TPC", 0);
556 fPID->AddDetector("EMCAL", 1);
557 }
558
559 fPID->SortDetectors();
560 fPIDqa = new AliHFEpidQAmanager();
561 fPIDqa->Initialize(fPID);
562
563 //--------Initialize correction Framework and Cuts
564 fCFM = new AliCFManager;
565 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
566 fCFM->SetNStepParticle(kNcutSteps);
567 for(Int_t istep = 0; istep < kNcutSteps; istep++)
568 fCFM->SetParticleCutsList(istep, NULL);
569
570 if(!fCuts){
571 AliWarning("Cuts not available. Default cuts will be used");
572 fCuts = new AliHFEcuts;
573 fCuts->CreateStandardCuts();
574 }
575 fCuts->Initialize(fCFM);
576
577 //---------Output Tlist
578 fOutputList = new TList();
579 fOutputList->SetOwner();
580 fOutputList->Add(fPIDqa->MakeList("PIDQA"));
581
582 fNoEvents = new TH1F("fNoEvents","",4,-0.5,3.5) ;
583 fOutputList->Add(fNoEvents);
584
fb87d707 585 Int_t binsE[5] = {250, 100, 1000, 100, 10};
586 Double_t xminE[5] = {1.0, -1, 0.0, 0, -0.5};
587 Double_t xmaxE[5] = {3.5, 1, 100.0, 100, 9.5};
588 fEMCAccE = new THnSparseD("fEMCAccE","EMC acceptance & E;#eta;#phi;Energy;Centrality;trugCondition;",5,binsE,xminE,xmaxE);
bfc7c23b 589 fOutputList->Add(fEMCAccE);
590
591 fTrkpt = new TH1F("fTrkpt","track pt",100,0,50);
592 fOutputList->Add(fTrkpt);
593
594 fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50);
595 fOutputList->Add(fTrackPtBefTrkCuts);
596
597 fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50);
598 fOutputList->Add(fTrackPtAftTrkCuts);
599
600 fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10);
601 fOutputList->Add(fTPCnsigma);
602
603 fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2);
604 fOutputList->Add(fTrkEovPBef);
605
606 fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2);
607 fOutputList->Add(fTrkEovPAft);
608
42c75692 609 fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,200,-10,10);
bfc7c23b 610 fOutputList->Add(fdEdxBef);
611
42c75692 612 fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,200,-10,10);
bfc7c23b 613 fOutputList->Add(fdEdxAft);
614
615 fIncpT = new TH2F("fIncpT","HFE pid electro vs. centrality",100,0,100,100,0,50);
616 fOutputList->Add(fIncpT);
617
fb87d707 618 fIncpTM20 = new TH2F("fIncpTM20","HFE pid electro vs. centrality with M20",100,0,100,100,0,50);
619 fOutputList->Add(fIncpTM20);
f09766dd 620
621 Int_t nBinspho[3] = { 100, 100, 500};
622 Double_t minpho[3] = { 0., 0., 0.};
623 Double_t maxpho[3] = {100., 50., 0.5};
624
625 fInvmassLS = new THnSparseD("fInvmassLS", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2);", 3, nBinspho,minpho, maxpho);
bfc7c23b 626 fOutputList->Add(fInvmassLS);
627
f09766dd 628 fInvmassULS = new THnSparseD("fInvmassULS", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2);", 3, nBinspho,minpho, maxpho);
bfc7c23b 629 fOutputList->Add(fInvmassULS);
630
631 fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
632 fOutputList->Add(fOpeningAngleLS);
633
634 fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
635 fOutputList->Add(fOpeningAngleULS);
636
637 fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50);
638 fOutputList->Add(fPhotoElecPt);
639
f09766dd 640 fPhoElecPt = new TH2F("fPhoElecPt", "Pho-inclusive electron pt",100,0,100,100,0,50);
bfc7c23b 641 fOutputList->Add(fPhoElecPt);
642
fb87d707 643 fPhoElecPtM20 = new TH2F("fPhoElecPtM20", "Pho-inclusive electron pt with M20",100,0,100,100,0,50);
644 fOutputList->Add(fPhoElecPtM20);
645
f09766dd 646 fSameElecPt = new TH2F("fSameElecPt", "Same-inclusive electron pt",100,0,100,100,0,50);
647 fOutputList->Add(fSameElecPt);
648
fb87d707 649 fSameElecPtM20 = new TH2F("fSameElecPtM20", "Same-inclusive electron pt with M20",100,0,100,100,0,50);
650 fOutputList->Add(fSameElecPtM20);
651
bfc7c23b 652 fCent = new TH1F("fCent","Centrality",100,0,100) ;
653 fOutputList->Add(fCent);
654
655 // Make common binning
f09766dd 656 const Double_t kMinP = 2.;
bfc7c23b 657 const Double_t kMaxP = 50.;
658 const Double_t kTPCSigMim = 40.;
659 const Double_t kTPCSigMax = 140.;
660
661 // 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select)
e7d87aef 662 Int_t nBins[17] = { 480, 200, 60, 20, 600, 300, 100, 40, 200, 200, 200, 100, 3, 3, 3, 10, 3};
663 Double_t min[17] = {kMinP, kTPCSigMim, 1.0, -1.0, -8.0, 0, 0, 0, 0.0, 0.0, 0.0, 0, -1.5, -0.5, -0.5, -0.5, -1.5};
664 Double_t max[17] = {kMaxP, kTPCSigMax, 4.0, 1.0, 4.0, 3.0, 0.1, 40, 2.0, 2.0, 2.0, 100, 1.5, 2.5, 2.5, 9.5, 1.5};
665 fEleInfo = new THnSparseD("fEleInfo", "Electron Info; pT [GeV/c]; TPC signal;phi;eta;nSig; E/p;Rmatch;Ncell;M02;M20;Disp;Centrality;charge;opp;same;trigCond;", 17, nBins, min, max);
bfc7c23b 666 fOutputList->Add(fEleInfo);
667
fb87d707 668 //<--- Trigger info
669
670 fClsEBftTrigCut = new TH1F("fClsEBftTrigCut","cluster E before trigger selection",1000,0,100);
671 fOutputList->Add(fClsEBftTrigCut);
672
673 fClsEAftTrigCut = new TH1F("fClsEAftTrigCut","cluster E if cls has 0 trigcut channel",1000,0,100);
674 fOutputList->Add(fClsEAftTrigCut);
675
676 fClsEAftTrigCut1 = new TH1F("fClsEAftTrigCut1","cluster E if cls with trig channel",1000,0,100);
677 fOutputList->Add(fClsEAftTrigCut1);
678
679 fClsEAftTrigCut2 = new TH1F("fClsEAftTrigCut2","cluster E if cls with trigcut channel",1000,0,100);
680 fOutputList->Add(fClsEAftTrigCut2);
681
682 fClsEAftTrigCut3 = new TH1F("fClsEAftTrigCut3","cluster E if cls with trigcut channel + nCell>Ecorrect",1000,0,100);
683 fOutputList->Add(fClsEAftTrigCut3);
684
685 fClsEAftTrigCut4 = new TH1F("fClsEAftTrigCut4","cluster E if cls with trigcut channel + nCell>Ecorrect + cls time cut",1000,0,100);
686 fOutputList->Add(fClsEAftTrigCut4);
687
688 fClsETime = new TH2F("fClsETime", "Cls time vs E; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
689 fOutputList->Add(fClsETime);
690
691 fClsETime1 = new TH2F("fClsETime1", "Cls time vs E if cls contains trigger channel; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
692 fOutputList->Add(fClsETime1);
693
694 fTrigTimes = new TH1F("fTrigTimes", "Trigger time; time; N;",25,0,25);
695 fOutputList->Add(fTrigTimes);
696
42c75692 697 fCellCheck = new TH2F("fCellCheck", "Cell vs E; E GeV; Cell ID",10,6,26,12000,0,12000);
698 fOutputList->Add(fCellCheck);
699
e7d87aef 700 //<---------- MC
701
702 fInputHFEMC = new TH2F("fInputHFEMC","Input MC HFE pid electro vs. centrality",100,0,100,100,0,50);
703 fOutputList->Add(fInputHFEMC);
704
705 fIncpTMChfe = new TH2F("fIncpTMChfe","MC HFE pid electro vs. centrality",100,0,100,100,0,50);
706 fOutputList->Add(fIncpTMChfe);
707
708 fIncpTMCM20hfe = new TH2F("fIncpTMCM20hfe","MC HFE pid electro vs. centrality with M20",100,0,100,100,0,50);
709 fOutputList->Add(fIncpTMCM20hfe);
710
711 fIncpTMCpho = new TH2F("fIncpTMCpho","MC Pho pid electro vs. centrality",100,0,100,100,0,50);
712 fOutputList->Add(fIncpTMCpho);
713
714 fIncpTMCM20pho = new TH2F("fIncpTMCM20pho","MC Pho pid electro vs. centrality with M20",100,0,100,100,0,50);
715 fOutputList->Add(fIncpTMCM20pho);
716
717 fPhoElecPtMC = new TH2F("fPhoElecPtMC", "MC Pho-inclusive electron pt",100,0,100,100,0,50);
718 fOutputList->Add(fPhoElecPtMC);
719
720 fPhoElecPtMCM20 = new TH2F("fPhoElecPtMCM20", "MC Pho-inclusive electron pt with M20",100,0,100,100,0,50);
721 fOutputList->Add(fPhoElecPtM20);
722
723 fSameElecPtMC = new TH2F("fSameElecPtMC", "MC Same-inclusive electron pt",100,0,100,100,0,50);
724 fOutputList->Add(fSameElecPtMC);
725
726 fSameElecPtMCM20 = new TH2F("fSameElecPtMCM20", "MC Same-inclusive electron pt with M20",100,0,100,100,0,50);
727 fOutputList->Add(fSameElecPtMCM20);
728
729
bfc7c23b 730 PostData(1,fOutputList);
731}
732
733//________________________________________________________________________
734void AliAnalysisTaskHFECal::Terminate(Option_t *)
735{
736 // Info("Terminate");
737 AliAnalysisTaskSE::Terminate();
738}
739
740//________________________________________________________________________
741Bool_t AliAnalysisTaskHFECal::ProcessCutStep(Int_t cutStep, AliVParticle *track)
742{
743 // Check single track cuts for a given cut step
744 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
745 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
746 return kTRUE;
747}
748//_________________________________________
f09766dd 749//void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec)
750void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec, Bool_t &fFlagConvinatElec)
bfc7c23b 751{
752 //Identify non-heavy flavour electrons using Invariant mass method
753
fb87d707 754 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
f09766dd 755 //fTrackCuts->SetRequireTPCRefit(kTRUE);
756 //fTrackCuts->SetEtaRange(-0.7,0.7);
757 //fTrackCuts->SetRequireSigmaToVertex(kTRUE);
bfc7c23b 758 fTrackCuts->SetMaxChi2PerClusterTPC(3.5);
f09766dd 759 fTrackCuts->SetMinNClustersTPC(70);
bfc7c23b 760
761 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
762
763 Bool_t flagPhotonicElec = kFALSE;
f09766dd 764 Bool_t flagConvinatElec = kFALSE;
bfc7c23b 765
5e65985c 766 //for(Int_t jTracks = itrack+1; jTracks<fESD->GetNumberOfTracks(); jTracks++){
767 for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){
bfc7c23b 768 AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
769 if (!trackAsso) {
770 printf("ERROR: Could not receive track %d\n", jTracks);
771 continue;
772 }
5e65985c 773 if(itrack==jTracks)continue;
774
f09766dd 775 Double_t dEdxAsso = -999., ptPrim=-999., ptAsso=-999., openingAngle = -999.;
bfc7c23b 776 Double_t mass=999., width = -999;
777 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
778
f09766dd 779 ptPrim = track->Pt();
780
bfc7c23b 781 dEdxAsso = trackAsso->GetTPCsignal();
782 ptAsso = trackAsso->Pt();
783 Int_t chargeAsso = trackAsso->Charge();
784 Int_t charge = track->Charge();
785
5e65985c 786 //if(ptAsso <0.3) continue;
787 if(ptAsso <0.5) continue;
bfc7c23b 788 if(!fTrackCuts->AcceptTrack(trackAsso)) continue;
f09766dd 789 if(dEdxAsso <65 || dEdxAsso>100) continue; //11a pass1
bfc7c23b 790
791 Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
792 if(charge>0) fPDGe1 = -11;
793 if(chargeAsso>0) fPDGe2 = -11;
794
795 if(charge == chargeAsso) fFlagLS = kTRUE;
796 if(charge != chargeAsso) fFlagULS = kTRUE;
797
798 AliKFParticle ge1(*track, fPDGe1);
799 AliKFParticle ge2(*trackAsso, fPDGe2);
800 AliKFParticle recg(ge1, ge2);
801
802 if(recg.GetNDF()<1) continue;
803 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
804 if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
805
806 AliKFVertex primV(*pVtx);
807 primV += recg;
808 recg.SetProductionVertex(primV);
809
810 recg.SetMassConstraint(0,0.0001);
811
812 openingAngle = ge1.GetAngle(ge2);
813 if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
814 if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
815
816 if(openingAngle > fOpeningAngleCut) continue;
817
818 recg.GetMass(mass,width);
819
f09766dd 820 double phoinfo[3];
821 phoinfo[0] = cent;
822 phoinfo[1] = ptPrim;
823 phoinfo[2] = mass;
824
825 if(fFlagLS) fInvmassLS->Fill(phoinfo);
826 if(fFlagULS) fInvmassULS->Fill(phoinfo);
bfc7c23b 827
828 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
829 flagPhotonicElec = kTRUE;
830 }
f09766dd 831 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec){
832 flagConvinatElec = kTRUE;
833 }
bfc7c23b 834
835 }
836 fFlagPhotonicElec = flagPhotonicElec;
f09766dd 837 fFlagConvinatElec = flagConvinatElec;
bfc7c23b 838
839}
840
fb87d707 841
842//_________________________________________
843void AliAnalysisTaskHFECal::FindTriggerClusters()
844{
845 // constants
846 const int nModuleCols = 2;
847 const int nModuleRows = 5;
848 const int nColsFeeModule = 48;
849 const int nRowsFeeModule = 24;
850 const int nColsFaltroModule = 24;
851 const int nRowsFaltroModule = 12;
852 //const int faltroWidthMax = 20;
853
854 // part 1, trigger extraction -------------------------------------
855 Int_t globCol, globRow;
856 //Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0, trigInCut=0;
857 Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0;
858
859 //Int_t trigtimes[faltroWidthMax];
860 Double_t cellTime[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
861 Double_t cellEnergy[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
862 //Double_t fTrigCutLow = 6;
863 //Double_t fTrigCutHigh = 10;
864 Double_t fTimeCutLow = 469e-09;
865 Double_t fTimeCutHigh = 715e-09;
866
867 AliESDCaloTrigger * fCaloTrigger = fESD->GetCaloTrigger( "EMCAL" );
868
869 // erase trigger maps
870 for(Int_t i = 0; i < nColsFaltroModule*nModuleCols; i++ )
871 {
872 for(Int_t j = 0; j < nRowsFaltroModule*nModuleRows; j++ )
873 {
874 ftriggersCut[i][j] = 0;
875 ftriggers[i][j] = 0;
876 ftriggersTime[i][j] = 0;
877 }
878 }
879
880 Int_t iglobCol=0, iglobRow=0;
881 // go through triggers
882 if( fCaloTrigger->GetEntries() > 0 )
883 {
884 // needs reset
885 fCaloTrigger->Reset();
886 while( fCaloTrigger->Next() )
887 {
888 fCaloTrigger->GetPosition( globCol, globRow );
889 fCaloTrigger->GetNL0Times( ntimes );
890 /*
891 // no L0s
892 if( ntimes < 1 ) continue;
893 // get precise timings
894 fCaloTrigger->GetL0Times( trigtimes );
895 trigInCut = 0;
896 for(Int_t i = 0; i < ntimes; i++ )
897 {
898 // save the first trigger time in channel
899 if( i == 0 || triggersTime[globCol][globRow] > trigtimes[i] )
900 triggersTime[globCol][globRow] = trigtimes[i];
901 //printf("trigger times: %d\n",trigtimes[i]);
902 // check if in cut
903 if(trigtimes[i] > fTrigCutLow && trigtimes[i] < fTrigCutHigh )
904 trigInCut = 1;
905
906 fTrigTimes->Fill(trigtimes[i]);
907 }
908 */
909
910 //L1 analysis from AliAnalysisTaskEMCALTriggerQA
911 Int_t bit = 0;
912 fCaloTrigger->GetTriggerBits(bit);
913
914 Int_t ts = 0;
915 fCaloTrigger->GetL1TimeSum(ts);
916 if (ts > 0)ftriggers[globCol][globRow] = 1;
917 // number of triggered channels in event
918 nTrigChannel++;
919 // ... inside cut
920 if(ts>0 && (bit >> 6 & 0x1))
921 {
922 iglobCol = globCol;
923 iglobRow = globRow;
924 nTrigChannelCut++;
925 ftriggersCut[globCol][globRow] = 1;
926 }
927
928 } // calo trigger entries
929 } // has calo trigger entries
930
931 // part 2 go through the clusters here -----------------------------------
932 Int_t nCluster=0, nCell=0, iCell=0, gCell=0;
933 Short_t cellAddr, nSACell, mclabel;
934 //Int_t nSACell, iSACell, mclabel;
935 Int_t iSACell;
42c75692 936 Double_t cellAmp=0, cellTimeT=0, clusterTime=0, efrac=0;
fb87d707 937 Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta, gphi, geta, feta, fphi;
938
939 TRefArray *fCaloClusters = new TRefArray();
940 fESD->GetEMCALClusters( fCaloClusters );
941 nCluster = fCaloClusters->GetEntries();
942
943
944 // save all cells times if there are clusters
945 if( nCluster > 0 ){
946 // erase time map
947 for(Int_t i = 0; i < nColsFeeModule*nModuleCols; i++ ){
948 for(Int_t j = 0; j < nRowsFeeModule*nModuleRows; j++ ){
949 cellTime[i][j] = 0.;
950 cellEnergy[i][j] = 0.;
951 }
952 }
953
954 // get cells
955 AliESDCaloCells *fCaloCells = fESD->GetEMCALCells();
956 //AliVCaloCells fCaloCells = fESD->GetEMCALCells();
957 nSACell = fCaloCells->GetNumberOfCells();
958 for(iSACell = 0; iSACell < nSACell; iSACell++ ){
959 // get the cell info *fCal
960 fCaloCells->GetCell( iSACell, cellAddr, cellAmp, cellTimeT , mclabel, efrac);
961
962 // get cell position
963 fGeom->GetCellIndex( cellAddr, nSupMod, nModule, nIphi, nIeta );
964 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
965
966 // convert co global phi eta
967 gphi = iphi + nRowsFeeModule*(nSupMod/2);
968 geta = ieta + nColsFeeModule*(nSupMod%2);
969
970 // save cell time and energy
971 cellTime[geta][gphi] = cellTimeT;
972 cellEnergy[geta][gphi] = cellAmp;
973
974 }
975 }
976
977 Int_t nClusterTrig, nClusterTrigCut;
978 UShort_t *cellAddrs;
979 Double_t clsE=-999, clsEta=-999, clsPhi=-999;
980 Float_t clsPos[3] = {0.,0.,0.};
981
982 for(Int_t icl=0; icl<fESD->GetNumberOfCaloClusters(); icl++)
983 {
984 AliESDCaloCluster *cluster = fESD->GetCaloCluster(icl);
985 if(!cluster || !cluster->IsEMCAL()) continue;
986
987 // get cluster cells
988 nCell = cluster->GetNCells();
989
990 // get cluster energy
991 clsE = cluster->E();
992
993 // get cluster position
994 cluster->GetPosition(clsPos);
995 TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
996 clsEta = clsPosVec.Eta();
997 clsPhi = clsPosVec.Phi();
998
999 // get the cell addresses
1000 cellAddrs = cluster->GetCellsAbsId();
1001
1002 // check if the cluster contains cell, that was marked as triggered
1003 nClusterTrig = 0;
1004 nClusterTrigCut = 0;
1005
1006 // loop the cells to check, if cluser in acceptance
1007 // any cluster with a cell outside acceptance is not considered
1008 for( iCell = 0; iCell < nCell; iCell++ )
1009 {
42c75692 1010 // check hot cell
1011 if(clsE>6.0)fCellCheck->Fill(clsE,cellAddrs[iCell]);
1012
fb87d707 1013 // get cell position
1014 fGeom->GetCellIndex( cellAddrs[iCell], nSupMod, nModule, nIphi, nIeta );
1015 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
1016
1017 // convert co global phi eta
1018 gphi = iphi + nRowsFeeModule*(nSupMod/2);
1019 geta = ieta + nColsFeeModule*(nSupMod%2);
1020
1021 if( cellTime[geta][gphi] > 0. ){
1022 clusterTime += cellTime[geta][gphi];
1023 gCell++;
1024 }
1025
1026 // get corresponding FALTRO
1027 fphi = gphi / 2;
1028 feta = geta / 2;
1029
1030 // try to match with a triggered
1031 if( ftriggers[feta][fphi]==1)
1032 { nClusterTrig++;
1033 }
1034 if( ftriggersCut[feta][fphi]==1)
1035 { nClusterTrigCut++;
1036 }
1037
1038 } // cells
1039
1040
1041 if( gCell > 0 )
1042 clusterTime = clusterTime / (Double_t)gCell;
1043 // fix the reconstriction code time 100ns jumps
1044 if( fESD->GetBunchCrossNumber() % 4 < 2 )
1045 clusterTime -= 0.0000001;
1046
1047 fClsETime->Fill(clsE,clusterTime);
1048 fClsEBftTrigCut->Fill(clsE);
1049
1050 if(nClusterTrig>0){
1051 fClsETime1->Fill(clsE,clusterTime);
1052 }
1053
1054 if(nClusterTrig>0){
1055 cluster->SetChi2(1);
1056 fClsEAftTrigCut1->Fill(clsE);
1057 }
1058
1059 if(nClusterTrigCut>0){
1060 cluster->SetChi2(2);
1061 fClsEAftTrigCut2->Fill(clsE);
1062 }
1063
1064 if(nClusterTrigCut>0 && ( nCell > (1 + clsE / 3)))
1065 {
1066 cluster->SetChi2(3);
1067 fClsEAftTrigCut3->Fill(clsE);
1068 }
1069
1070 if(nClusterTrigCut>0 && (nCell > (1 + clsE / 3) )&&( clusterTime > fTimeCutLow && clusterTime < fTimeCutHigh ))
1071 {
1072 // cluster->SetChi2(4);
1073 fClsEAftTrigCut4->Fill(clsE);
1074 }
1075 if(nClusterTrigCut<1)
1076 {
1077 cluster->SetChi2(0);
1078
1079 fClsEAftTrigCut->Fill(clsE);
1080 }
1081
1082 } // clusters
1083}
1084
1085
1086
1087
1088
1089
1090