1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: Dmitri Peressounko (RRC KI) *
6 * Permission to use, copy, modify and distribute this software and its *
7 * documentation strictly for non-commercial purposes is hereby granted *
8 * without fee, provided that the above copyright notice appears in all *
9 * copies and that both the copyright notice and this permission notice *
10 * appear in the supporting documentation. The authors make no claims *
11 * about the suitability of this software for any purpose. It is *
12 * provided "as is" without express or implied warranty. *
13 **************************************************************************/
15 ////////////////////////////////////////////////
16 //---------------------------------------------
17 // Class used to do analysis on conversion + calorimeter pairs
18 // A lot of code cut-and-pasted from AliV0Reader and GammaConversion classes
20 //---------------------------------------------
21 ////////////////////////////////////////////////
27 #include "TDirectory.h"
28 #include "TLorentzVector.h"
32 #include "AliAnalysisManager.h"
33 #include "AliESDInputHandler.h"
34 #include "AliAnalysisTaskCaloConv.h"
37 #include "AliESDEvent.h"
38 #include "AliESDpid.h"
39 #include "AliESDtrackCuts.h"
40 #include "AliESDtrackCuts.h"
41 #include "AliCFContainer.h" // for CF
42 #include "AliESDCaloCluster.h"
43 #include "AliPHOSGeoUtils.h"
44 #include "AliEMCALGeometry.h"
45 #include "AliCFContainer.h"
46 #include "AliMCEventHandler.h"
47 #include "AliMCEvent.h"
49 #include "AliKFParticle.h"
50 #include "AliKFVertex.h"
55 ClassImp(AliAnalysisTaskCaloConv)
58 AliAnalysisTaskCaloConv::AliAnalysisTaskCaloConv():
64 fOutputContainer(NULL),
65 fCFOutputContainer(NULL),
71 fTriggerCINT1B(kFALSE),
73 fMinOpeningAngleGhostCut(0.),
86 fnSigmaAboveElectronLine(5.),
87 fnSigmaBelowElectronLine(-3.),
88 fnSigmaAbovePionLine(0.),
89 fpnSigmaAbovePionLine(1.),
95 fchi2CutConversion(30.)
97 // Default constructor
99 for(Int_t i=0;i<nBin;i++){
105 for(Int_t i=0; i<6; i++){
106 snprintf(key,55,"PHOS_BadMap_mod%d",i) ;
107 fPHOSBadMap[i]=new TH2I(key,"Bad Modules map",64,0.,64.,56,0.,56.) ;
109 for(Int_t i=0; i<10; i++){
110 snprintf(key,55,"EMCAL_BadMap_mod%d",i) ;
111 fEMCALBadMap[i] = new TH2I(key,"Bad Modules map",24,0.,24.,48,0.,48.) ;
114 AliAnalysisTaskCaloConv::AliAnalysisTaskCaloConv(const char* name):
115 AliAnalysisTaskSE(name),
120 fOutputContainer(NULL),
121 fCFOutputContainer(NULL),
127 fTriggerCINT1B(kFALSE),
129 fMinOpeningAngleGhostCut(0.),
134 fBadDistCutPHOS(3.3),
135 fBadDistCutEMCAL(6.),
142 fnSigmaAboveElectronLine(5.),
143 fnSigmaBelowElectronLine(-3.),
144 fnSigmaAbovePionLine(0.),
145 fpnSigmaAbovePionLine(1.),
151 fchi2CutConversion(30.)
153 // Common I/O in slot 0
154 DefineInput (0, TChain::Class());
155 DefineOutput(0, TTree::Class());
157 // Your private output
158 DefineOutput(1, TList::Class());
159 DefineOutput(2, TList::Class()); // for CF
162 for(Int_t i=0;i<nBin;i++){
168 for(Int_t i=0; i<6; i++){
169 snprintf(key,55,"PHOS_BadMap_mod%d",i) ;
170 fPHOSBadMap[i]=new TH2I(key,"Bad Modules map",64,0.,64.,56,0.,56.) ;
172 for(Int_t i=0; i<10; i++){
173 snprintf(key,55,"EMCAL_BadMap_mod%d",i) ;
174 fEMCALBadMap[i] = new TH2I(key,"Bad Modules map",24,0.,24.,48,0.,48.) ;
176 // fESDpid = new AliESDpid;
178 //_____________________________________________________
179 AliAnalysisTaskCaloConv::~AliAnalysisTaskCaloConv()
181 // Remove all pointers
183 if (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
184 AliAnalysisManager::kProofAnalysis) {
186 if(fOutputContainer){
187 fOutputContainer->Clear() ;
188 delete fOutputContainer ;
190 if(fCFOutputContainer){
191 fCFOutputContainer->Clear() ;
192 delete fCFOutputContainer ;
205 for(Int_t ivtx=0; ivtx<10; ivtx++){
206 if(fPHOSEvents[ivtx]){
207 delete fPHOSEvents[ivtx] ;
208 fPHOSEvents[ivtx]=0x0 ;
210 if(fEMCALEvents[ivtx]){
211 delete fEMCALEvents[ivtx] ;
212 fEMCALEvents[ivtx]=0x0 ;
214 if(fConvEvents[ivtx]){
215 delete fConvEvents[ivtx] ;
216 fConvEvents[ivtx]=0x0 ;
219 for(Int_t i=0; i<6; i++)
221 delete fPHOSBadMap[i] ;
224 for(Int_t i=0; i<10; i++)
226 delete fEMCALBadMap[i];
231 //_____________________________________________________
232 void AliAnalysisTaskCaloConv::Init()
235 // AliLog::SetGlobalLogLevel(AliLog::kError);
237 //_____________________________________________________
238 void AliAnalysisTaskCaloConv::UserExec(Option_t */*option*/)
240 // Execute analysis for current event
241 // First select conversion and calorimeter photons
242 // then construct inv. mass distributions
243 //First try to find Stack information.
244 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()){
245 if(static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent())
246 fStack = static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent()->Stack();
250 AliESDInputHandler *esdHandler=dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
252 if( esdHandler && esdHandler->GetESDpid()){
253 fESDpid=new AliESDpid(*(esdHandler->GetESDpid())) ;
256 fESDpid=new AliESDpid;
257 Double_t alephParameters[5];
258 if(fStack){// simulation
259 alephParameters[0] = 2.15898e+00/50.;
260 alephParameters[1] = 1.75295e+01;
261 alephParameters[2] = 3.40030e-09;
262 alephParameters[3] = 1.96178e+00;
263 alephParameters[4] = 3.91720e+00;
264 fESDpid->GetTOFResponse().SetTimeResolution(80.);
267 alephParameters[0] = 0.0283086;
268 alephParameters[1] = 2.63394e+01;
269 alephParameters[2] = 5.04114e-11;
270 alephParameters[3] = 2.12543e+00;
271 alephParameters[4] = 4.88663e+00;
272 fESDpid->GetTOFResponse().SetTimeResolution(130.);
273 fESDpid->GetTPCResponse().SetMip(47.9);
276 fESDpid->GetTPCResponse().SetBetheBlochParameters(
277 alephParameters[0],alephParameters[1],alephParameters[2],
278 alephParameters[3],alephParameters[4]);
279 fESDpid->GetTPCResponse().SetSigma(3.79301e-03, 2.21280e+04);
284 // fESDtrackCuts= AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
285 fESDtrackCuts = new AliESDtrackCuts;
288 fESDtrackCuts->SetMinNClustersTPC(70);
289 fESDtrackCuts->SetMaxChi2PerClusterTPC(4);
290 fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
291 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
293 fESDtrackCuts->SetRequireITSRefit(kTRUE);
294 fESDtrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
295 AliESDtrackCuts::kAny);
296 // if(selPrimaries) {
297 // 7*(0.0026+0.0050/pt^1.01)
298 fESDtrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
300 fESDtrackCuts->SetMaxDCAToVertexZ(2);
301 fESDtrackCuts->SetDCAToVertex2D(kFALSE);
302 fESDtrackCuts->SetRequireSigmaToVertex(kFALSE);
307 fESDtrackCuts = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
309 fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
310 fESDtrackCuts->SetMinNClustersTPC(70);
311 fESDtrackCuts->SetMaxChi2PerClusterTPC(4);
312 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
313 fESDtrackCuts->SetRequireITSRefit(kTRUE);
314 // fESDtrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny); //TEMPORARY <-> REMOVE
319 fESDEvent=(AliESDEvent*)InputEvent();
320 FillHistogram("hEventsTrig",0.5) ;
321 Bool_t isSelected = esdHandler && ((esdHandler->IsEventSelected()& AliVEvent::kMB) == AliVEvent::kMB);
323 printf("Not selected !!!!! \n") ;
324 PostData(1, fOutputContainer);
327 FillHistogram("hEventsTrig",1.5) ;
329 //Take Only events with proper trigger
330 //No trigger in MC data => no check
331 // if(!fStack && !fESDEvent->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL")){
333 if(!fStack && !fESDEvent->IsTriggerClassFired("CINT1-B-NOPF-ALLNOTRD")){
334 PostData(1, fOutputContainer);
337 FillHistogram("hEventsTrig",2.5) ;
339 //checks if we have a prim vertex
340 if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) {
341 PostData(1, fOutputContainer);
344 FillHistogram("hEventsTrig",3.5) ;
346 if(TMath::Abs(fESDEvent->GetPrimaryVertex()->GetZ())>10.){
347 PostData(1, fOutputContainer);
350 FillHistogram("hEventsTrig",4.5) ;
353 //Calculate charged multiplicity
354 Int_t trackCounter = 0;
355 for (Int_t i=0;i<fESDEvent->GetNumberOfTracks();++i) {
356 AliESDtrack *track = new AliESDtrack(*fESDEvent->GetTrack(i)) ;
357 if(fESDtrackCuts->AcceptTrack(track) && TMath::Abs(track->Eta())< 0.9)
361 fCentr=trackCounter+0.5 ;
362 FillHistogram("hMult",fCentr) ;
364 //Init geometry if not done yet
367 //Select conversion and calorimeter photons
368 //Conversion photons should go first since they are used in calibration in SelectCALOPhotons()
369 SelectConvPhotons() ;
370 SelectPHOSPhotons() ;
371 SelectEMCALPhotons() ;
372 //Fill MC histograms if MC is present
376 PostData(1, fOutputContainer);
378 PostData(2, fCFOutputContainer); // for CF
382 //____________________________________________________________
383 void AliAnalysisTaskCaloConv::ConnectInputData(Option_t *option){
384 // see header file for documentation
386 AliAnalysisTaskSE::ConnectInputData(option);
389 //____________________________________________________________
390 void AliAnalysisTaskCaloConv::UserCreateOutputObjects()
392 //UserCreateOutputObjects
393 if(fDebug)gDirectory->Print() ;
394 // Create the output container
395 if(fOutputContainer != NULL){
396 delete fOutputContainer;
397 fOutputContainer = NULL;
399 fOutputContainer = new TList();
400 fOutputContainer->SetOwner(kTRUE);
402 if(fCFOutputContainer != NULL){
403 delete fCFOutputContainer;
404 fCFOutputContainer = NULL;
406 //===========Correction Framework ======================
408 fCFOutputContainer = new TList();
409 fCFOutputContainer->SetOwner(kTRUE);
412 Int_t iBin[3]={500,40,100};
413 fConvCFCont = new AliCFContainer("ConvContainer","container for converted photons", 23,3,iBin);
414 fConvCFCont->SetBinLimits(0,0.,50.);
415 fConvCFCont->SetBinLimits(1,-2.,2.) ;
416 fConvCFCont->SetBinLimits(2,0.,1.);
417 fCFOutputContainer->Add(fConvCFCont) ;
419 fPHOSCFCont = new AliCFContainer("PHOSContainer","container for PHOS photons", 10,2,iBin);
420 fPHOSCFCont->SetBinLimits(0,0.,50.);
421 fPHOSCFCont->SetBinLimits(1,-2.,2.) ;
422 fCFOutputContainer->Add(fPHOSCFCont) ;
424 fEMCALCFCont = new AliCFContainer("EMCALContainer","container for EMCAL photons", 10,2,iBin);
425 fEMCALCFCont->SetBinLimits(0,0.,50.);
426 fEMCALCFCont->SetBinLimits(1,-2.,2.) ;
427 fCFOutputContainer->Add(fEMCALCFCont) ;
429 fPi0CFCont = new AliCFContainer("Pi0Container","container for EMCAL photons", 10,2,iBin);
430 fPi0CFCont->SetBinLimits(0,0.,50.);
431 fPi0CFCont->SetBinLimits(1,-2.,2.) ;
432 fCFOutputContainer->Add(fPi0CFCont) ;
435 //========================================================
437 //Adding the histograms to the output container
438 Int_t firstRun= 125000 ;
439 Int_t lastRun = 135000 ;
440 Int_t nRuns =lastRun-firstRun+1 ;
443 fOutputContainer->Add(new TH2F("hRunTrigger","Triggers fired",nRuns,float(firstRun),float(lastRun),2,0.,2.)) ;
444 fOutputContainer->Add(new TH1F("hRunEvents","Events per run",nRuns,float(firstRun),float(lastRun))) ;
445 fOutputContainer->Add(new TH1F("hRunConvs","Conversion photons per run",nRuns,float(firstRun),float(lastRun))) ;
446 fOutputContainer->Add(new TH1F("hRunPHOS","PHOS photons per run",nRuns,float(firstRun),float(lastRun))) ;
447 fOutputContainer->Add(new TH1F("hRunEMCAL","EMCAL photons per run",nRuns,float(firstRun),float(lastRun))) ;
448 fOutputContainer->Add(new TH1F("hVtxBin","Vtx distribution",10,0.,10.)) ;
449 fOutputContainer->Add(new TH1F("hEvents","Events processed",1,0.,1.)) ;
450 fOutputContainer->Add(new TH1F("hEventsTrig","Events processed",10,0.,10.)) ;
452 fOutputContainer->Add(new TH2F("hQA_PHOS_mod1_soft","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
453 fOutputContainer->Add(new TH2F("hQA_PHOS_mod2_soft","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
454 fOutputContainer->Add(new TH2F("hQA_PHOS_mod3_soft","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
455 fOutputContainer->Add(new TH2F("hQA_PHOS_mod4_soft","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
456 fOutputContainer->Add(new TH2F("hQA_PHOS_mod5_soft","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
457 fOutputContainer->Add(new TH2F("hQA_PHOS_mod1_hard","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
458 fOutputContainer->Add(new TH2F("hQA_PHOS_mod2_hard","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
459 fOutputContainer->Add(new TH2F("hQA_PHOS_mod3_hard","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
460 fOutputContainer->Add(new TH2F("hQA_PHOS_mod4_hard","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
461 fOutputContainer->Add(new TH2F("hQA_PHOS_mod5_hard","number of clusters per cell",64,0.,64.,56,0.,56.)) ;
463 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM0_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
464 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM1_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
465 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM2_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
466 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM3_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
467 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM4_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
468 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM5_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
469 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM6_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
470 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM7_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
471 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM8_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
472 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM9_soft","number of clusters per cell",24,0.,24,48,0.,48.)) ;
473 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM0_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
474 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM1_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
475 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM2_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
476 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM3_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
477 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM4_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
478 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM5_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
479 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM6_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
480 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM7_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
481 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM8_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
482 fOutputContainer->Add(new TH2F("hQA_EMCAL_SM9_hard","number of clusters per cell",24,0.,24,48,0.,48.)) ;
484 fOutputContainer->Add(new TH2F("hQA_ConvPhiEta","Number of V0s phi eta",100,0.,TMath::TwoPi(),40,-1.5,1.5)) ;
486 fOutputContainer->Add(new TH2F("hdEdx","dEdx of acceptaed electrons",1000,0.,10.,150,0.,150.)) ;
488 fOutputContainer->Add(new TH1F("hMult","Multiplicity",200,0.,200.)) ;
492 //Calibration of PHOS
493 fOutputContainer->Add(new TH3F("PHOS_mod1_th1","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
494 fOutputContainer->Add(new TH3F("PHOS_mod2_th1","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
495 fOutputContainer->Add(new TH3F("PHOS_mod3_th1","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
496 fOutputContainer->Add(new TH3F("PHOS_mod4_th1","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
497 fOutputContainer->Add(new TH3F("PHOS_mod5_th1","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
499 fOutputContainer->Add(new TH3F("PHOS_mod1_th2","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
500 fOutputContainer->Add(new TH3F("PHOS_mod2_th2","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
501 fOutputContainer->Add(new TH3F("PHOS_mod3_th2","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
502 fOutputContainer->Add(new TH3F("PHOS_mod4_th2","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
503 fOutputContainer->Add(new TH3F("PHOS_mod5_th2","Inv.Mass distr. per channel",64,0.,64,56,0.,56,100,0.,0.5)) ;
506 //Vary Conversion cuts
507 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_OnFly","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
508 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Offline","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
509 fOutputContainer->Add(new TH3F("PHOS_Re_mvsPt_OnFly_mult","Mass vs pt",400,0.,1.,npt,0.,ptmax,150,0.,150.)) ;
510 fOutputContainer->Add(new TH3F("PHOS_Re_mvsPt_Offline_mult","Mass vs pt",400,0.,1.,npt,0.,ptmax,150,0.,150.)) ;
511 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
512 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
513 // fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
514 // fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
515 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
516 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
517 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
518 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
519 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
520 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
521 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
522 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
523 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
524 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_Wcut_Neu","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
525 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
526 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_Wcut_Neu","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
527 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_On_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
528 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Off_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
530 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_OnFly","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
531 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Offline","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
532 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
533 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
534 // fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
535 // fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
536 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
537 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
538 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
539 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
540 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
541 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
542 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
543 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
544 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
545 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
546 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_On_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
547 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Off_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
549 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_OnFly","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
550 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Offline","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
551 fOutputContainer->Add(new TH3F("EMCAL_Re_mvsPt_OnFly_mult","Mass vs pt",400,0.,1.,npt,0.,ptmax,30,0.,60.)) ;
552 fOutputContainer->Add(new TH3F("EMCAL_Re_mvsPt_Offline_mult","Mass vs pt",400,0.,1.,npt,0.,ptmax,30,0.,60.)) ;
553 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
554 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
555 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
556 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
557 // fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
558 // fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
559 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
560 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
561 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
562 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
563 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
564 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
565 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
566 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
567 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_On_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
568 fOutputContainer->Add(new TH2F("EMCAL_Re_mvsPt_Off_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
570 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_OnFly","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
571 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Offline","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
572 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
573 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_ArmQt","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
574 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
575 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_dEdx","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
576 // fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
577 // fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_Prob","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
578 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
579 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_R120","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
580 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
581 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_Z","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
582 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
583 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_chi","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
584 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
585 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_Eta","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
586 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_On_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
587 fOutputContainer->Add(new TH2F("EMCAL_Mi_mvsPt_Off_Wcut","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
589 //PHOS PID variations
590 fOutputContainer->Add(new TH3F("PHOS_Re_mvsPt_alpha","Mass vs pt vs PHOS E",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
591 fOutputContainer->Add(new TH3F("PHOS_Re_mvsPt_E","Mass vs pt vs PHOS E",400,0.,1.,npt,0.,ptmax,100,0.,10.)) ;
592 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_all","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
593 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_all_dist","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
594 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Disp","Mass vs pt, disp cut",400,0.,1.,npt,0.,ptmax)) ;
595 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_TOF","Mass vs pt, TOF cut",400,0.,1.,npt,0.,ptmax)) ;
596 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_Neutral","Mass vs pt, Neutral cut",400,0.,1.,npt,0.,ptmax)) ;
597 fOutputContainer->Add(new TH2F("PHOS_Re_mvsPt_DispNeutral","Mass vs pt, Disp & neutral cut",400,0.,1.,npt,0.,ptmax)) ;
599 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_all","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
600 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_all_dist","Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
601 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Disp","Mass vs pt, disp cut",400,0.,1.,npt,0.,ptmax)) ;
602 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_TOF","Mass vs pt, TOF cut",400,0.,1.,npt,0.,ptmax)) ;
603 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_Neutral","Mass vs pt, Neutral cut",400,0.,1.,npt,0.,ptmax)) ;
604 fOutputContainer->Add(new TH2F("PHOS_Mi_mvsPt_DispNeutral","Mass vs pt, Disp & neutral cut",400,0.,1.,npt,0.,ptmax)) ;
607 for(Int_t mod=1; mod<=5;mod++){
608 snprintf(key,155,"PHOS_Re_mvsPt_mod%d_single",mod) ;
609 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
610 snprintf(key,155,"PHOS_Re_mvsPt_mod%d_all",mod) ;
611 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
614 //Single photon spectrum
616 fOutputContainer->Add(new TH1F("Single_conv_OnFly","Single photon spectrum",npt,0.,ptmax)) ;
617 fOutputContainer->Add(new TH1F("Single_conv_Offline","Single photon spectrum",npt,0.,ptmax)) ;
618 fOutputContainer->Add(new TH1F("Single_conv_On_ArmQt","Single photon spectrum",npt,0.,ptmax)) ;
619 fOutputContainer->Add(new TH1F("Single_conv_Off_ArmQt","Single photon spectrum",npt,0.,ptmax)) ;
620 fOutputContainer->Add(new TH1F("Single_conv_On_dEdx","Single photon spectrum",npt,0.,ptmax)) ;
621 fOutputContainer->Add(new TH1F("Single_conv_Off_dEdx","Single photon spectrum",npt,0.,ptmax)) ;
622 // fOutputContainer->Add(new TH1F("Single_conv_On_Prob","Single photon spectrum",npt,0.,ptmax)) ;
623 // fOutputContainer->Add(new TH1F("Single_conv_Off_Prob","Single photon spectrum",npt,0.,ptmax)) ;
624 fOutputContainer->Add(new TH1F("Single_conv_On_R120","Single photon spectrum",npt,0.,ptmax)) ;
625 fOutputContainer->Add(new TH1F("Single_conv_Off_R120","Single photon spectrum",npt,0.,ptmax)) ;
626 fOutputContainer->Add(new TH1F("Single_conv_On_Z","Single photon spectrum",npt,0.,ptmax)) ;
627 fOutputContainer->Add(new TH1F("Single_conv_Off_Z","Single photon spectrum",npt,0.,ptmax)) ;
628 fOutputContainer->Add(new TH1F("Single_conv_On_chi","Single photon spectrum",npt,0.,ptmax)) ;
629 fOutputContainer->Add(new TH1F("Single_conv_Off_chi","Single photon spectrum",npt,0.,ptmax)) ;
630 fOutputContainer->Add(new TH1F("Single_conv_On_Eta","Single photon spectrum",npt,0.,ptmax)) ;
631 fOutputContainer->Add(new TH1F("Single_conv_Off_Eta","Single photon spectrum",npt,0.,ptmax)) ;
632 fOutputContainer->Add(new TH1F("Single_conv_On_Wcut","Single photon spectrum",npt,0.,ptmax)) ;
633 fOutputContainer->Add(new TH1F("Single_conv_Off_Wcut","Single photon spectrum",npt,0.,ptmax)) ;
636 fOutputContainer->Add(new TH2F("PHOS_single_all_mult","Single photon spectrum",npt,0.,ptmax,150,0.,150.)) ;
637 fOutputContainer->Add(new TH2F("PHOS_single_disp_mult","Single photon spectrum",npt,0.,ptmax,150,0.,150.)) ;
638 fOutputContainer->Add(new TH2F("PHOS_single_neu_mult","Single photon spectrum",npt,0.,ptmax,150,0.,150.)) ;
640 for(Int_t mod=1; mod<=5;mod++){
641 snprintf(key,155,"PHOS_single_mod%d_all",mod) ;
642 fOutputContainer->Add(new TH1F(key,"Single photon spectrum",npt,0.,ptmax)) ;
643 snprintf(key,155,"PHOS_single_mod%d_disp",mod) ;
644 fOutputContainer->Add(new TH1F(key,"Single photon spectrum",npt,0.,ptmax)) ;
645 snprintf(key,155,"PHOS_single_mod%d_neutral",mod) ;
646 fOutputContainer->Add(new TH1F(key,"Single photon spectrum",npt,0.,ptmax)) ;
647 snprintf(key,155,"PHOS_single_mod%d_dispneutral",mod) ;
648 fOutputContainer->Add(new TH1F(key,"Single photon spectrum",npt,0.,ptmax)) ;
649 snprintf(key,155,"PHOS_single_mod%d_dist1",mod) ;
650 fOutputContainer->Add(new TH1F(key,"Single photon spectrum",npt,0.,ptmax)) ;
651 snprintf(key,155,"PHOS_single_mod%d_dist2",mod) ;
652 fOutputContainer->Add(new TH1F(key,"Single photon spectrum",npt,0.,ptmax)) ;
655 for(Int_t mod=0; mod<4;mod++){
656 snprintf(key,155,"EMCAL_mod%d_th1",mod) ;
657 fOutputContainer->Add(new TH3F(key,"Inv.Mass distr. per channel",24,0.,24,48,0.,48,100,0.,0.5)) ;
658 snprintf(key,155,"EMCAL_mod%d_th2",mod) ;
659 fOutputContainer->Add(new TH3F(key,"Inv.Mass distr. per channel",24,0.,24,48,0.,48,100,0.,0.5)) ;
661 snprintf(key,155,"EMCAL_Re_mvsPt_mod%d_single",mod) ;
662 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
664 snprintf(key,155,"EMCAL_Re_mvsPt_mod%d_all",mod) ;
665 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
666 snprintf(key,155,"EMCAL_Re_mvsPt_mod%d_Disp",mod) ;
667 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
668 snprintf(key,155,"EMCAL_Re_mvsPt_mod%d_TOF",mod) ;
669 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
670 snprintf(key,155,"EMCAL_Re_mvsPt_mod%d_Neutral",mod) ;
671 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
672 snprintf(key,155,"EMCAL_Re_mvsPt_mod%d_DispNeutral",mod) ;
673 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
675 snprintf(key,155,"EMCAL_Mi_mvsPt_mod%d_all",mod) ;
676 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
677 snprintf(key,155,"EMCAL_Mi_mvsPt_mod%d_Disp",mod) ;
678 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
679 snprintf(key,155,"EMCAL_Mi_mvsPt_mod%d_TOF",mod) ;
680 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
681 snprintf(key,155,"EMCAL_Mi_mvsPt_mod%d_Neutral",mod) ;
682 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
683 snprintf(key,155,"EMCAL_Mi_mvsPt_mod%d_DispNeutral",mod) ;
684 fOutputContainer->Add(new TH2F(key,"Mass vs pt",400,0.,1.,npt,0.,ptmax)) ;
688 fOutputContainer->Add(new TH2F("hMC_CaloConv_pi0_unitEta","Primary #pi^{0}",npt,0.,ptmax,150,0.,150.)) ;
689 fOutputContainer->Add(new TH2F("hMC_CaloConv_eta_unitEta","Primary #pi^{0}",npt,0.,ptmax,150,0.,150.)) ;
690 fOutputContainer->Add(new TH1F("hMC_CaloConv_allpi0","Primary #pi^{0}",npt,0.,ptmax)) ;
691 fOutputContainer->Add(new TH1F("hMC_CaloConv_alleta","Primary #pi^{0}",npt,0.,ptmax)) ;
692 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0PHOSacc","#pi^{0} decayed in PHOS acc",npt,0.,ptmax)) ;
693 fOutputContainer->Add(new TH1F("hMC_CaloConv_etaPHOSacc","#pi^{0} decayed in PHOS acc",npt,0.,ptmax)) ;
694 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0EMCALacc","#pi^{0} decayed in EMCAL acc",npt,0.,ptmax)) ;
695 fOutputContainer->Add(new TH1F("hMC_CaloConv_etaEMCALacc","#pi^{0} decayed in EMCAL acc",npt,0.,ptmax)) ;
696 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_PHOS_conv","#pi^{0} decayed in PHOS acc asnd conv. photon",npt,0.,ptmax)) ;
697 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_PHOS_conv","#pi^{0} decayed in PHOS acc asnd conv. photon",npt,0.,ptmax)) ;
698 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_EMCAL_conv","#pi^{0} decayed in EMCAL acc asnd conv. photon",npt,0.,ptmax)) ;
699 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_EMCAL_conv","#pi^{0} decayed in EMCAL acc asnd conv. photon",npt,0.,ptmax)) ;
700 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_bothphot_conv","#pi^{0} both photons converted",npt,0.,ptmax)) ;
701 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_bothphot_conv","#pi^{0} both photons converted",npt,0.,ptmax)) ;
702 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_convPhotInCalo","#pi^{0} photon in calo converted",npt,0.,ptmax)) ;
703 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_convPhotInCalo","#pi^{0} photon in calo converted",npt,0.,ptmax)) ;
704 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
705 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
706 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
707 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
708 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_EMCALacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
709 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_EMCALacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
710 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_EMCALacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
711 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_EMCALacc","#pi^{0} photon converted and V0 found",npt,0.,ptmax)) ;
712 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
713 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
714 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_pid","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
715 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_pid","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
716 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
717 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
718 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_pid","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
719 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_pid","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
720 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_good","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
721 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_good","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
722 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_good","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
723 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_good","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
724 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_mod1","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
725 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_mod1","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
726 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_mod2","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
727 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_mod2","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
728 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_mod3","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
729 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_mod3","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
730 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_mod4","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
731 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_mod4","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
732 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_PHOSclu_mod5","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
733 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0onfly_PHOSclu_mod5","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
734 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_mod1","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
735 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_mod1","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
736 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_mod2","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
737 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_mod2","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
738 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_mod3","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
739 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_mod3","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
740 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_mod4","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
741 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_mod4","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
742 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_PHOSclu_mod5","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
743 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0offline_PHOSclu_mod5","#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
745 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0on_PHOSclu_ptRec","#pi^{0} V0 and cluster in PHOS found(rec pt)",npt,0.,ptmax)) ;
746 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0on_PHOSclu_ptRec","#pi^{0} V0 and cluster in PHOS found(rec pt)",npt,0.,ptmax)) ;
747 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0off_PHOSclu_ptRec","#pi^{0} V0 and cluster in PHOS found(rec pt)",npt,0.,ptmax)) ;
748 fOutputContainer->Add(new TH1F("hMC_CaloConv_eta_v0off_PHOSclu_ptRec","#pi^{0} V0 and cluster in PHOS found(rec pt)",npt,0.,ptmax)) ;
749 fOutputContainer->Add(new TH2F("hMC_CaloConv_pi0_v0on_PHOSclu_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax)) ;
750 fOutputContainer->Add(new TH2F("hMC_CaloConv_eta_v0on_PHOSclu_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax)) ;
751 fOutputContainer->Add(new TH2F("hMC_CaloConv_pi0_v0off_PHOSclu_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax)) ;
752 fOutputContainer->Add(new TH2F("hMC_CaloConv_eta_v0off_PHOSclu_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax)) ;
754 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0off_EMCALclu_ptRec","#pi^{0} V0 and cluster in EMCAL found (rec pt)",npt,0.,ptmax)) ;
755 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0on_EMCALclu_ptRec","#pi^{0} V0 and cluster in EMCAL found (rec pt)",npt,0.,ptmax)) ;
756 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_EMCALclu","#pi^{0} V0 and cluster in EMCAL found",npt,0.,ptmax)) ;
757 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_EMCALclu_pid","#pi^{0} V0 and cluster in EMCAL found",npt,0.,ptmax)) ;
758 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_EMCALclu","#pi^{0} V0 and cluster in EMCAL found",npt,0.,ptmax)) ;
759 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_EMCALclu_pid","#pi^{0} V0 and cluster in EMCAL found",npt,0.,ptmax)) ;
760 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0onfly_EMCALclu_good","#pi^{0} V0 and cluster in EMCAL found",npt,0.,ptmax)) ;
761 fOutputContainer->Add(new TH1F("hMC_CaloConv_pi0_v0offline_EMCALclu_good","#pi^{0} V0 and cluster in EMCAL found",npt,0.,ptmax)) ;
762 fOutputContainer->Add(new TH2F("hMC_CaloConv_pi0_v0on_EMCALclu_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax)) ;
763 fOutputContainer->Add(new TH2F("hMC_CaloConv_pi0_v0off_EMCALclu_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax)) ;
765 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_Phot_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
766 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_Pi0_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
767 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_eta_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
768 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_K_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
769 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_pi_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
770 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_pbar_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
771 fOutputContainer->Add(new TH3F("hMC_Resid_PHOS_other_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
772 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_Phot_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
773 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_Pi0_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
774 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_eta_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
775 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_K_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
776 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_pi_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
777 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_pbar_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
778 fOutputContainer->Add(new TH3F("hMC_Resid_EMCAL_other_mvsPt","m vs pt for rec pi0s",400,0.,1.,npt,0.,ptmax,10,0.,1.)) ;
781 fOutputContainer->Add(new TH1F("hMC_CaloConv_phot","Primary photons",npt,0.,ptmax)) ;
782 fOutputContainer->Add(new TH1F("hMC_CaloConv_gammaPHOSacc","Photons in PHOS acc",npt,0.,ptmax)) ;
783 fOutputContainer->Add(new TH1F("hMC_CaloConv_gammaEMCALacc","Photons in EMCAL acc",npt,0.,ptmax)) ;
784 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_conv","Converted photons",npt,0.,ptmax)) ;
785 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_v0","Converted photons with V0",npt,0.,ptmax)) ;
786 fOutputContainer->Add(new TH2F("hMC_CaloConv_gamma_v0_devsE","Converted photons with V0",200,-1.,1.,npt,0.,ptmax)) ;
787 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_PHOSclu","Photons with cluster in PHOS",npt,0.,ptmax)) ;
788 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_PHOSclu_dist1","Photons with cluster in PHOS",npt,0.,ptmax)) ;
789 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_PHOSclu_dist2","Photons with cluster in PHOS",npt,0.,ptmax)) ;
790 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_EMCALclu","Photons with cluster in EMCAL",npt,0.,ptmax)) ;
791 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_EMCALclu_dist1","Photons with cluster in EMCAL",npt,0.,ptmax)) ;
792 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_EMCALclu_dist2","Photons with cluster in EMCAL",npt,0.,ptmax)) ;
793 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_PHOSclu_recE","Photons with cluster in PHOS",npt,0.,ptmax)) ;
794 fOutputContainer->Add(new TH1F("hMC_CaloConv_gamma_EMCALclu_recE","Photons with cluster in EMCAL",npt,0.,ptmax)) ;
795 fOutputContainer->Add(new TH2F("hMC_CaloConv_gamma_PHOSclu_devsE","Photons with cluster in PHOS",200,-1.,1.,npt,0.,ptmax)) ;
796 fOutputContainer->Add(new TH2F("hMC_CaloConv_gamma_EMCALclu_devsE","Photons with cluster in EMCAL",200,-1.,1.,npt,0.,ptmax)) ;
801 for(Int_t iw=0;iw<10;iw++){ //resolution
802 for(Int_t in=0;in<10;in++){
803 snprintf(keym,55,"hMC_nonlinearity_w%d_n%d",iw,in) ;
804 fOutputContainer->Add(new TH2F(keym,"m vs pt, nonlinearity test" ,200,0.,0.5,npt,0.,ptmax)) ;
805 snprintf(keym,55,"hMC_nonlinearity_ConvPHOS_w%d_n%d",iw,in) ;
806 fOutputContainer->Add(new TH2F(keym,"m vs pt, nonlinearity test" ,200,0.,0.5,npt,0.,ptmax)) ;
807 snprintf(keym,55,"hMC_nonlinearity_EMCAL_w%d_n%d",iw,in) ;
808 fOutputContainer->Add(new TH2F(keym,"m vs pt, nonlinearity test" ,200,0.,0.5,npt,0.,ptmax)) ;
809 snprintf(keym,55,"hMC_CaloConv_pi0_v0onfly_PHOSclu_w%d_n%d",iw,in) ;
810 fOutputContainer->Add(new TH1F(keym,"#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
811 snprintf(keym,55,"hMC_CaloConv_pi0_v0onfly_ConvPHOSclu_w%d_n%d",iw,in) ;
812 fOutputContainer->Add(new TH1F(keym,"#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
813 snprintf(keym,55,"hMC_CaloConv_pi0_v0onfly_EMCALclu_w%d_n%d",iw,in) ;
814 fOutputContainer->Add(new TH1F(keym,"#pi^{0} V0 and cluster in PHOS found",npt,0.,ptmax)) ;
819 fOutputContainer->Add(new TH3F("All_chi2_eta_pt","MC chi2 vs eta vs phi",100,0.,100.,200,-2.,2.,npt,0.,ptmax)) ;
820 fOutputContainer->Add(new TH2F("All_w_vs_m","MC w vs m",300,0.,TMath::Pi(),400,0.,1.)) ;
821 fOutputContainer->Add(new TH3F("MC_V0_pt_eta_phi","MC pt vs eta vs phi",npt,0.,ptmax,200,-2.,2.,200,0.,TMath::TwoPi())) ;
822 fOutputContainer->Add(new TH3F("MC_V0_m_eta_pt","MC m vs eta vs phi",400,0.,1.,200,-2.,2.,npt,0.,ptmax)) ;
823 fOutputContainer->Add(new TH3F("MC_V0_chi2_eta_pt","MC chi2 vs eta vs phi",100,0.,100.,200,-2.,2.,npt,0.,ptmax)) ;
824 fOutputContainer->Add(new TH2F("MC_V0_w_vs_m","MC w vs m",300,0.,TMath::Pi(),400,0.,1.)) ;
826 fOutputContainer->SetName(GetName());
828 //______________________________________________________________________
829 void AliAnalysisTaskCaloConv::InitGeometry()
831 //If not done yet, create Geometry for PHOS and EMCAL
832 //and read misalignment matrixes from ESD/AOD (AOD not implemented yet)
835 if(fPHOSgeom && fEMCALgeom){ //already initialized
839 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent()) ;
841 AliFatal("Can not read geometry matrixes from ESD/AOD: NO ESD") ;
842 if(!fPHOSgeom){//reading PHOS matrixes
843 fPHOSgeom = new AliPHOSGeoUtils("IHEP","");
844 for(Int_t mod=0; mod<5; mod++){
846 const TGeoHMatrix* m=esd->GetPHOSMatrix(mod) ;
848 fPHOSgeom->SetMisalMatrix(m, mod) ;
853 fEMCALgeom = AliEMCALGeometry::GetInstance("EMCAL_FIRSTYEARV1");
854 for(Int_t mod=0; mod < (fEMCALgeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){ //<---Gustavo, could you check???
856 const TGeoHMatrix* m=esd->GetEMCALMatrix(mod) ;
858 fEMCALgeom->SetMisalMatrix(m, mod) ;
863 //________________________________________________________________
864 void AliAnalysisTaskCaloConv::SelectPHOSPhotons(){
866 // Loop over all CaloClusters
868 fPHOSEvent->Clear() ;
870 fPHOSEvent = new TClonesArray("TLorentzVector",10) ;
876 vtx[0] = fESDEvent->GetPrimaryVertex()->GetX();
877 vtx[1] = fESDEvent->GetPrimaryVertex()->GetY();
878 vtx[2] = fESDEvent->GetPrimaryVertex()->GetZ();
879 for (Int_t i=0; i<fESDEvent->GetNumberOfCaloClusters(); i++) {
880 AliESDCaloCluster * clu = fESDEvent->GetCaloCluster(i);
884 clu ->GetMomentum(p ,vtx);
888 if(clu->GetNCells()<=2){
892 Bool_t isNeutral = kTRUE ;
893 Bool_t isDispOK = kTRUE ;
894 Bool_t isTOFOK = kTRUE ;
897 isNeutral = clu->GetEmcCpvDistance()>5. ; //To be improved
899 Double_t l0=clu->GetM02(),l1=clu->GetM20() ;
900 if(l1>= 0 && l0>= 0 && l1 < 0.1 && l0 < 0.1) isDispOK=kFALSE ;
901 if(l1>= 0 && l0 > 0.5 && l1 < 0.1 && l0 < 1.5) isDispOK=kTRUE ;
902 if(l1>= 0 && l0 > 2.0 && l1 < 0.1 && l0 < 2.7) isDispOK=kFALSE ;
903 if(l1>= 0 && l0 > 2.7 && l1 < 0.1 && l0 < 4.0) isDispOK=kFALSE ;
904 if(l1 > 0.1 && l1 < 0.7 && l0 > 0.7 && l0 < 2.1) isDispOK=kTRUE ;
905 if(l1 > 0.1 && l1 < 0.3 && l0 > 3.0 && l0 < 5.0) isDispOK=kFALSE ;
906 if(l1 > 0.3 && l1 < 0.7 && l0 > 2.5 && l0 < 4.0) isDispOK=kFALSE ;
907 if(l1 > 0.7 && l1 < 1.3 && l0 > 1.0 && l0 < 1.6) isDispOK=kTRUE ;
908 if(l1 > 0.7 && l1 < 1.3 && l0 > 1.6 && l0 < 3.5) isDispOK=kTRUE ;
909 if(l1 > 1.3 && l1 < 3.5 && l0 > 1.3 && l0 < 3.5) isDispOK=kTRUE ;
911 Float_t xyz[3] = {0,0,0};
912 clu->GetPosition(xyz); //Global position in ALICE system
913 TVector3 global(xyz) ;
915 if(!fPHOSgeom->GlobalPos2RelId(global,relid)){
916 printf("PHOS_beyond: x=%f, y=%f, z=%f \n",xyz[0],xyz[1],xyz[2]) ;
922 if(!IsGoodChannel("PHOS",iMod,iX,iZ))
925 Bool_t closeToBad=(clu->GetDistanceToBadChannel()>fBadDistCutPHOS) ;
927 p.SetBit(kCaloPIDdisp,isDispOK) ;
928 p.SetBit(kCaloPIDtof,isTOFOK) ;
929 p.SetBit(kCaloPIDneutral,isNeutral) ;
930 p.SetBit(BIT(17+iMod),kTRUE) ;
931 p.SetBit(kCaloDistBad,closeToBad) ;
932 new((*fPHOSEvent)[inPHOS]) TLorentzVector(p) ;
933 fGammaPHOS[inPHOS] = i ;
937 //Single photon spectra
938 Double_t pt= p.Pt() ;
939 TString skey="PHOS_single_"; skey+="mod" ; skey+=iMod ; skey+="_all" ;
940 FillHistogram(skey,pt) ;
941 FillHistogram("PHOS_single_all_mult",pt,fCentr) ;
943 skey="PHOS_single_"; skey+="mod" ; skey+=iMod ; skey+="_disp" ;
944 FillHistogram("PHOS_single_disp_mult",pt,fCentr) ;
945 FillHistogram(skey,pt) ;
948 skey="PHOS_single_"; skey+="mod" ; skey+=iMod ; skey+="_neutral" ;
949 FillHistogram(skey,pt) ;
950 FillHistogram("PHOS_single_neu_mult",pt,fCentr) ;
952 if(isNeutral && isDispOK){
953 skey="PHOS_single_"; skey+="mod" ; skey+=iMod ; skey+="_dispneutral" ;
954 FillHistogram(skey,pt) ;
956 //Distance to bad channel
957 if(clu->GetDistanceToBadChannel()>fBadDistCutPHOS){
958 skey="PHOS_single_"; skey+="mod" ; skey+=iMod ; skey+="_dist1" ;
959 FillHistogram(skey,pt) ;
961 if(clu->GetDistanceToBadChannel()>2.*fBadDistCutPHOS){
962 skey="PHOS_single_"; skey+="mod" ; skey+=iMod ; skey+="_dist2" ;
963 FillHistogram(skey,pt) ;
968 skey="hQA_PHOS_mod"; skey+=iMod; skey+="_soft" ;
969 FillHistogram(skey,iX-0.5, iZ-0.5,1.) ;
971 skey="hQA_PHOS_mod"; skey+=iMod; skey+="_hard" ;
972 FillHistogram(skey,iX-0.5, iZ-0.5,1.) ;
976 //Fill histogams for calibration
977 if(clu->E()<fPi0Thresh1 ) continue;
978 for(Int_t iconv=0;iconv<fConvEvent->GetEntriesFast();iconv++){
979 TLorentzVector *gammaConv=static_cast<TLorentzVector*>(fConvEvent->At(iconv)) ;
981 skey="PHOS_"; skey+="mod" ; skey+=iMod ; skey+="_th1" ;
982 FillHistogram(skey,iX-0.5, iZ-0.5,pi0.M());
984 skey="PHOS_"; skey+="mod"; skey+=iMod ; skey+="_th2" ;
985 FillHistogram(skey,iX-0.5, iZ-0.5,pi0.M());
990 //____________________________________________________________
991 void AliAnalysisTaskCaloConv::SelectEMCALPhotons(){
993 // Loop over all CaloClusters
995 fEMCALEvent->Clear() ;
997 fEMCALEvent = new TClonesArray("TLorentzVector",10) ;
998 Int_t inEMCAL = 0 ; //, inEMCALRecal=0;
1002 Double_t vtx[3]={0.,0.,0.};
1003 vtx[0] = fESDEvent->GetPrimaryVertex()->GetX();
1004 vtx[1] = fESDEvent->GetPrimaryVertex()->GetY();
1005 vtx[2] = fESDEvent->GetPrimaryVertex()->GetZ();
1006 Int_t nEMCAL=0 ; //For QA
1007 for (Int_t i=0; i<fESDEvent->GetNumberOfCaloClusters(); i++) {
1008 AliESDCaloCluster * clu = fESDEvent->GetCaloCluster(i);
1012 TLorentzVector pRecal ;
1013 clu ->GetMomentum(p ,vtx);
1014 Bool_t isNeutral = kTRUE ;
1015 Bool_t isDispOK = kTRUE ;
1016 Bool_t isTOFOK = kTRUE ;
1021 isNeutral = clu->GetEmcCpvDistance()>10. ; //To be improved
1022 if(clu->GetTOF()>550.e-9 && clu->GetTOF()<750.e-9)
1026 Float_t phi = p.Phi();
1027 if(phi < 0) phi+=TMath::TwoPi();
1029 fEMCALgeom->GetAbsCellIdFromEtaPhi(p.Eta(),phi, absId);
1030 iMod=fEMCALgeom->GetSuperModuleNumber(absId) ;
1031 Int_t imod = -1, iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1;
1032 fEMCALgeom->GetCellIndex(absId,imod,iTower,iIphi,iIeta);
1033 fEMCALgeom->GetCellPhiEtaIndexInSModule(imod,iTower, iIphi, iIeta,iphi,ieta);
1034 if(imod<0 || imod>5){
1035 printf("EMCAL: Beyond the geometry!\n") ;
1036 printf("phi=%f, eta=%f, absId=%d, SM=%d \n",p.Eta(),phi, absId, imod) ;
1042 if(!IsGoodChannel("EMCAL",iMod,iX,iZ))
1044 p.SetBit(kCaloPIDdisp,isDispOK) ;
1045 p.SetBit(kCaloPIDtof,isTOFOK) ;
1046 p.SetBit(kCaloPIDneutral,isNeutral) ;
1047 p.SetBit(BIT(17+imod),kTRUE) ;
1048 new((*fEMCALEvent)[inEMCAL]) TLorentzVector(p) ;
1049 fGammaEMCAL[inEMCAL] = i ;
1053 //Fill QA histograms
1054 if(clu->E()>0.5 && iMod>=0){ //Sometimes modules is negative not found??
1055 TString skey="hQA_EMCAL_SM";skey+=iMod ; skey+="_soft" ;
1056 FillHistogram(skey,iX-0.5, iZ-0.5,1.) ;
1058 skey="hQA_EMCAL_SM";skey+=iMod ; skey+="_hard" ;
1059 FillHistogram(skey,iX-0.5, iZ-0.5,1.) ;
1063 //Fill histograms for recalibration
1064 if(clu->E()<fPi0Thresh1) continue ;
1065 for(Int_t iconv=0;iconv<fConvEvent->GetEntriesFast();iconv++){
1066 TLorentzVector *gammaConv=static_cast<TLorentzVector*>(fConvEvent->At(iconv)) ;
1068 TString skey="EMCAL_"; skey+="mod" ; skey+=iMod ; skey+="_th1" ;
1069 FillHistogram(skey,iX-0.5, iZ-0.5,pi0.M());
1070 if(clu->E()>fPi0Thresh2){
1071 skey="EMCAL_"; skey+="mod"; skey+=iMod ; skey+="_th2" ;
1072 FillHistogram(skey,iX-0.5, iZ-0.5,pi0.M());
1078 //______________________________________________________________________
1079 void AliAnalysisTaskCaloConv::SelectConvPhotons(){
1080 //Fill list of conversion photons
1081 //that is scan v0s and select photon-like
1083 //set some constants
1084 const Double_t cutSigmaMass=0.0001; //Constraint on photon mass
1085 const Bool_t useImprovedVertex=kTRUE ; //Use verted with converted photon?
1086 // const Double_t zrSlope = TMath::Tan(2*TMath::ATan(TMath::Exp(-fetaCut)));
1087 const Double_t zrSlope12 = TMath::Tan(2*TMath::ATan(TMath::Exp(-1.2)));
1088 const Double_t zrSlope09 = TMath::Tan(2*TMath::ATan(TMath::Exp(-0.9)));
1089 const Double_t zOffset = 7.;
1092 fConvEvent = new TClonesArray("TLorentzVector",10) ;
1094 fConvEvent->Clear() ;
1096 //No primary vertex in event => scip
1097 if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) {
1102 for(Int_t iv0=0; iv0<fESDEvent->GetNumberOfV0s();iv0++){
1103 AliESDv0 * v0 = fESDEvent->GetV0(iv0) ;
1105 AliESDtrack * pos = fESDEvent->GetTrack(v0->GetPindex()) ;
1106 AliESDtrack * neg = fESDEvent->GetTrack(v0->GetNindex()) ;
1107 const AliExternalTrackParam * paramPos = v0->GetParamP() ;
1108 const AliExternalTrackParam * paramNeg = v0->GetParamN() ;
1109 if(pos->GetSign() <0){//change tracks
1111 neg=fESDEvent->GetTrack(v0->GetPindex()) ;
1113 paramNeg=v0->GetParamP() ;
1115 AliKFParticle negKF(*paramNeg,11);
1116 AliKFParticle posKF(*paramPos,-11);
1117 AliKFParticle photKF(negKF,posKF) ;
1118 photKF.SetMassConstraint(0,cutSigmaMass);
1120 if(useImprovedVertex){
1121 AliKFVertex primaryVertexImproved(*(fESDEvent->GetPrimaryVertex()));
1123 if(primaryVertexImproved.GetNContributors()>1){
1124 primaryVertexImproved+=photKF;
1125 photKF.SetProductionVertex(primaryVertexImproved);
1128 Double_t m=0., width=0. ;
1129 photKF.GetMass(m,width);
1131 TLorentzVector photLV;
1132 // photLV.SetXYZM(negKF.Px()+posKF.Px(),negKF.Py()+posKF.Py(),negKF.Pz()+negKF.Pz(),0.) ;
1133 // photLV.SetXYZT(photKF.GetPx(),photKF.GetPy(),photKF.GetPz(),photKF.GetE()) ;
1134 photLV.SetXYZM(photKF.GetPx(),photKF.GetPy(),photKF.GetPz(),0.) ; //Produces slightly better pi0 width
1136 //Parameters for correction function
1137 Double_t a[3]={photLV.Pt(),photLV.Eta(),m} ;
1139 fConvCFCont->Fill(a,0) ;
1142 Bool_t isOnFly=kTRUE ;
1144 if (v0->GetOnFlyStatus()){
1146 fConvCFCont->Fill(a,1) ;
1151 fConvCFCont->Fill(a,2) ;
1154 //Number of TPC clusters
1155 if(neg->GetNcls(1) <2 || pos->GetNcls(1) <2){
1159 //remove like sign pairs
1160 if(pos->GetSign() == neg->GetSign()){
1165 fConvCFCont->Fill(a,3) ;
1167 fConvCFCont->Fill(a,4) ;
1170 if( !(pos->GetStatus() & AliESDtrack::kTPCrefit) ||
1171 !(neg->GetStatus() & AliESDtrack::kTPCrefit) ){
1176 fConvCFCont->Fill(a,5) ;
1178 fConvCFCont->Fill(a,6) ;
1181 if( neg->GetKinkIndex(0) > 0 ||
1182 pos->GetKinkIndex(0) > 0) {
1187 if( fESDpid->NumberOfSigmasTPC(pos,AliPID::kElectron)<fnSigmaBelowElectronLine ||
1188 fESDpid->NumberOfSigmasTPC(pos,AliPID::kElectron)>fnSigmaAboveElectronLine ||
1189 fESDpid->NumberOfSigmasTPC(neg,AliPID::kElectron)<fnSigmaBelowElectronLine ||
1190 fESDpid->NumberOfSigmasTPC(neg,AliPID::kElectron)>fnSigmaAboveElectronLine ){
1193 const Double_t minPnSigmaAbovePionLine = 1. ;
1194 const Double_t maxPnSigmaAbovePionLine = 3. ;
1195 const Double_t nSigmaAbovePionLine = 0 ;
1196 if(pos->P()>minPnSigmaAbovePionLine && pos->P()<maxPnSigmaAbovePionLine ){
1197 if(fESDpid->NumberOfSigmasTPC(pos,AliPID::kPion)<nSigmaAbovePionLine){
1201 if(neg->P()>minPnSigmaAbovePionLine && neg->P()<maxPnSigmaAbovePionLine){
1202 if(fESDpid->NumberOfSigmasTPC(neg,AliPID::kPion)<nSigmaAbovePionLine){
1207 Bool_t isdEdx=kTRUE;
1208 if(pos->P()>minPnSigmaAbovePionLine && pos->P()<maxPnSigmaAbovePionLine ){
1209 if(fESDpid->NumberOfSigmasTPC(pos,AliPID::kPion)<2.){
1213 if(neg->P()>minPnSigmaAbovePionLine && neg->P()<maxPnSigmaAbovePionLine){
1214 if(fESDpid->NumberOfSigmasTPC(neg,AliPID::kPion)<2.){
1221 const Double_t minPKaonRejection=1.5 ;
1222 const Double_t sigmaAroundLine=1. ;
1223 if(neg->P()<minPKaonRejection ){
1224 if(TMath::Abs(fESDpid->NumberOfSigmasTPC(neg,AliPID::kKaon))<sigmaAroundLine){
1228 if(pos->P()<minPKaonRejection ){
1229 if(TMath::Abs(fESDpid->NumberOfSigmasTPC(pos,AliPID::kKaon))<sigmaAroundLine){
1235 const Double_t minPProtonRejection=2. ;
1236 if(neg->P()<minPProtonRejection){
1237 if(TMath::Abs(fESDpid->NumberOfSigmasTPC(neg,AliPID::kProton))<sigmaAroundLine){
1241 if(pos->P()<minPProtonRejection ){
1242 if(TMath::Abs(fESDpid->NumberOfSigmasTPC(pos,AliPID::kProton))<sigmaAroundLine){
1247 const Double_t minPPionRejection=0.5 ;
1248 if(neg->P()<minPPionRejection ){
1249 if(TMath::Abs(fESDpid->NumberOfSigmasTPC(neg,AliPID::kPion))<sigmaAroundLine){
1253 if(pos->P()<minPPionRejection ){
1254 if( TMath::Abs(fESDpid->NumberOfSigmasTPC(pos,AliPID::kPion))<sigmaAroundLine){
1261 FillHistogram("hdEdx",paramPos->GetP(),pos->GetTPCsignal()) ;
1262 FillHistogram("hdEdx",paramNeg->GetP(),neg->GetTPCsignal()) ;
1266 //Check the pid probability
1267 Bool_t isProb=kTRUE ;
1268 Double_t posProbArray[10];
1269 Double_t negProbArray[10];
1270 neg->GetTPCpid(negProbArray);
1271 pos->GetTPCpid(posProbArray);
1272 if(negProbArray[AliPID::kElectron]<fprobCut || posProbArray[AliPID::kElectron]<fprobCut){
1277 fConvCFCont->Fill(a,9) ;
1279 fConvCFCont->Fill(a,10) ;
1282 Double_t v0x=0.,v0y=0.,v0z=0.;
1283 v0->GetXYZ(v0x,v0y,v0z) ;
1284 Double_t r=TMath::Sqrt(v0x*v0x + v0y*v0y) ;
1286 const Double_t rMin=2.8 ;
1289 if(r>fmaxR){ // cuts on distance from collision point
1292 Bool_t isStrictR=kFALSE ;
1297 fConvCFCont->Fill(a,11) ;
1299 fConvCFCont->Fill(a,12) ;
1303 if((TMath::Abs(v0z)*zrSlope12)-zOffset > r ){ // cuts out regions where we do not reconstruct
1308 fConvCFCont->Fill(a,13) ;
1310 fConvCFCont->Fill(a,14) ;
1313 if(TMath::Abs(v0z) > fmaxZ ){ // cuts out regions where we do not reconstruct
1316 Bool_t isStrictZ=kFALSE ;
1317 if((TMath::Abs(v0z)*zrSlope09)-zOffset < r )
1322 fConvCFCont->Fill(a,15) ;
1324 fConvCFCont->Fill(a,16) ;
1327 if(photKF.GetNDF()<=0){
1332 fConvCFCont->Fill(a,17) ;
1334 fConvCFCont->Fill(a,18) ;
1337 Double_t chi2V0 = photKF.GetChi2()/photKF.GetNDF();
1338 FillHistogram("All_chi2_eta_pt",chi2V0,photLV.Eta(),photLV.Pt()) ;
1340 if(chi2V0 > fchi2CutConversion || chi2V0 <=0){
1343 Bool_t isStrictChi=kFALSE ;
1344 if(chi2V0 < 0.7*fchi2CutConversion && chi2V0 >0){
1350 fConvCFCont->Fill(a,19) ;
1352 fConvCFCont->Fill(a,20) ;
1355 const Double_t wideEtaCut=1.2 ;
1356 if(TMath::Abs(photLV.Eta())> wideEtaCut){
1359 if(TMath::Abs(paramPos->Eta())> wideEtaCut ||
1360 TMath::Abs(paramNeg->Eta())> wideEtaCut ){
1364 Bool_t isWideEta=kTRUE ;
1365 if(TMath::Abs(photLV.Eta())< fetaCut &&
1366 TMath::Abs(paramPos->Eta())<fetaCut &&
1367 TMath::Abs(paramNeg->Eta()) < fetaCut){
1372 if(photLV.Pt()<fptCut){
1377 fConvCFCont->Fill(a,21) ;
1379 fConvCFCont->Fill(a,22) ;
1384 if(photLV.Pt()>0.5){
1385 Double_t phi=photLV.Phi() ;
1386 while(phi<0.)phi+=TMath::TwoPi() ;
1387 while(phi>TMath::TwoPi())phi-=TMath::TwoPi() ;
1388 FillHistogram("hQA_ConvPhiEta",phi,photLV.Eta()) ;
1391 Double_t w=PlanarityAngle(paramPos,paramNeg) ;
1392 Bool_t isPlanarityCut = (0.08-0.22*w > m || 0.15*(w-2.4)>m) ;
1393 FillHistogram("All_w_vs_m",w,m) ;
1395 const Double_t armenterosAlphaCut=0.05 ;
1396 Double_t armenterosQtAlfa[2]={0.,0.} ;
1397 GetArmenterosQtAlfa(&posKF, &negKF, &photKF, armenterosQtAlfa ) ;
1398 Bool_t isArmQt=(armenterosQtAlfa[1]<armenterosAlphaCut) ;
1400 photLV.SetBit(kConvOnFly,isOnFly) ;
1401 photLV.SetBit(kConvArmQt,isArmQt) ;
1402 photLV.SetBit(kConvdEdx,isdEdx) ;
1403 photLV.SetBit(kConvProb,isProb) ;
1404 photLV.SetBit(kConvR,isStrictR) ;
1405 photLV.SetBit(kConvZR,isStrictZ) ;
1406 photLV.SetBit(kConvNDF,isStrictChi) ;
1407 photLV.SetBit(kConvEta,isWideEta) ;
1408 photLV.SetBit(kConvPlan,isPlanarityCut) ;
1410 new((*fConvEvent)[inConv]) TLorentzVector(photLV) ;
1411 fGammaV0s[inConv] = iv0 ;
1414 //Single photon spectrum
1415 Double_t pt=photLV.Pt() ;
1419 FillHistogram("Single_conv_OnFly",pt) ;
1420 if(isdEdx && !isWideEta)
1421 FillHistogram("Single_conv_On_dEdx",pt) ;
1422 if(!isWideEta && isStrictR)
1423 FillHistogram("Single_conv_On_R120",pt) ;
1424 if( !isWideEta && isStrictZ)
1425 FillHistogram("Single_conv_On_Z",pt) ;
1426 if(!isWideEta && isStrictChi)
1427 FillHistogram("Single_conv_On_chi",pt) ;
1429 FillHistogram("Single_conv_On_Eta",pt) ;
1430 if(!isWideEta && isPlanarityCut)
1431 FillHistogram("Single_conv_On_Wcut",pt) ;
1432 if(!isWideEta && isArmQt)
1433 FillHistogram("Single_conv_On_ArmQt",pt) ;
1437 FillHistogram("Single_conv_Offline",pt) ;
1438 if(isdEdx && !isWideEta)
1439 FillHistogram("Single_conv_Off_dEdx",pt) ;
1440 if(!isWideEta && isStrictR)
1441 FillHistogram("Single_conv_Off_R120",pt) ;
1442 if(!isWideEta && isStrictZ)
1443 FillHistogram("Single_conv_Off_Z",pt) ;
1444 if(!isWideEta && isStrictChi)
1445 FillHistogram("Single_conv_Off_chi",pt) ;
1447 FillHistogram("Single_conv_Off_Eta",pt) ;
1448 if(!isWideEta && isPlanarityCut)
1449 FillHistogram("Single_conv_Off_Wcut",pt) ;
1450 if(!isWideEta && isArmQt)
1451 FillHistogram("Single_conv_Off_ArmQt",pt) ;
1454 //Fill MC information
1456 TParticle * negativeMC = fStack->Particle(TMath::Abs(neg->GetLabel()));
1457 TParticle * positiveMC = fStack->Particle(TMath::Abs(pos->GetLabel()));
1459 if(!negativeMC || !positiveMC)
1462 if(negativeMC->GetMother(0) != positiveMC->GetMother(0))
1465 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1468 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
1472 TParticle * v0Gamma = fStack->Particle(negativeMC->GetMother(0));
1474 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){ // id 5 is conversion
1477 if(v0Gamma->GetPdgCode() == 22){
1478 FillHistogram("MC_V0_pt_eta_phi",v0Gamma->Pt(),v0Gamma->Eta(),v0Gamma->Phi()) ;
1479 FillHistogram("MC_V0_m_eta_pt",m,v0Gamma->Eta(),v0Gamma->Pt()) ;
1480 FillHistogram("MC_V0_chi2_eta_pt",chi2V0,v0Gamma->Eta(),v0Gamma->Pt()) ;
1481 FillHistogram("MC_V0_w_vs_m",w,m) ;
1486 //______________________________________________________________________
1487 void AliAnalysisTaskCaloConv::FillRealMixed(){
1488 // Fills real (same event) and Mixed (different events) inv.mass dsitributions
1489 // Moves current event to the list of stored events at the end
1492 vtx[0] = fESDEvent->GetPrimaryVertex()->GetX();
1493 vtx[1] = fESDEvent->GetPrimaryVertex()->GetY();
1494 vtx[2] = fESDEvent->GetPrimaryVertex()->GetZ();
1497 Int_t zvtx = (Int_t)((vtx[2]+10.)/2.) ;
1501 Double_t run = fESDEvent->GetRunNumber()+0.5;
1502 TString trigClasses = fESDEvent->GetFiredTriggerClasses();
1503 if(trigClasses.Contains("CINT1B-ABCE-NOPF-ALL"))
1504 FillHistogram("hRunTrigger",run,0.5) ;
1506 FillHistogram("hRunTrigger",run,1.5) ;
1508 FillHistogram("hEvents",0.5) ;
1509 FillHistogram("hVtxBin",zvtx-0.5) ;
1510 FillHistogram("hRunEvents",run) ;
1512 //check if containers for mixed events exist
1513 //create new if necessary
1514 if(!fPHOSEvents[zvtx]) fPHOSEvents[zvtx]=new TList() ;
1515 if(!fEMCALEvents[zvtx]) fEMCALEvents[zvtx]=new TList() ;
1516 if(!fConvEvents[zvtx]) fConvEvents[zvtx]=new TList() ;
1518 Int_t nPHOS=fPHOSEvent->GetEntriesFast() ;
1519 Int_t nEMCAL=fEMCALEvent->GetEntriesFast() ;
1520 Int_t nConv = fConvEvent->GetEntriesFast() ;
1521 //Some QA histograms
1522 //Calculate number of good converion photons
1524 for(Int_t iConv = 0; iConv<nConv; iConv++){
1525 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1526 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1530 FillHistogram("hRunConvs",run,double(nConvGood)) ;
1531 FillHistogram("hRunPHOS", run,double(nPHOS)) ;
1532 FillHistogram("hRunEMCAL",run,double(nEMCAL)) ;
1534 //Fill Real distributions
1535 for(Int_t iPHOS=0; iPHOS<nPHOS;iPHOS++){
1536 TLorentzVector * cal = static_cast<TLorentzVector*>(fPHOSEvent->At(iPHOS)) ;
1537 for(Int_t iConv = 0; iConv<nConv; iConv++){
1538 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1539 TLorentzVector pi=*cal + *cnv ;
1540 Double_t alpha=TMath::Abs(cal->Energy()-cnv->Energy())/(cal->Energy()+cnv->Energy()) ;
1541 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1542 FillHistogram("PHOS_Re_mvsPt_all",pi.M(),pi.Pt()) ;
1544 for(Int_t iw=0;iw<10;iw++){ //resolution
1545 for(Int_t in=0;in<10;in++){
1546 snprintf(keym,55,"hMC_nonlinearity_w%d_n%d",iw,in) ;
1547 Double_t mMod=0.,ptMod=0. ;
1548 Recalibrate(mMod, ptMod, cal, cnv, iw, in) ;
1549 FillHistogram(keym,mMod,ptMod) ;
1550 snprintf(keym,55,"hMC_nonlinearity_ConvPHOS_w%d_n%d",iw,in) ;
1551 RecalibrateConvPHOS(mMod, ptMod, cal, cnv, iw, in) ;
1552 FillHistogram(keym,mMod,ptMod) ;
1556 if(cal->TestBit(kCaloDistBad))
1557 FillHistogram("PHOS_Re_mvsPt_all_dist",pi.M(),pi.Pt()) ;
1558 FillHistogram("PHOS_Re_mvsPt_E",pi.M(),pi.Pt(),cal->Energy()) ;
1559 FillHistogram("PHOS_Re_mvsPt_alpha",pi.M(),pi.Pt(),alpha) ;
1560 if(cal->TestBit(kCaloPIDdisp))
1561 FillHistogram("PHOS_Re_mvsPt_Disp",pi.M(),pi.Pt()) ;
1562 if(cal->TestBit(kCaloPIDtof))
1563 FillHistogram("PHOS_Re_mvsPt_TOF",pi.M(),pi.Pt()) ;
1564 if(cal->TestBit(kCaloPIDneutral))
1565 FillHistogram("PHOS_Re_mvsPt_Neutral",pi.M(),pi.Pt()) ;
1566 if(cal->TestBit(kCaloPIDneutral) && cal->TestBit(kCaloPIDdisp))
1567 FillHistogram("PHOS_Re_mvsPt_DispNeutral",pi.M(),pi.Pt()) ;
1569 //Vary Conversion cuts
1570 if(cnv->TestBit(kConvOnFly)){
1572 if(!cnv->TestBit(kConvEta)){
1573 FillHistogram("PHOS_Re_mvsPt_OnFly",pi.M(),pi.Pt()) ;
1574 FillHistogram("PHOS_Re_mvsPt_OnFly_mult",pi.M(),pi.Pt(),fCentr) ;
1576 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1577 FillHistogram("PHOS_Re_mvsPt_On_dEdx",pi.M(),pi.Pt()) ;
1578 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1579 FillHistogram("PHOS_Re_mvsPt_On_R120",pi.M(),pi.Pt()) ;
1580 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1581 FillHistogram("PHOS_Re_mvsPt_On_Z",pi.M(),pi.Pt()) ;
1582 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1583 FillHistogram("PHOS_Re_mvsPt_On_chi",pi.M(),pi.Pt()) ;
1585 FillHistogram("PHOS_Re_mvsPt_On_Eta",pi.M(),pi.Pt()) ;
1586 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan)){
1587 FillHistogram("PHOS_Re_mvsPt_On_Wcut",pi.M(),pi.Pt()) ;
1588 if(cal->TestBit(kCaloPIDneutral))
1589 FillHistogram("PHOS_Re_mvsPt_On_Wcut_Neu",pi.M(),pi.Pt()) ;
1591 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1592 FillHistogram("PHOS_Re_mvsPt_On_ArmQt",pi.M(),pi.Pt()) ;
1596 if(!cnv->TestBit(kConvEta)){
1597 FillHistogram("PHOS_Re_mvsPt_Offline",pi.M(),pi.Pt()) ;
1598 FillHistogram("PHOS_Re_mvsPt_Offline_mult",pi.M(),pi.Pt(),fCentr) ;
1600 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1601 FillHistogram("PHOS_Re_mvsPt_Off_dEdx",pi.M(),pi.Pt()) ;
1602 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1603 FillHistogram("PHOS_Re_mvsPt_Off_R120",pi.M(),pi.Pt()) ;
1604 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1605 FillHistogram("PHOS_Re_mvsPt_Off_Z",pi.M(),pi.Pt()) ;
1606 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1607 FillHistogram("PHOS_Re_mvsPt_Off_chi",pi.M(),pi.Pt()) ;
1609 FillHistogram("PHOS_Re_mvsPt_Off_Eta",pi.M(),pi.Pt()) ;
1610 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan)){
1611 FillHistogram("PHOS_Re_mvsPt_Off_Wcut",pi.M(),pi.Pt()) ;
1612 if(cal->TestBit(kCaloPIDneutral))
1613 FillHistogram("PHOS_Re_mvsPt_Off_Wcut_Neu",pi.M(),pi.Pt()) ;
1615 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1616 FillHistogram("PHOS_Re_mvsPt_Off_ArmQt",pi.M(),pi.Pt()) ;
1620 //PHOS module-dependent histograms
1621 for(Int_t iPHOS=0; iPHOS<nPHOS;iPHOS++){
1622 TLorentzVector * cal = static_cast<TLorentzVector*>(fPHOSEvent->At(iPHOS)) ;
1624 while(!cal->TestBit(BIT(17+mod)) && mod<5)mod++ ;
1625 TString base("PHOS_Re_mvsPt_mod") ; base+=mod ;
1627 for(Int_t iConv = 0; iConv<nConv; iConv++){
1628 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1629 TLorentzVector pi=*cal + *cnv ;
1630 full=base ; full+="_single" ;
1631 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1632 FillHistogram(full,pi.M(),cal->Pt()) ;
1633 full=base ; full+="_all" ;
1634 FillHistogram(full,pi.M(),pi.Pt()) ;
1638 for(Int_t iEMCAL=0; iEMCAL<nEMCAL;iEMCAL++){
1639 TLorentzVector * cal = static_cast<TLorentzVector*>(fEMCALEvent->At(iEMCAL)) ;
1641 while(!cal->TestBit(BIT(17+mod)) && mod<6)mod++ ;
1642 TString base("EMCAL_Re_mvsPt_mod") ; base+=mod ;
1644 for(Int_t iConv = 0; iConv<nConv; iConv++){
1645 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1646 TLorentzVector pi=*cal + *cnv ;
1647 // Double_t alpha=TMath::Abs(cal->Energy()-cnv->Energy())/(cal->Energy()+cnv->Energy()) ;
1648 full=base+"_single" ;
1649 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1650 FillHistogram(full,pi.M(),cal->Pt()) ;
1652 FillHistogram(full,pi.M(),pi.Pt()) ;
1654 for(Int_t iw=0;iw<10;iw++){ //resolution
1655 for(Int_t in=0;in<10;in++){
1656 snprintf(keym,55,"hMC_nonlinearity_EMCAL_w%d_n%d",iw,in) ;
1657 Double_t mMod=0.,ptMod=0. ;
1658 RecalibrateEMCAL(mMod, ptMod, cal, cnv, iw, in) ;
1659 FillHistogram(keym,mMod,ptMod) ;
1663 //Vary Conversion cuts
1664 if(cnv->TestBit(kConvOnFly)){
1666 if(!cnv->TestBit(kConvEta)){
1667 FillHistogram("EMCAL_Re_mvsPt_OnFly",pi.M(),pi.Pt()) ;
1668 FillHistogram("EMCAL_Re_mvsPt_OnFly_mult",pi.M(),pi.Pt(),fCentr) ;
1670 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1671 FillHistogram("EMCAL_Re_mvsPt_On_dEdx",pi.M(),pi.Pt()) ;
1672 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1673 FillHistogram("EMCAL_Re_mvsPt_On_R120",pi.M(),pi.Pt()) ;
1674 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1675 FillHistogram("EMCAL_Re_mvsPt_On_Z",pi.M(),pi.Pt()) ;
1676 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1677 FillHistogram("EMCAL_Re_mvsPt_On_chi",pi.M(),pi.Pt()) ;
1679 FillHistogram("EMCAL_Re_mvsPt_On_Eta",pi.M(),pi.Pt()) ;
1680 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1681 FillHistogram("EMCAL_Re_mvsPt_On_Wcut",pi.M(),pi.Pt()) ;
1682 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1683 FillHistogram("EMCAL_Re_mvsPt_On_ArmQt",pi.M(),pi.Pt()) ;
1687 if(!cnv->TestBit(kConvEta)){
1688 FillHistogram("EMCAL_Re_mvsPt_Offline",pi.M(),pi.Pt()) ;
1689 FillHistogram("EMCAL_Re_mvsPt_Offline_mult",pi.M(),pi.Pt(),fCentr) ;
1691 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1692 FillHistogram("EMCAL_Re_mvsPt_Off_dEdx",pi.M(),pi.Pt()) ;
1693 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1694 FillHistogram("EMCAL_Re_mvsPt_Off_R120",pi.M(),pi.Pt()) ;
1695 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1696 FillHistogram("EMCAL_Re_mvsPt_Off_Z",pi.M(),pi.Pt()) ;
1697 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1698 FillHistogram("EMCAL_Re_mvsPt_Off_chi",pi.M(),pi.Pt()) ;
1700 FillHistogram("EMCAL_Re_mvsPt_Off_Eta",pi.M(),pi.Pt()) ;
1701 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1702 FillHistogram("EMCAL_Re_mvsPt_Off_Wcut",pi.M(),pi.Pt()) ;
1703 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1704 FillHistogram("EMCAL_Re_mvsPt_Off_ArmQt",pi.M(),pi.Pt()) ;
1709 TList * prevPHOS = fPHOSEvents[zvtx] ;
1710 TList * prevEMCAL = fEMCALEvents[zvtx] ;
1711 TList * prevConv = fConvEvents[zvtx] ;
1714 for(Int_t iPHOS=0; iPHOS<nPHOS;iPHOS++){
1715 TLorentzVector * cal = static_cast<TLorentzVector*>(fPHOSEvent->At(iPHOS)) ;
1716 for(Int_t ev=0; ev<prevConv->GetSize();ev++){
1717 TClonesArray * mixConv = static_cast<TClonesArray*>(prevConv->At(ev)) ;
1718 for(Int_t iConv = 0; iConv<mixConv->GetEntriesFast(); iConv++){
1719 TLorentzVector * cnv = static_cast<TLorentzVector*>(mixConv->At(iConv)) ;
1720 TLorentzVector pi=*cal + *cnv ;
1721 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1722 FillHistogram("PHOS_Mi_mvsPt_all",pi.M(),pi.Pt()) ;
1723 if(cal->TestBit(kCaloPIDdisp))
1724 FillHistogram("PHOS_Mi_mvsPt_Disp",pi.M(),pi.Pt()) ;
1725 if(cal->TestBit(kCaloPIDtof))
1726 FillHistogram("PHOS_Mi_mvsPt_TOF",pi.M(),pi.Pt()) ;
1727 if(cal->TestBit(kCaloPIDneutral))
1728 FillHistogram("PHOS_Mi_mvsPt_Neutral",pi.M(),pi.Pt()) ;
1729 if(cal->TestBit(kCaloPIDneutral) && cal->TestBit(kCaloPIDdisp))
1730 FillHistogram("PHOS_Mi_mvsPt_DispNeutral",pi.M(),pi.Pt()) ;
1732 //Vary Conversion cuts
1733 if(cnv->TestBit(kConvOnFly)){
1735 if(!cnv->TestBit(kConvEta))
1736 FillHistogram("PHOS_Mi_mvsPt_OnFly",pi.M(),pi.Pt()) ;
1737 if( cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1738 FillHistogram("PHOS_Mi_mvsPt_On_dEdx",pi.M(),pi.Pt()) ;
1739 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1740 FillHistogram("PHOS_Mi_mvsPt_On_R120",pi.M(),pi.Pt()) ;
1741 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1742 FillHistogram("PHOS_Mi_mvsPt_On_Z",pi.M(),pi.Pt()) ;
1743 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1744 FillHistogram("PHOS_Mi_mvsPt_On_chi",pi.M(),pi.Pt()) ;
1746 FillHistogram("PHOS_Mi_mvsPt_On_Eta",pi.M(),pi.Pt()) ;
1747 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1748 FillHistogram("PHOS_Mi_mvsPt_On_Wcut",pi.M(),pi.Pt()) ;
1749 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1750 FillHistogram("PHOS_Mi_mvsPt_On_ArmQt",pi.M(),pi.Pt()) ;
1754 if(!cnv->TestBit(kConvEta))
1755 FillHistogram("PHOS_Mi_mvsPt_Offline",pi.M(),pi.Pt()) ;
1756 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1757 FillHistogram("PHOS_Mi_mvsPt_Off_dEdx",pi.M(),pi.Pt()) ;
1758 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1759 FillHistogram("PHOS_Mi_mvsPt_Off_R120",pi.M(),pi.Pt()) ;
1760 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1761 FillHistogram("PHOS_Mi_mvsPt_Off_Z",pi.M(),pi.Pt()) ;
1762 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1763 FillHistogram("PHOS_Mi_mvsPt_Off_chi",pi.M(),pi.Pt()) ;
1765 FillHistogram("PHOS_Mi_mvsPt_Off_Eta",pi.M(),pi.Pt()) ;
1766 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1767 FillHistogram("PHOS_Mi_mvsPt_Off_Wcut",pi.M(),pi.Pt()) ;
1768 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1769 FillHistogram("PHOS_Mi_mvsPt_Off_ArmQt",pi.M(),pi.Pt()) ;
1774 for(Int_t iConv = 0; iConv<nConv; iConv++){
1775 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1776 for(Int_t ev=0; ev<prevPHOS->GetSize();ev++){
1777 TClonesArray * mixPHOS = static_cast<TClonesArray*>(prevPHOS->At(ev)) ;
1778 for(Int_t iPHOS=0; iPHOS<mixPHOS->GetEntriesFast();iPHOS++){
1779 TLorentzVector * cal = static_cast<TLorentzVector*>(mixPHOS->At(iPHOS)) ;
1780 TLorentzVector pi=*cal + *cnv ;
1781 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1782 FillHistogram("PHOS_Mi_mvsPt_all",pi.M(),pi.Pt()) ;
1783 if(cal->TestBit(kCaloPIDdisp))
1784 FillHistogram("PHOS_Mi_mvsPt_Disp",pi.M(),pi.Pt()) ;
1785 if(cal->TestBit(kCaloPIDtof))
1786 FillHistogram("PHOS_Mi_mvsPt_TOF",pi.M(),pi.Pt()) ;
1787 if(cal->TestBit(kCaloPIDneutral))
1788 FillHistogram("PHOS_Mi_mvsPt_Neutral",pi.M(),pi.Pt()) ;
1789 if(cal->TestBit(kCaloPIDneutral) && cal->TestBit(kCaloPIDdisp))
1790 FillHistogram("PHOS_Mi_mvsPt_DispNeutral",pi.M(),pi.Pt()) ;
1792 //Vary Conversion cuts
1793 if(cnv->TestBit(kConvOnFly)){
1795 if(!cnv->TestBit(kConvEta))
1796 FillHistogram("PHOS_Mi_mvsPt_OnFly",pi.M(),pi.Pt()) ;
1797 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1798 FillHistogram("PHOS_Mi_mvsPt_On_dEdx",pi.M(),pi.Pt()) ;
1799 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1800 FillHistogram("PHOS_Mi_mvsPt_On_R120",pi.M(),pi.Pt()) ;
1801 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1802 FillHistogram("PHOS_Mi_mvsPt_On_Z",pi.M(),pi.Pt()) ;
1803 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1804 FillHistogram("PHOS_Mi_mvsPt_On_chi",pi.M(),pi.Pt()) ;
1806 FillHistogram("PHOS_Mi_mvsPt_On_Eta",pi.M(),pi.Pt()) ;
1807 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1808 FillHistogram("PHOS_Mi_mvsPt_On_Wcut",pi.M(),pi.Pt()) ;
1809 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1810 FillHistogram("PHOS_Mi_mvsPt_On_ArmQt",pi.M(),pi.Pt()) ;
1814 if(!cnv->TestBit(kConvEta))
1815 FillHistogram("PHOS_Mi_mvsPt_Offline",pi.M(),pi.Pt()) ;
1816 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1817 FillHistogram("PHOS_Mi_mvsPt_Off_dEdx",pi.M(),pi.Pt()) ;
1818 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1819 FillHistogram("PHOS_Mi_mvsPt_Off_R120",pi.M(),pi.Pt()) ;
1820 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1821 FillHistogram("PHOS_Mi_mvsPt_Off_Z",pi.M(),pi.Pt()) ;
1822 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1823 FillHistogram("PHOS_Mi_mvsPt_Off_chi",pi.M(),pi.Pt()) ;
1825 FillHistogram("PHOS_Mi_mvsPt_Off_Eta",pi.M(),pi.Pt()) ;
1826 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1827 FillHistogram("PHOS_Mi_mvsPt_Off_Wcut",pi.M(),pi.Pt()) ;
1828 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1829 FillHistogram("PHOS_Mi_mvsPt_Off_ArmQt",pi.M(),pi.Pt()) ;
1836 //PHOS module dependent
1837 for(Int_t iPHOS=0; iPHOS<nPHOS;iPHOS++){
1838 TLorentzVector * cal = static_cast<TLorentzVector*>(fPHOSEvent->At(iPHOS)) ;
1840 while(!cal->TestBit(BIT(17+mod)) && mod<5)mod++ ;
1841 TString base("PHOS_Mi_mvsPt_mod") ; base+=mod ;
1843 for(Int_t ev=0; ev<prevConv->GetSize();ev++){
1844 TClonesArray * mixConv = static_cast<TClonesArray*>(prevConv->At(ev)) ;
1845 for(Int_t iConv = 0; iConv<mixConv->GetEntriesFast(); iConv++){
1846 TLorentzVector * cnv = static_cast<TLorentzVector*>(mixConv->At(iConv)) ;
1847 TLorentzVector pi=*cal + *cnv ;
1849 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1850 FillHistogram(full,pi.M(),pi.Pt()) ;
1855 for(Int_t iConv = 0; iConv<nConv; iConv++){
1856 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1857 for(Int_t ev=0; ev<prevPHOS->GetSize();ev++){
1858 TClonesArray * mixPHOS = static_cast<TClonesArray*>(prevPHOS->At(ev)) ;
1859 for(Int_t iPHOS=0; iPHOS<mixPHOS->GetEntriesFast();iPHOS++){
1860 TLorentzVector * cal = static_cast<TLorentzVector*>(mixPHOS->At(iPHOS)) ;
1862 while(!cal->TestBit(BIT(17+mod)) && mod<5)mod++ ;
1863 TString base("PHOS_Mi_mvsPt_mod") ; base+=mod ;
1865 TLorentzVector pi=*cal + *cnv ;
1867 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1868 FillHistogram(full,pi.M(),pi.Pt()) ;
1876 for(Int_t iEMCAL=0; iEMCAL<nEMCAL;iEMCAL++){
1877 TLorentzVector * cal = static_cast<TLorentzVector*>(fEMCALEvent->At(iEMCAL)) ;
1879 while(!cal->TestBit(BIT(17+mod)) && mod<6)mod++ ;
1880 TString base("EMCAL_Mi_mvsPt_mod") ; base+=mod ;
1882 for(Int_t ev=0; ev<prevConv->GetSize();ev++){
1883 TClonesArray * mixConv = static_cast<TClonesArray*>(prevConv->At(ev)) ;
1884 for(Int_t iConv = 0; iConv<mixConv->GetEntriesFast(); iConv++){
1885 TLorentzVector * cnv = static_cast<TLorentzVector*>(mixConv->At(iConv)) ;
1886 TLorentzVector pi=*cal + *cnv ;
1888 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1889 FillHistogram(full,pi.M(),pi.Pt()) ;
1891 if(cnv->TestBit(kConvOnFly)){
1893 if(!cnv->TestBit(kConvEta))
1894 FillHistogram("EMCAL_Mi_mvsPt_OnFly",pi.M(),pi.Pt()) ;
1895 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1896 FillHistogram("EMCAL_Mi_mvsPt_On_dEdx",pi.M(),pi.Pt()) ;
1897 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1898 FillHistogram("EMCAL_Mi_mvsPt_On_R120",pi.M(),pi.Pt()) ;
1899 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1900 FillHistogram("EMCAL_Mi_mvsPt_On_Z",pi.M(),pi.Pt()) ;
1901 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1902 FillHistogram("EMCAL_Mi_mvsPt_On_chi",pi.M(),pi.Pt()) ;
1904 FillHistogram("EMCAL_Mi_mvsPt_On_Eta",pi.M(),pi.Pt()) ;
1905 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1906 FillHistogram("EMCAL_Mi_mvsPt_On_Wcut",pi.M(),pi.Pt()) ;
1907 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1908 FillHistogram("EMCAL_Mi_mvsPt_On_ArmQt",pi.M(),pi.Pt()) ;
1912 if(!cnv->TestBit(kConvEta))
1913 FillHistogram("EMCAL_Mi_mvsPt_Offline",pi.M(),pi.Pt()) ;
1914 if( cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1915 FillHistogram("EMCAL_Mi_mvsPt_Off_dEdx",pi.M(),pi.Pt()) ;
1916 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1917 FillHistogram("EMCAL_Mi_mvsPt_Off_R120",pi.M(),pi.Pt()) ;
1918 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1919 FillHistogram("EMCAL_Mi_mvsPt_Off_Z",pi.M(),pi.Pt()) ;
1920 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1921 FillHistogram("EMCAL_Mi_mvsPt_Off_chi",pi.M(),pi.Pt()) ;
1923 FillHistogram("EMCAL_Mi_mvsPt_Off_Eta",pi.M(),pi.Pt()) ;
1924 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1925 FillHistogram("EMCAL_Mi_mvsPt_Off_Wcut",pi.M(),pi.Pt()) ;
1926 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1927 FillHistogram("EMCAL_Mi_mvsPt_Off_ArmQt",pi.M(),pi.Pt()) ;
1933 for(Int_t iConv = 0; iConv<nConv; iConv++){
1934 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
1935 for(Int_t ev=0; ev<prevEMCAL->GetSize();ev++){
1936 TClonesArray * mixEMCAL = static_cast<TClonesArray*>(prevEMCAL->At(ev)) ;
1937 for(Int_t iEMCAL=0; iEMCAL<mixEMCAL->GetEntriesFast();iEMCAL++){
1938 TLorentzVector * cal = static_cast<TLorentzVector*>(mixEMCAL->At(iEMCAL)) ;
1940 while(!cal->TestBit(BIT(17+mod)) && mod<6)mod++ ;
1941 TString base("EMCAL_Mi_mvsPt_mod") ; base+=mod ;
1943 TLorentzVector pi=*cal + *cnv ;
1945 if(cnv->TestBit(kConvOnFly) && !cnv->TestBit(kConvEta)){
1946 FillHistogram(full,pi.M(),pi.Pt()) ;
1947 if(cal->TestBit(kCaloPIDdisp)){
1949 FillHistogram(full,pi.M(),pi.Pt()) ;
1951 if(cal->TestBit(kCaloPIDtof)){
1953 FillHistogram(full,pi.M(),pi.Pt()) ;
1955 if(cal->TestBit(kCaloPIDneutral)){
1956 full=base+"_Neutral" ;
1957 FillHistogram(full,pi.M(),pi.Pt()) ;
1959 if(cal->TestBit(kCaloPIDneutral) && cal->TestBit(kCaloPIDdisp)){
1960 full=base+"_DispNeutral" ;
1961 FillHistogram(full,pi.M(),pi.Pt()) ;
1964 if(cnv->TestBit(kConvOnFly)){
1966 if(!cnv->TestBit(kConvEta))
1967 FillHistogram("EMCAL_Mi_mvsPt_OnFly",pi.M(),pi.Pt()) ;
1968 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1969 FillHistogram("EMCAL_Mi_mvsPt_On_dEdx",pi.M(),pi.Pt()) ;
1970 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1971 FillHistogram("EMCAL_Mi_mvsPt_On_R120",pi.M(),pi.Pt()) ;
1972 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1973 FillHistogram("EMCAL_Mi_mvsPt_On_Z",pi.M(),pi.Pt()) ;
1974 if( !cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1975 FillHistogram("EMCAL_Mi_mvsPt_On_chi",pi.M(),pi.Pt()) ;
1977 FillHistogram("EMCAL_Mi_mvsPt_On_Eta",pi.M(),pi.Pt()) ;
1978 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1979 FillHistogram("EMCAL_Mi_mvsPt_On_Wcut",pi.M(),pi.Pt()) ;
1980 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
1981 FillHistogram("EMCAL_Mi_mvsPt_On_ArmQt",pi.M(),pi.Pt()) ;
1985 if(!cnv->TestBit(kConvEta))
1986 FillHistogram("EMCAL_Mi_mvsPt_Offline",pi.M(),pi.Pt()) ;
1987 if(cnv->TestBit(kConvdEdx) && !cnv->TestBit(kConvEta))
1988 FillHistogram("EMCAL_Mi_mvsPt_Off_dEdx",pi.M(),pi.Pt()) ;
1989 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvR))
1990 FillHistogram("EMCAL_Mi_mvsPt_Off_R120",pi.M(),pi.Pt()) ;
1991 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvZR))
1992 FillHistogram("EMCAL_Mi_mvsPt_Off_Z",pi.M(),pi.Pt()) ;
1993 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvNDF))
1994 FillHistogram("EMCAL_Mi_mvsPt_Off_chi",pi.M(),pi.Pt()) ;
1996 FillHistogram("EMCAL_Mi_mvsPt_Off_Eta",pi.M(),pi.Pt()) ;
1997 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvPlan))
1998 FillHistogram("EMCAL_Mi_mvsPt_Off_Wcut",pi.M(),pi.Pt()) ;
1999 if(!cnv->TestBit(kConvEta) && cnv->TestBit(kConvArmQt))
2000 FillHistogram("EMCAL_Mi_mvsPt_Off_ArmQt",pi.M(),pi.Pt()) ;
2009 //Now we either add current events to stack or remove
2010 //Here we have some difficulty: conversion and calorimeter photons have different spectra.
2011 //So to correctly reproduce combinatorial background we have to preserve average number of
2012 //photons of both kinds per event. Therefore we should not reject empty PHOS/EMCAL events
2013 //though it will cost some memory. Reject only those events where no photons anywhere
2015 if((fPHOSEvent->GetEntriesFast()>0 || fEMCALEvent->GetEntriesFast()>0) && fConvEvent->GetEntriesFast()>0){
2016 prevPHOS->AddFirst(fPHOSEvent) ;
2018 prevEMCAL->AddFirst(fEMCALEvent) ;
2020 prevConv->AddFirst(fConvEvent) ;
2022 if(prevPHOS->GetSize()>100){//Remove redundant events
2023 TClonesArray * tmp = static_cast<TClonesArray*>(prevPHOS->Last()) ;
2024 prevPHOS->RemoveLast() ;
2026 tmp = static_cast<TClonesArray*>(prevEMCAL->Last()) ;
2027 prevEMCAL->RemoveLast() ;
2029 tmp = static_cast<TClonesArray*>(prevConv->Last()) ;
2030 prevConv->RemoveLast() ;
2036 //___________________________________________________________________________
2037 void AliAnalysisTaskCaloConv::ProcessMC(){
2039 //fill histograms for efficiensy etc. calculation
2040 if(!fStack) return ;
2042 const Double_t rcut = 1. ; //cut for primary particles
2044 vtx[0] = fESDEvent->GetPrimaryVertex()->GetX();
2045 vtx[1] = fESDEvent->GetPrimaryVertex()->GetY();
2046 vtx[2] = fESDEvent->GetPrimaryVertex()->GetZ();
2048 Int_t nPHOS=fPHOSEvent->GetEntriesFast() ;
2049 Int_t nEMCAL=fEMCALEvent->GetEntriesFast() ;
2050 Int_t nConv = fConvEvent->GetEntriesFast() ;
2052 //---------First pi0/eta-----------------------------
2055 for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
2056 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
2057 if(particle->GetPdgCode() == 111)
2058 snprintf(partName,10,"pi0") ;
2060 if(particle->GetPdgCode() == 221)
2061 snprintf(partName,10,"eta") ;
2066 if(particle->R() >rcut)
2069 Double_t pt = particle->Pt() ;
2070 //Total number of pi0 with creation radius <1 cm
2071 snprintf(hkey,55,"hMC_CaloConv_all%s",partName) ;
2072 FillHistogram(hkey,pt) ;
2073 if(TMath::Abs(particle->Y())<1.){
2074 snprintf(hkey,55,"hMC_CaloConv_%s_unitEta",partName) ;
2075 FillHistogram(hkey,pt,fCentr) ;
2078 //Check if one of photons converted
2079 if(particle->GetNDaughters()!=2)
2080 continue ; //Do not account Dalitz decays
2082 TParticle * gamma1 = fStack->Particle(particle->GetFirstDaughter());
2083 TParticle * gamma2 = fStack->Particle(particle->GetLastDaughter());
2084 //Number of pi0s decayed into acceptance
2085 Bool_t inAcc1 = (TMath::Abs(gamma1->Eta())<0.9) ;
2086 Bool_t inAcc2 = (TMath::Abs(gamma2->Eta())<0.9) ;
2088 Double_t x=0.,z=0. ;
2089 Bool_t hitPHOS1 = fPHOSgeom->ImpactOnEmc(gamma1, mod1, z,x) ;
2090 Bool_t hitPHOS2 = fPHOSgeom->ImpactOnEmc(gamma2, mod2, z,x) ;
2091 Bool_t hitEMCAL1= fEMCALgeom->Impact(gamma1) ;
2092 Bool_t hitEMCAL2= fEMCALgeom->Impact(gamma2) ;
2094 Bool_t goodPair=kFALSE ;
2095 if((inAcc1 && hitPHOS2) || (inAcc2 && hitPHOS1)){
2096 snprintf(hkey,55,"hMC_CaloConv_%sPHOSacc",partName) ;
2097 FillHistogram(hkey,pt) ;
2100 if((inAcc1 && hitEMCAL2) || (inAcc2 && hitEMCAL1)){
2101 snprintf(hkey,55,"hMC_CaloConv_%sEMCALacc",partName) ;
2102 FillHistogram(hkey,pt) ;
2109 Bool_t converted1 = kFALSE ;
2110 if(gamma1->GetNDaughters()==2){
2111 TParticle * e1=fStack->Particle(gamma1->GetFirstDaughter()) ;
2112 TParticle * e2=fStack->Particle(gamma1->GetLastDaughter()) ;
2113 if(TMath::Abs(e1->GetPdgCode())==11 && TMath::Abs(e2->GetPdgCode())==11){ //conversion
2115 converted1 = kTRUE ;
2118 Bool_t converted2 = kFALSE ;
2119 if(gamma2->GetNDaughters()==2){
2120 TParticle * e1=fStack->Particle(gamma2->GetFirstDaughter()) ;
2121 TParticle * e2=fStack->Particle(gamma2->GetLastDaughter()) ;
2122 if(TMath::Abs(e1->GetPdgCode())==11 && TMath::Abs(e2->GetPdgCode())==11){ //conversion
2124 converted2 = kTRUE ;
2128 //Number of pi0s with one photon converted
2129 if((converted1 && !converted2 && hitPHOS2) || (!converted1 && hitPHOS1 && converted2)) {
2130 snprintf(hkey,55,"hMC_CaloConv_%s_PHOS_conv",partName) ;
2131 FillHistogram(hkey,pt) ;
2134 if((converted1 && !converted2 && hitEMCAL2) || (!converted1 && hitEMCAL1 && converted2)) {
2135 snprintf(hkey,55,"hMC_CaloConv_%s_EMCAL_conv",partName) ;
2136 FillHistogram(hkey,pt) ;
2140 if(converted1 && converted2) {
2141 snprintf(hkey,55,"hMC_CaloConv_%s_bothphot_conv",partName) ;
2142 FillHistogram(hkey,pt) ;
2146 //photon pointing calorimeter converted
2147 if((converted1 && hitPHOS1 && !hitEMCAL2) || (converted2 && hitPHOS2 && !hitEMCAL1) ||
2148 (converted1 && hitEMCAL1 && !hitPHOS2) || (converted2 && hitEMCAL2 && !hitPHOS1)){
2149 snprintf(hkey,55,"hMC_CaloConv_%s_convPhotInCalo",partName) ;
2150 FillHistogram(hkey,pt) ;
2154 //Converted pi0 with v0 and photon PHOS or EMCAL
2155 Bool_t foundV01onfly=kFALSE, foundV01offline=kFALSE, foundV02onfly=kFALSE, foundV02offline=kFALSE ;
2156 Bool_t foundV01onflyPID=kFALSE, foundV01offlinePID=kFALSE, foundV02onflyPID=kFALSE, foundV02offlinePID=kFALSE ;
2157 TLorentzVector pConvOn,pConvOff ;
2158 for(Int_t iv0=0; iv0<fESDEvent->GetNumberOfV0s();iv0++){
2159 AliESDv0 * v0 = fESDEvent->GetV0(iv0) ;
2161 TParticle * negativeMC = fStack->Particle(TMath::Abs(fESDEvent->GetTrack(v0->GetNindex())->GetLabel()));
2162 TParticle * positiveMC = fStack->Particle(TMath::Abs(fESDEvent->GetTrack(v0->GetPindex())->GetLabel()));
2164 if(negativeMC->GetMother(0) != positiveMC->GetMother(0))
2167 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
2170 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
2174 TParticle * v0Gamma = fStack->Particle(negativeMC->GetMother(0));
2175 Bool_t same = (v0Gamma == gamma1) ;
2176 TParticle * tmp = v0Gamma ;
2177 while(!same && tmp->GetFirstMother()>=0){
2178 tmp = fStack->Particle(tmp->GetFirstMother());
2179 same = (tmp == gamma1) ;
2182 if(v0->GetOnFlyStatus())
2183 foundV01onfly = kTRUE ;
2185 foundV01offline= kTRUE ;
2186 for(Int_t iconv=0; iconv<nConv;iconv++){
2187 if(fGammaV0s[iconv] == iv0){
2188 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iconv)) ;
2190 if(!cnv->TestBit(kConvEta)){
2191 if(v0->GetOnFlyStatus()){
2193 foundV01onflyPID = kTRUE ;
2197 foundV01offlinePID = kTRUE ;
2205 same = (v0Gamma == gamma2) ;
2207 while(!same && tmp->GetFirstMother()>=0){
2208 tmp = fStack->Particle(tmp->GetFirstMother());
2209 same = (tmp == gamma2) ;
2212 if(v0->GetOnFlyStatus())
2213 foundV02onfly = kTRUE ;
2215 foundV02offline = kTRUE ;
2216 for(Int_t iconv=0; iconv<nConv;iconv++){
2217 if(fGammaV0s[iconv] == iv0){
2218 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iconv)) ;
2220 if(!cnv->TestBit(kConvEta)){
2221 if(v0->GetOnFlyStatus()){
2223 foundV02onflyPID = kTRUE ;
2227 foundV02offlinePID = kTRUE ;
2237 if((foundV01onfly && hitPHOS2) || (foundV02onfly && hitPHOS1)){
2238 snprintf(hkey,55,"hMC_CaloConv_%s_v0onfly_PHOSacc",partName) ;
2239 FillHistogram(hkey,pt) ;
2242 if((foundV01offline && hitPHOS2) || (foundV02offline && hitPHOS1)){
2243 snprintf(hkey,55,"hMC_CaloConv_%s_v0offline_PHOSacc",partName) ;
2244 FillHistogram(hkey,pt) ;
2247 if((foundV01onfly && hitEMCAL2) || (foundV02onfly && hitEMCAL1)){
2248 snprintf(hkey,55,"hMC_CaloConv_%s_v0onfly_EMCALacc",partName) ;
2249 FillHistogram(hkey,pt) ;
2252 if((foundV01offline && hitEMCAL2) || (foundV02offline && hitEMCAL1)){
2253 snprintf(hkey,55,"hMC_CaloConv_%s_v0offline_EMCALacc",partName) ;
2254 FillHistogram(hkey,pt) ;
2261 //Converted pi0 with v0 and cluster in PHOS/EMCAL
2262 Bool_t cluInPHOS = kFALSE,cluInEMCAL=kFALSE ;
2263 Bool_t cluInPHOSpid = kFALSE,cluInEMCALpid=kFALSE ;
2264 Bool_t closeToBad= kFALSE ;
2265 TLorentzVector pCalo ;
2266 for (Int_t i=0; i<fESDEvent->GetNumberOfCaloClusters(); i++) {
2267 AliESDCaloCluster * clu = fESDEvent->GetCaloCluster(i);
2268 Int_t iprim = clu->GetLabel() ; //# of particle hit PHOS/EMCAL
2269 Bool_t matched = kFALSE ;
2271 if(iprim==particle->GetFirstDaughter() || iprim==particle->GetLastDaughter()){
2276 iprim=fStack->Particle(iprim)->GetFirstMother() ;
2281 if(clu->IsPHOS() && (hitPHOS1 || hitPHOS2)){
2283 //Check if cluster passed PID
2284 for(Int_t inPHOS=0; inPHOS<nPHOS;inPHOS++){
2285 if(fGammaPHOS[inPHOS] == i){
2286 cluInPHOSpid=kTRUE ;
2290 clu->GetMomentum(pCalo ,vtx);
2291 if(clu->GetDistanceToBadChannel()<fBadDistCutPHOS)
2295 if(!clu->IsPHOS() && (hitEMCAL1 || hitEMCAL2)){
2297 //Check if cluster passed PID
2298 for(Int_t inEMCAL=0; inEMCAL<nEMCAL;inEMCAL++){
2299 if(fGammaEMCAL[inEMCAL] == i){
2300 cluInPHOSpid=kTRUE ;
2304 clu->GetMomentum(pCalo ,vtx);
2305 if(clu->GetDistanceToBadChannel()<fBadDistCutEMCAL)
2313 if(foundV01onfly ||foundV02onfly){
2314 snprintf(hkey,55,"hMC_CaloConv_%s_v0onfly_PHOSclu",partName) ;
2315 FillHistogram(hkey,pt) ;
2317 if((foundV01onflyPID ||foundV02onflyPID) && cluInPHOSpid){
2318 snprintf(hkey,55,"hMC_CaloConv_%s_v0onfly_PHOSclu_pid",partName) ;
2319 FillHistogram(hkey,pt) ;
2320 for(Int_t iw=0;iw<10;iw++){ //resolution
2321 for(Int_t in=0;in<10;in++){
2323 snprintf(keym,55,"hMC_CaloConv_%s_v0onfly_PHOSclu_w%d_n%d",partName,iw,in) ;
2324 Double_t mMod=0.,ptMod=0. ;
2325 Recalibrate(mMod, ptMod, &pCalo, &pConvOn, iw, in) ;
2326 FillHistogram(keym,ptMod) ;
2327 snprintf(keym,55,"hMC_CaloConv_%s_v0onfly_ConvPHOSclu_w%d_n%d",partName,iw,in) ;
2328 RecalibrateConvPHOS(mMod, ptMod, &pCalo, &pConvOn, iw, in) ;
2329 FillHistogram(keym,ptMod) ;
2333 snprintf(hkey,55,"hMC_CaloConv_%s_v0onfly_PHOSclu_good",partName) ;
2334 FillHistogram(hkey,pt) ;
2336 Double_t m=(pCalo+pConvOn).M() ;
2337 Double_t ptm=(pCalo+pConvOn).Pt() ;
2338 snprintf(hkey,55,"hMC_CaloConv_%s_v0on_PHOSclu_ptRec",partName) ;
2339 FillHistogram(hkey,ptm) ;
2340 snprintf(hkey,55,"hMC_CaloConv_%s_v0on_PHOSclu_mvsPt",partName) ;
2341 FillHistogram(hkey,m,ptm) ;
2346 if(foundV01offline ||foundV02offline){
2347 snprintf(hkey,55,"hMC_CaloConv_%s_v0offline_PHOSclu",partName) ;
2348 FillHistogram(hkey,pt) ;
2349 if((foundV01offlinePID ||foundV02offlinePID) && cluInPHOSpid){
2350 snprintf(hkey,55,"hMC_CaloConv_%s_v0offline_PHOSclu_pid",partName) ;
2351 FillHistogram(hkey,pt) ;
2352 Double_t m=(pCalo+pConvOff).M() ;
2353 Double_t ptm=(pCalo+pConvOff).Pt() ;
2354 snprintf(hkey,55,"hMC_CaloConv_%s_v0off_PHOSclu_ptRec",partName) ;
2355 FillHistogram(hkey,ptm) ;
2356 snprintf(hkey,55,"hMC_CaloConv_%s_v0off_PHOSclu_mvsPt",partName) ;
2357 FillHistogram(hkey,m,ptm) ;
2359 snprintf(hkey,55,"hMC_CaloConv_%s_v0offline_PHOSclu_good",partName) ;
2360 FillHistogram(hkey,pt) ;
2365 if((foundV01onflyPID ||foundV02onflyPID) && cluInPHOSpid){
2366 TString base("hMC_CaloConv_") ; base+=partName; base+="_v0onfly_PHOSclu_mod" ;
2371 FillHistogram(base.Data(),pt) ;
2374 if((foundV01offlinePID ||foundV02offlinePID) && cluInPHOSpid){
2375 TString base("hMC_CaloConv_") ; base+=partName; base+="_v0offline_PHOSclu_mod" ;
2380 FillHistogram(base.Data(),pt) ;
2383 if(cluInEMCAL && strcmp(partName,"pi0")==0){
2385 if(foundV01onfly ||foundV02onfly){
2386 FillHistogram("hMC_CaloConv_pi0_v0onfly_EMCALclu",pt) ;
2388 if((foundV01onflyPID ||foundV02onflyPID) && cluInEMCALpid){
2389 FillHistogram("hMC_CaloConv_pi0_v0onfly_EMCALclu_pid",pt) ;
2390 for(Int_t iw=0;iw<10;iw++){ //resolution
2391 for(Int_t in=0;in<10;in++){
2393 snprintf(keym,55,"hMC_CaloConv_pi0_v0onfly_EMCALclu_w%d_n%d",iw,in) ;
2394 Double_t mMod=0.,ptMod=0. ;
2395 RecalibrateEMCAL(mMod, ptMod, &pCalo, &pConvOn, iw, in) ;
2396 FillHistogram(keym,ptMod) ;
2400 FillHistogram("hMC_CaloConv_pi0_v0onfly_EMCALclu_good",pt) ;
2401 Double_t m=(pCalo+pConvOn).M() ;
2402 Double_t ptm=(pCalo+pConvOn).Pt() ;
2403 FillHistogram("hMC_CaloConv_pi0_v0on_EMCALclu_ptRec",ptm) ;
2404 FillHistogram("hMC_CaloConv_pi0_v0on_EMCALclu_mvsPt",m,ptm) ;
2409 if(foundV01offline ||foundV02offline){
2410 FillHistogram("hMC_CaloConv_pi0_v0offline_EMCALclu",pt) ;
2411 if((foundV01offlinePID ||foundV02offlinePID) && cluInEMCALpid){
2412 FillHistogram("hMC_CaloConv_pi0_v0offline_EMCALclu_pid",pt) ;
2413 Double_t m=(pCalo+pConvOff).M() ;
2414 Double_t ptm=(pCalo+pConvOff).Pt() ;
2415 FillHistogram("hMC_CaloConv_pi0_v0off_EMCALclu_ptRec",ptm) ;
2416 FillHistogram("hMC_CaloConv_pi0_v0off_EMCALclu_mvsPt",m,ptm) ;
2418 FillHistogram("hMC_CaloConv_pi0_v0offline_EMCALclu_good",pt) ;
2424 //Construct Inv mass distributions for residual correlations
2425 if(fPHOSEvent && fConvEvent){
2426 for(Int_t iPHOS=0; iPHOS<fPHOSEvent->GetEntriesFast();iPHOS++){
2427 TLorentzVector * cal = static_cast<TLorentzVector*>(fPHOSEvent->At(iPHOS)) ;
2428 Int_t iclu=fGammaPHOS[iPHOS] ;
2429 AliESDCaloCluster * clu = fESDEvent->GetCaloCluster(iclu);
2430 Int_t iprimPHOS = clu->GetLabel() ; //# of particle hit PHOS/EMCAL
2431 for(Int_t iConv = 0; iConv<fConvEvent->GetEntriesFast(); iConv++){
2432 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
2433 if(!cnv->TestBit(kConvOnFly) || cnv->TestBit(kConvEta))
2436 Int_t iv0=fGammaV0s[iConv] ;
2437 AliESDv0 * v0 = fESDEvent->GetV0(iv0) ;
2438 Int_t iprimNeg = TMath::Abs(fESDEvent->GetTrack(v0->GetNindex())->GetLabel()) ;
2439 Int_t iprimPos = TMath::Abs(fESDEvent->GetTrack(v0->GetPindex())->GetLabel()) ;
2441 //Check if there was a common ancistor
2442 Bool_t found = kFALSE ;
2443 Int_t curPHOS=iprimPHOS ;
2445 while(!found && curPHOS>-1){
2446 Int_t curNeg=iprimNeg ;
2447 while(!found && curNeg>-1){
2448 if(curNeg==curPHOS){
2453 curNeg=fStack->Particle(curNeg)->GetFirstMother() ;
2456 curPHOS=fStack->Particle(curPHOS)->GetFirstMother() ;
2461 while(!found && curPHOS>-1){
2462 Int_t curPos=iprimPos ;
2463 while(!found && curPos>-1){
2464 if(curPos==curPHOS){
2469 curPos=fStack->Particle(curPos)->GetFirstMother() ;
2472 curPHOS=fStack->Particle(curPHOS)->GetFirstMother() ;
2474 if(commonA != commonB){
2476 AliInfo(Form("CommonA=%d, commonB=%d",commonA,commonB)) ;
2478 if(commonA>-1){//There was common particles
2479 Int_t pdg = fStack->Particle(commonA)->GetPdgCode() ;
2480 TLorentzVector pi=*cal + *cnv ;
2482 Double_t pt=pi.Pt() ;
2483 Double_t alpha=TMath::Abs(cal->Energy()-cnv->Energy())/(cal->Energy()+cnv->Energy()) ;
2487 case 22: //conversion
2488 FillHistogram("hMC_Resid_PHOS_Phot_mvsPt",m,pt,alpha) ;
2491 FillHistogram("hMC_Resid_PHOS_Pi0_mvsPt",m,pt,alpha) ;
2494 FillHistogram("hMC_Resid_PHOS_eta_mvsPt",m,pt,alpha) ;
2500 FillHistogram("hMC_Resid_PHOS_K_mvsPt",m,pt,alpha) ;
2504 FillHistogram("hMC_Resid_PHOS_pi_mvsPt",m,pt,alpha) ;
2508 FillHistogram("hMC_Resid_PHOS_pbar_mvsPt",m,pt,alpha) ;
2511 FillHistogram("hMC_Resid_PHOS_other_mvsPt",m,pt,alpha) ;
2520 if(fEMCALEvent && fConvEvent){
2521 for(Int_t iEMCAL=0; iEMCAL<fEMCALEvent->GetEntriesFast();iEMCAL++){
2522 TLorentzVector * cal = static_cast<TLorentzVector*>(fEMCALEvent->At(iEMCAL)) ;
2523 Int_t iclu=fGammaEMCAL[iEMCAL] ;
2524 AliESDCaloCluster * clu = fESDEvent->GetCaloCluster(iclu);
2525 Int_t iprimEMCAL = clu->GetLabel() ; //# of particle hit EMCAL
2526 for(Int_t iConv = 0; iConv<fConvEvent->GetEntriesFast(); iConv++){
2527 TLorentzVector * cnv = static_cast<TLorentzVector*>(fConvEvent->At(iConv)) ;
2528 if(!cnv->TestBit(kConvOnFly) || cnv->TestBit(kConvEta))
2530 Int_t iv0=fGammaV0s[iConv] ;
2531 AliESDv0 * v0 = fESDEvent->GetV0(iv0) ;
2532 Int_t iprimNeg = TMath::Abs(fESDEvent->GetTrack(v0->GetNindex())->GetLabel()) ;
2533 Int_t iprimPos = TMath::Abs(fESDEvent->GetTrack(v0->GetPindex())->GetLabel()) ;
2535 //Check if there was a common ancistor
2536 Bool_t found = kFALSE ;
2537 Int_t curEMCAL=iprimEMCAL ;
2539 while(!found && curEMCAL>-1){
2540 Int_t curNeg=iprimNeg ;
2541 while(!found && curNeg>-1){
2542 if(curNeg==curEMCAL){
2547 curNeg=fStack->Particle(curNeg)->GetFirstMother() ;
2550 curEMCAL=fStack->Particle(curEMCAL)->GetFirstMother() ;
2553 curEMCAL=iprimEMCAL ;
2555 while(!found && curEMCAL>-1){
2556 Int_t curPos=iprimPos ;
2557 while(!found && curPos>-1){
2558 if(curPos==curEMCAL){
2563 curPos=fStack->Particle(curPos)->GetFirstMother() ;
2566 curEMCAL=fStack->Particle(curEMCAL)->GetFirstMother() ;
2569 if(commonA != commonB){
2571 AliInfo(Form("CommonA=%d, commonB=%d",commonA,commonB)) ;
2573 if(commonA>-1){//There was common particles
2574 Int_t pdg = fStack->Particle(commonA)->GetPdgCode() ;
2575 TLorentzVector pi=*cal + *cnv ;
2577 Double_t pt=pi.Pt() ;
2578 Double_t alpha=TMath::Abs(cal->Energy()-cnv->Energy())/(cal->Energy()+cnv->Energy()) ;
2582 case 22: //conversion
2583 FillHistogram("hMC_Resid_EMCAL_Phot_mvsPt",m,pt,alpha) ;
2586 FillHistogram("hMC_Resid_EMCAL_Pi0_mvsPt",m,pt,alpha) ;
2589 FillHistogram("hMC_Resid_EMCAL_eta_mvsPt",m,pt,alpha) ;
2595 FillHistogram("hMC_Resid_EMCAL_K_mvsPt",m,pt,alpha) ;
2599 FillHistogram("hMC_Resid_EMCAL_pi_mvsPt",m,pt,alpha) ;
2603 FillHistogram("hMC_Resid_EMCAL_pbar_mvsPt",m,pt,alpha) ;
2606 FillHistogram("hMC_Resid_EMCAL_other_mvsPt",m,pt,alpha) ;
2615 //------------- now photons ----------------
2616 for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
2617 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
2618 if(particle->GetPdgCode() != 22)
2621 if(particle->R() >rcut)
2624 if(TMath::Abs(particle->Eta())>0.9)
2627 Double_t pt = particle->Pt() ;
2628 //Total number of pi0 with creation radius <1 cm
2629 FillHistogram("hMC_CaloConv_phot",pt) ;
2632 Double_t x=0.,z=0. ;
2633 Bool_t hitPHOS = fPHOSgeom->ImpactOnEmc(particle, mod, z,x) ;
2634 Bool_t hitEMCAL= fEMCALgeom->Impact(particle) ;
2636 //Photons in PHOS/EMCAL acceptance
2638 FillHistogram("hMC_CaloConv_gammaPHOSacc",pt) ;
2640 FillHistogram("hMC_CaloConv_gammaEMCALacc",pt) ;
2642 //number of photons converted
2643 Bool_t converted = kFALSE ;
2644 if(particle->GetNDaughters()==2){
2645 TParticle * e1=fStack->Particle(particle->GetFirstDaughter()) ;
2646 TParticle * e2=fStack->Particle(particle->GetLastDaughter()) ;
2647 if(TMath::Abs(e1->GetPdgCode())==11 && TMath::Abs(e2->GetPdgCode())==11){ //conversion
2653 FillHistogram("hMC_CaloConv_gamma_conv",pt) ;
2655 //Converted photons with V0
2656 TLorentzVector pConv ;
2657 Bool_t foundV0=kFALSE ;
2658 for(Int_t iv0=0; iv0<fESDEvent->GetNumberOfV0s();iv0++){
2659 AliESDv0 * v0 = fESDEvent->GetV0(iv0) ;
2661 TParticle * negativeMC = fStack->Particle(TMath::Abs(fESDEvent->GetTrack(v0->GetNindex())->GetLabel()));
2662 TParticle * positiveMC = fStack->Particle(TMath::Abs(fESDEvent->GetTrack(v0->GetPindex())->GetLabel()));
2664 if(negativeMC->GetMother(0) != positiveMC->GetMother(0))
2667 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
2670 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
2674 TParticle * v0Gamma = fStack->Particle(negativeMC->GetMother(0));
2675 Bool_t same = (v0Gamma == particle) ;
2676 TParticle * tmp = v0Gamma ;
2677 while(!same && tmp->GetFirstMother()>=0){
2678 tmp = fStack->Particle(tmp->GetFirstMother());
2679 same = (tmp == particle) ;
2683 const AliExternalTrackParam * paramPos = v0->GetParamP() ;
2684 const AliExternalTrackParam * paramNeg = v0->GetParamN() ;
2685 AliKFParticle negKF(*paramNeg,11);
2686 AliKFParticle posKF(*paramPos,-11);
2687 pConv.SetXYZM(negKF.Px()+posKF.Px(),negKF.Py()+posKF.Py(),negKF.Pz()+negKF.Pz(),0.) ;
2692 FillHistogram("hMC_CaloConv_gamma_v0",pt) ;
2693 FillHistogram("hMC_CaloConv_gamma_v0_devsE",(particle->Energy()-pConv.E())/particle->Energy(),particle->Energy()) ;
2696 //Registered in PHOS/EMCAL
2697 Bool_t cluInPHOS = kFALSE,cluInEMCAL=kFALSE ;
2698 TLorentzVector pCalo ;
2699 Bool_t dist1=kFALSE, dist2=kFALSE ;
2700 for (Int_t i=0; i<fESDEvent->GetNumberOfCaloClusters(); i++) {
2701 AliESDCaloCluster * clu = fESDEvent->GetCaloCluster(i);
2702 Int_t iprim = clu->GetLabel() ; //# of particle hit PHOS/EMCAL
2703 Bool_t matched = kFALSE ;
2710 iprim=fStack->Particle(iprim)->GetFirstMother() ;
2715 if(clu->IsPHOS() && hitPHOS){
2717 clu->GetMomentum(pCalo ,vtx);
2718 if(clu->GetDistanceToBadChannel()<fBadDistCutPHOS)
2720 if(clu->GetDistanceToBadChannel()<2.*fBadDistCutPHOS)
2724 if(!clu->IsPHOS() && hitEMCAL){
2726 clu->GetMomentum(pCalo ,vtx);
2727 if(clu->GetDistanceToBadChannel()<fBadDistCutEMCAL)
2729 if(clu->GetDistanceToBadChannel()<2.*fBadDistCutEMCAL)
2736 FillHistogram("hMC_CaloConv_gamma_PHOSclu",pt) ;
2738 FillHistogram("hMC_CaloConv_gamma_PHOSclu_dist1",pt) ;
2740 FillHistogram("hMC_CaloConv_gamma_PHOSclu_dist2",pt) ;
2741 FillHistogram("hMC_CaloConv_gamma_PHOSclu_recE",pCalo.E()) ;
2742 FillHistogram("hMC_CaloConv_gamma_PHOSclu_devsE",(particle->Energy()-pCalo.E())/particle->Energy(),particle->Energy()) ;
2745 FillHistogram("hMC_CaloConv_gamma_EMCALclu",pt) ;
2747 FillHistogram("hMC_CaloConv_gamma_EMCALclu_dist1",pt) ;
2749 FillHistogram("hMC_CaloConv_gamma_EMCALclu_dist2",pt) ;
2750 FillHistogram("hMC_CaloConv_gamma_EMCALclu_recE",pCalo.E()) ;
2751 FillHistogram("hMC_CaloConv_gamma_EMCALclu_devsE",(particle->Energy()-pCalo.E())/particle->Energy(),particle->Energy()) ;
2755 //_____________________________________________________________________________
2756 void AliAnalysisTaskCaloConv::FillHistogram(const char * key,Double_t x)const{
2758 TH1F * tmp = dynamic_cast<TH1F*>(fOutputContainer->FindObject(key)) ;
2760 AliInfo(Form("can not find histogram <%s> ",key)) ;
2765 //_____________________________________________________________________________
2766 void AliAnalysisTaskCaloConv::FillHistogram(const char * key,Double_t x,Double_t y)const{
2768 TObject * tmp = fOutputContainer->FindObject(key) ;
2770 AliInfo(Form("can not find histogram <%s> ",key)) ;
2773 if(tmp->IsA() == TClass::GetClass("TH1F")){
2774 ((TH1F*)tmp)->Fill(x,y) ;
2777 if(tmp->IsA() == TClass::GetClass("TH2F")){
2778 ((TH2F*)tmp)->Fill(x,y) ;
2781 AliError(Form("Calling FillHistogram with 2 parameters for histo <%s> of type %s",key,tmp->IsA()->GetName())) ;
2784 //_____________________________________________________________________________
2785 void AliAnalysisTaskCaloConv::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
2786 //Fills 1D histograms with key
2787 TObject * tmp = fOutputContainer->FindObject(key) ;
2789 AliInfo(Form("can not find histogram <%s> ",key)) ;
2792 if(tmp->IsA() == TClass::GetClass("TH2F")){
2793 ((TH2F*)tmp)->Fill(x,y,z) ;
2796 if(tmp->IsA() == TClass::GetClass("TH3F")){
2797 ((TH3F*)tmp)->Fill(x,y,z) ;
2801 //______________________________________________________________________________
2802 Double_t AliAnalysisTaskCaloConv::PlanarityAngle(const AliExternalTrackParam * paramPos,const AliExternalTrackParam * paramNeg)const {
2803 //calculate angle between e+e- plain and perpendicular to MF
2804 //We need sign of MagField to calculate orienation
2806 TVector3 u(paramPos->Px()+paramNeg->Px(),paramPos->Py()+paramNeg->Py(),paramPos->Pz()+paramNeg->Pz()) ;
2808 TVector3 vPos(paramPos->Px(),paramPos->Py(),paramPos->Pz()) ;
2809 TVector3 vNeg(paramNeg->Px(),paramNeg->Py(),paramNeg->Pz()) ;
2810 TVector3 v=vPos.Cross(vNeg) ;
2811 TVector3 w = u.Cross(v);
2813 TVector3 ua=u.Cross(z);
2814 Double_t wa = w.Angle(ua);
2815 Double_t mfield=fESDEvent->GetMagneticField() ;
2817 return wa; //normal field
2819 return TMath::Pi()-wa ; //reverse field
2822 //______________________________________________________________________________
2823 Bool_t AliAnalysisTaskCaloConv::IsGoodChannel(const char * det, Int_t mod, Int_t ix, Int_t iz){
2824 //Check if this channel belogs to the good ones
2826 if(strcmp(det,"PHOS")==0){
2828 AliError(Form("No bad map for PHOS module %d ",mod)) ;
2831 if(!fPHOSBadMap[mod]){
2832 AliError(Form("No Bad map for PHOS module %d",mod)) ;
2835 if(fPHOSBadMap[mod]->GetBinContent(ix,iz)>0)
2841 if(strcmp(det,"EMCAL")==0){
2843 AliError(Form("No bad map for EMCAL module %d ",mod)) ;
2846 if(!fEMCALBadMap[mod]){
2847 AliError(Form("No bad map for EMCAL module %d ",mod)) ;
2850 if(fEMCALBadMap[mod]->GetBinContent(ix,iz)>0)
2856 AliError(Form("Can not find bad channels for detector %s ",det)) ;
2862 //______________________________________________________________________________
2863 void AliAnalysisTaskCaloConv::Recalibrate(Double_t &m, Double_t &pt, const TLorentzVector *calo, const TLorentzVector * conv, Int_t iw, Int_t in) {
2864 //Apply decalibration and non-linearity
2865 TLorentzVector calo2(*calo) ;
2866 Double_t en=calo2.E() ;
2868 Double_t sigma=0.06+0.005*iw ; //additional smearing
2870 Double_t a=0.02*(in%6-2.5) ;
2871 Double_t b=0.5+1.*((Int_t)in/6) ;
2872 Double_t enNew=1.-a*TMath::Exp(-en/b) ;
2873 Double_t corr=gRandom->Gaus(enNew,sigma) ;
2876 m=(calo2+ *conv).M() ;
2877 pt=(calo2+ *conv).Pt() ;
2880 //______________________________________________________________________________
2881 void AliAnalysisTaskCaloConv::RecalibrateEMCAL(Double_t &m, Double_t &pt, const TLorentzVector *calo, const TLorentzVector * conv, Int_t iw, Int_t in) {
2882 //Apply decalibration and non-linearity
2883 TLorentzVector calo2(*calo) ;
2884 Double_t en=calo2.E() ;
2886 Double_t sigma=0.04+0.005*iw ; //additional smearing
2888 Double_t a=0.02*(in%6-2.5) ;
2889 Double_t b=0.25+0.5*((Int_t)in/6) ;
2890 Double_t enNew=1.-a*TMath::Exp(-en/b) ;
2891 Double_t corr=gRandom->Gaus(enNew,sigma) ;
2894 m=(calo2+ *conv).M() ;
2895 pt=(calo2+ *conv).Pt() ;
2898 //______________________________________________________________________________
2899 void AliAnalysisTaskCaloConv::RecalibrateConvPHOS(Double_t &m, Double_t &pt, const TLorentzVector *calo, const TLorentzVector * conv, Int_t iw, Int_t in) {
2900 //Apply decalibration and non-linearity
2902 //First default PHOS smearing
2903 TLorentzVector calo2(*calo) ;
2904 Double_t en=calo2.E() ;
2906 Double_t sigma=0.065 ; //additional smearing
2910 Double_t enNew=1.+a*TMath::Exp(-en/b) ;
2911 Double_t corr=gRandom->Gaus(enNew,sigma) ;
2914 //Now conversion photon
2915 TLorentzVector conv2(*conv) ;
2916 //linear offset in z:
2917 Double_t eta=conv2.Eta() ;
2918 Double_t c=1.e-3*iw ;
2919 eta+= c *TMath::Sign(0.9-TMath::Abs(eta),eta) ;
2921 //Smear energy and add nonlinearity
2923 Double_t enConv=conv2.E() ;
2924 Double_t ac=0.02*(in%5) ;
2925 Double_t bc=0.25+0.5*((Int_t)in/5) ;
2926 Double_t enNewc=1.+ac*TMath::Exp(-enConv/bc) ;
2927 corr=gRandom->Gaus(enNewc,0.01) ;
2928 Double_t ptc=conv2.Pt()*corr ;
2929 conv2.SetPtEtaPhiM(ptc,eta,conv2.Phi(),0.) ;
2931 m =(calo2 + conv2).M() ;
2932 pt=(calo2 + conv2).Pt() ;
2935 //______________________________________________________________________________
2936 void AliAnalysisTaskCaloConv::GetArmenterosQtAlfa(AliKFParticle* positiveKFParticle, AliKFParticle * negativeKFParticle, AliKFParticle * gammaKFCandidate, Double_t armenterosQtAlfa[2] ){
2937 //see header file for documentation
2939 TVector3 momentumVectorPositiveKF(positiveKFParticle->GetPx(),positiveKFParticle->GetPy(),positiveKFParticle->GetPz());
2940 TVector3 momentumVectorNegativeKF(negativeKFParticle->GetPx(),negativeKFParticle->GetPy(),negativeKFParticle->GetPz());
2941 TVector3 vecV0(gammaKFCandidate->GetPx(),gammaKFCandidate->GetPy(),gammaKFCandidate->GetPz());
2943 Float_t thetaV0pos=TMath::ACos(( momentumVectorPositiveKF* vecV0)/(momentumVectorPositiveKF.Mag() * vecV0.Mag()));
2944 Float_t thetaV0neg=TMath::ACos(( momentumVectorNegativeKF* vecV0)/(momentumVectorNegativeKF.Mag() * vecV0.Mag()));
2946 Float_t alfa =((momentumVectorPositiveKF.Mag())*TMath::Cos(thetaV0pos)-(momentumVectorNegativeKF.Mag())*TMath::Cos(thetaV0neg))/
2947 ((momentumVectorPositiveKF.Mag())*TMath::Cos(thetaV0pos)+(momentumVectorNegativeKF.Mag())*TMath::Cos(thetaV0neg)) ;
2950 Float_t qt = momentumVectorPositiveKF.Mag()*TMath::Sin(thetaV0pos);
2952 armenterosQtAlfa[0]=qt;
2953 armenterosQtAlfa[1]=alfa;