--- /dev/null
+# -*- mode: cmake -*-
+
+set(SRCS
+ GammaConv/AliV0Reader.cxx GammaConv/AliAnalysisTaskGammaConversion.cxx GammaConv/AliGammaConversionHistograms.cxx
+)
+
+# fill list of header files from list of source files
+# by exchanging the file extension
+String(REPLACE ".cxx" ".h" HDRS "${SRCS}")
+
+AddLibrary(PWG4PartCorr "${SRCS}" "${HDRS}")
+
--- /dev/null
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *\r
+ * Version 1.0 *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+\r
+////////////////////////////////////////////////\r
+//--------------------------------------------- \r
+// Class used to do analysis on conversion pairs\r
+//---------------------------------------------\r
+////////////////////////////////////////////////\r
+\r
+// root\r
+#include <TChain.h>\r
+\r
+// analysis\r
+#include "AliAnalysisTaskGammaConversion.h"\r
+#include "AliStack.h"\r
+#include "AliLog.h"\r
+#include <vector>\r
+\r
+class AliKFVertex;\r
+class AliAODHandler;\r
+class AliAODEvent;\r
+class ALiESDEvent;\r
+class AliMCEvent;\r
+class AliMCEventHandler;\r
+class AliESDInputHandler;\r
+class AliAnalysisManager;\r
+class Riostream;\r
+class TFile;\r
+class TInterpreter;\r
+class TSystem;\r
+class TROOT;\r
+\r
+ClassImp(AliAnalysisTaskGammaConversion)\r
+\r
+\r
+AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion():\r
+ AliAnalysisTaskSE(),\r
+ fV0Reader(NULL),\r
+ fStack(NULL),\r
+ fOutputContainer(NULL),\r
+ fHistograms(NULL),\r
+ fDoMCTruth(kFALSE),\r
+ fMCAllGammas(),\r
+ fMCPi0s(),\r
+ fMCEtas(),\r
+ fMCGammaChic(),\r
+ fKFReconstructedGammas(),\r
+ fElectronMass(-1),\r
+ fGammaMass(-1),\r
+ fPi0Mass(-1),\r
+ fEtaMass(-1),\r
+ fGammaWidth(-1),\r
+ fPi0Width(-1),\r
+ fEtaWidth(-1),\r
+ fCalculateBackground(kFALSE)\r
+{\r
+ // Default constructor\r
+ // Common I/O in slot 0\r
+ DefineInput (0, TChain::Class());\r
+ DefineOutput(0, TTree::Class());\r
+\r
+ // Your private output\r
+ DefineOutput(1, TList::Class());\r
+}\r
+\r
+AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name):\r
+ AliAnalysisTaskSE(name),\r
+ fV0Reader(NULL),\r
+ fStack(NULL),\r
+ fOutputContainer(0x0),\r
+ fHistograms(NULL),\r
+ fDoMCTruth(kFALSE),\r
+ fMCAllGammas(),\r
+ fMCPi0s(),\r
+ fMCEtas(),\r
+ fMCGammaChic(),\r
+ fKFReconstructedGammas(),\r
+ fElectronMass(-1),\r
+ fGammaMass(-1),\r
+ fPi0Mass(-1),\r
+ fEtaMass(-1),\r
+ fGammaWidth(-1),\r
+ fPi0Width(-1),\r
+ fEtaWidth(-1),\r
+ fCalculateBackground(kFALSE)\r
+{\r
+ // Common I/O in slot 0\r
+ DefineInput (0, TChain::Class());\r
+ DefineOutput(0, TTree::Class());\r
+ \r
+ // Your private output\r
+ DefineOutput(1, TList::Class());\r
+}\r
+\r
+AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion() \r
+{\r
+ // Remove all pointers\r
+ \r
+ if(fOutputContainer){\r
+ fOutputContainer->Clear() ; \r
+ delete fOutputContainer ;\r
+ }\r
+ if(fHistograms){\r
+ delete fHistograms;\r
+ }\r
+ if(fV0Reader){\r
+ delete fV0Reader;\r
+ }\r
+}\r
+\r
+\r
+void AliAnalysisTaskGammaConversion::Init()\r
+{\r
+ // Initialization\r
+ AliLog::SetGlobalLogLevel(AliLog::kError);\r
+}\r
+\r
+\r
+void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/)\r
+{\r
+ // Execute analysis for current event\r
+ \r
+ ConnectInputData("");\r
+ \r
+ //clear vectors\r
+ fMCAllGammas.clear();\r
+ fMCPi0s.clear();\r
+ fMCEtas.clear();\r
+ fMCGammaChic.clear();\r
+\r
+ fKFReconstructedGammas.clear();\r
+\r
+ //Clear the data in the v0Reader\r
+ fV0Reader->UpdateEventByEventData();\r
+\r
+ // Process the MC information\r
+ if(fDoMCTruth){\r
+ ProcessMCData();\r
+ }\r
+\r
+ // Process the v0 information\r
+ ProcessV0s();\r
+\r
+ //calculate background if flag is set\r
+ if(fCalculateBackground){\r
+ CalculateBackground();\r
+ }\r
+\r
+ // Process reconstructed gammas\r
+ ProcessGammasForNeutralMesonAnalysis();\r
+\r
+ PostData(1, fOutputContainer);\r
+ \r
+}\r
+\r
+void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t */*option*/){\r
+ // see header file for documentation\r
+\r
+ if(fV0Reader == NULL){\r
+ // Write warning here cuts and so on are default if this ever happens\r
+ }\r
+ fV0Reader->Initialize();\r
+}\r
+\r
+void AliAnalysisTaskGammaConversion::ProcessMCData(){\r
+ // see header file for documentation\r
+ \r
+ fStack = fV0Reader->GetMCStack();\r
+\r
+ for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {\r
+ TParticle* particle = (TParticle *)fStack->Particle(iTracks);\r
+\r
+ if (!particle) {\r
+ //print warning here\r
+ continue;\r
+ }\r
+ \r
+ if(particle->Pt()<fV0Reader->GetPtCut()){\r
+ continue;\r
+ }\r
+ \r
+ if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut()){\r
+ continue;\r
+ }\r
+\r
+ if(particle->R()>fV0Reader->GetMaxRCut()){ // cuts on distance from collision point\r
+ continue;\r
+ }\r
+\r
+ Double_t tmpPhi=particle->Phi();\r
+ if(particle->Phi()> TMath::Pi()){\r
+ tmpPhi = particle->Phi()-(2*TMath::Pi());\r
+ }\r
+\r
+ \r
+ //process the gammas\r
+ if (particle->GetPdgCode()== 22){\r
+ fMCAllGammas.push_back(particle);\r
+ if(particle->GetMother(0)>-1){ //Means we have a mother\r
+ if( fStack->Particle(particle->GetMother(0))->GetPdgCode() != 22 ){//Checks for a non gamma mother.\r
+ fHistograms->FillHistogram("MC_Gamma_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_Gamma_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_Gamma_Eta", particle->Eta());\r
+ \r
+ fHistograms->FillHistogram("MC_Gamma_Phi", tmpPhi);\r
+\r
+ //adding the conversion points from all gammas with e+e- daughters\r
+ if(particle->GetNDaughters() >= 2){\r
+ TParticle* daughter0 = NULL;\r
+ TParticle* daughter1 = NULL;\r
+ \r
+ for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){\r
+ TParticle *tmpDaughter = fStack->Particle(daughterIndex);\r
+ if(tmpDaughter->GetUniqueID() == 5){\r
+ if(tmpDaughter->GetPdgCode() == 11){\r
+ daughter0 = tmpDaughter;\r
+ }\r
+ else if(tmpDaughter->GetPdgCode() == -11){\r
+ daughter1 = tmpDaughter;\r
+ }\r
+ }\r
+ }\r
+\r
+ if(daughter0 == NULL || daughter1 == NULL){ // means we do not have two daughters from pair production\r
+ continue;\r
+ }\r
+\r
+ if(daughter0->R()>fV0Reader->GetMaxRCut() || daughter1->R()>fV0Reader->GetMaxRCut()){\r
+ continue;\r
+ }\r
+ \r
+ if((daughter0->GetPdgCode() == -11 && daughter1->GetPdgCode()) == 11 ||\r
+ (daughter0->GetPdgCode() == 11 && daughter1->GetPdgCode() == -11)){\r
+\r
+ // begin Mapping \r
+ Int_t rBin = fHistograms->GetRBin(daughter0->R());\r
+ Int_t phiBin = fHistograms->GetPhiBin(daughter0->Phi());\r
+ \r
+ TString nameMCMappingPhiR="";\r
+ nameMCMappingPhiR.Form("MC_EP_Mapping-Phi%02d-R%02d",phiBin,rBin);\r
+ fHistograms->FillHistogram(nameMCMappingPhiR, daughter0->Vz(), particle->Eta());\r
+ \r
+ TString nameMCMappingPhi="";\r
+ nameMCMappingPhi.Form("MC_EP_Mapping-Phi%02d",phiBin);\r
+ fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());\r
+ \r
+ TString nameMCMappingR="";\r
+ nameMCMappingR.Form("MC_EP_Mapping-R%02d",rBin);\r
+ fHistograms->FillHistogram(nameMCMappingR, particle->Eta());\r
+ \r
+ TString nameMCMappingPhiInR="";\r
+ nameMCMappingPhiInR.Form("MC_EP_Mapping_Phi_R-%02d",rBin);\r
+ fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);\r
+ //end mapping\r
+\r
+ fHistograms->FillHistogram("MC_EP_R",daughter0->R());\r
+ fHistograms->FillHistogram("MC_EP_ZR",daughter0->Vz(),daughter0->R());\r
+ fHistograms->FillHistogram("MC_EP_XY",daughter0->Vx(),daughter0->Vy());\r
+ fHistograms->FillHistogram("MC_EP_OpeningAngle",GetMCOpeningAngle(daughter0, daughter1));\r
+ }// end if((daughter0->GetPdgCode() == -11 && daughter1->GetPdgCode()) == 11 ||....... approx 20 lines above\r
+ }// end if(particle->GetNDaughters() >= 2){\r
+ } // end if( fStack->Particle(particle->GetMother(0))->GetPdgCode() != 22 )\r
+ if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0 \r
+ fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S\r
+ fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2\r
+ ){ \r
+ fMCGammaChic.push_back(particle);\r
+ }\r
+ }// end if(particle->GetMother(0)>-1)\r
+ else{//means we have a primary particle\r
+ fHistograms->FillHistogram("MC_DirectGamma_Energy",particle->Energy());\r
+ fHistograms->FillHistogram("MC_DirectGamma_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_DirectGamma_Eta", particle->Eta());\r
+ fHistograms->FillHistogram("MCDirectGammaPhi", tmpPhi);\r
+\r
+ //adding the conversion points from all gammas with e+e- daughters\r
+ if(particle->GetNDaughters() == 2){\r
+ TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());\r
+ TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());\r
+ if((daughter0->GetPdgCode() == -11 && daughter1->GetPdgCode() == 11) ||\r
+ (daughter0->GetPdgCode() == 11 && daughter1->GetPdgCode() == -11)){\r
+ \r
+ fHistograms->FillHistogram("MC_EP_R",daughter0->R());\r
+ fHistograms->FillHistogram("MC_EP_ZR",daughter0->Vz(),daughter0->R());\r
+ fHistograms->FillHistogram("MC_EP_XY",daughter0->Vx(),daughter0->Vy());\r
+ fHistograms->FillHistogram("MC_EP_OpeningAngle",GetMCOpeningAngle(daughter0, daughter1));\r
+\r
+ }\r
+ }\r
+ }// end else\r
+ }// end if (particle->GetPdgCode()== 22){\r
+ else if (TMath::Abs(particle->GetPdgCode())== 11){ // Means we have an electron or a positron\r
+ if(particle->GetMother(0)>-1){ // means we have a mother\r
+ if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==22 ){ // Means we have a gamma mother\r
+ if(particle->GetPdgCode() == 11){//electron \r
+ fHistograms->FillHistogram("MC_E_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_E_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_E_Eta", particle->Eta());\r
+ fHistograms->FillHistogram("MC_E_Phi", tmpPhi);\r
+ }\r
+ if(particle->GetPdgCode() == -11){//positron \r
+ fHistograms->FillHistogram("MC_P_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_P_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_P_Eta", particle->Eta());\r
+ fHistograms->FillHistogram("MC_P_Phi", tmpPhi);\r
+ }\r
+ }\r
+ }\r
+ } // end else if (TMath::Abs(particle->GetPdgCode())== 11)\r
+ else if(particle->GetNDaughters() == 2){\r
+\r
+ TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());\r
+ TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());\r
+ if(daughter0->GetPdgCode() == 22 && daughter1->GetPdgCode() == 22){//check for gamma gamma daughters\r
+ \r
+ if(particle->GetPdgCode()==111){//Pi0\r
+ \r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);\r
+\r
+ if( iTracks >= fStack->GetNprimary()){\r
+ \r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());\r
+\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_R", particle->R());\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_ZR", particle->Vz(),particle->R());\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
+ fHistograms->FillHistogram("MC_Pi0_Secondaries_XY", particle->Vx(),particle->Vy());//only fill from one daughter to avoid multiple filling\r
+ }\r
+ else{\r
+ fHistograms->FillHistogram("MC_Pi0_Eta", particle->Eta());\r
+\r
+ fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);\r
+ fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_Pi0_R", particle->R());\r
+ fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());\r
+ fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
+ fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling\r
+ }\r
+ }\r
+ else if(particle->GetPdgCode()==221){//Eta\r
+ fHistograms->FillHistogram("MC_Eta_Eta", particle->Eta());\r
+\r
+ fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);\r
+ fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_Eta_R", particle->R());\r
+ fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());\r
+ fHistograms->FillHistogram("MC_Eta_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
+ fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling\r
+ }\r
+ \r
+ //the match data should be filled no matter which mother the gamma-gamma comes from\r
+ fHistograms->FillHistogram("MC_Match_Gamma_R", particle->R());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_ZR", particle->Vz(),particle->R());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_XY", particle->Vx(),particle->Vy());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_Mass", particle->GetCalcMass());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
+ fHistograms->FillHistogram("MC_Match_Gamma_Energy", particle->Energy());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_Pt", particle->Pt());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_Eta", particle->Eta());\r
+ fHistograms->FillHistogram("MC_Match_Gamma_Phi",tmpPhi);\r
+ }\r
+ }// end else if(particle->GetNDaughters() == 2)\r
+ }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++)\r
+} // end ProcessMCData\r
+\r
+void AliAnalysisTaskGammaConversion::ProcessV0s(){\r
+ // see header file for documentation\r
+\r
+ Int_t nSurvivingV0s=0;\r
+ while(fV0Reader->NextV0()){\r
+ nSurvivingV0s++;\r
+ //-------------------------- filling v0 information -------------------------------------\r
+ fHistograms->FillHistogram("ESD_EP_OpeningAngle", fV0Reader->GetOpeningAngle()); \r
+ fHistograms->FillHistogram("ESD_EP_R", fV0Reader->GetXYRadius());\r
+ fHistograms->FillHistogram("ESD_EP_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());\r
+ fHistograms->FillHistogram("ESD_EP_XY", fV0Reader->GetX(),fV0Reader->GetY());\r
+ \r
+ \r
+ fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());\r
+ fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());\r
+ fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());\r
+ fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());\r
+ \r
+ fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());\r
+ fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());\r
+ fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());\r
+ fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());\r
+ \r
+ fHistograms->FillHistogram("ESD_Gamma_Energy", fV0Reader->GetMotherCandidateEnergy());\r
+ fHistograms->FillHistogram("ESD_Gamma_Pt", fV0Reader->GetMotherCandidatePt());\r
+ fHistograms->FillHistogram("ESD_Gamma_Eta", fV0Reader->GetMotherCandidateEta());\r
+ fHistograms->FillHistogram("ESD_Gamma_Phi", fV0Reader->GetMotherCandidatePhi());\r
+\r
+\r
+ // begin mapping\r
+ Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());\r
+ Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());\r
+ Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();\r
+ \r
+ TString nameESDMappingPhiR="";\r
+ nameESDMappingPhiR.Form("ESD_EP_Mapping-Phi%02d-R%02d",phiBin,rBin);\r
+ fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);\r
+\r
+ TString nameESDMappingPhi="";\r
+ nameESDMappingPhi.Form("ESD_EP_Mapping-Phi%02d",phiBin);\r
+ fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);\r
+\r
+ TString nameESDMappingR="";\r
+ nameESDMappingR.Form("ESD_EP_Mapping-R%02d",rBin);\r
+ fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta); \r
+\r
+ TString nameESDMappingPhiInR="";\r
+ nameESDMappingPhiInR.Form("ESD_EP_Mapping_Phi_R-%02d",rBin);\r
+ fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());\r
+ // end mapping\r
+ \r
+ fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());\r
+\r
+ //----------------------------------- checking for "real" conversions (MC match) --------------------------------------\r
+ if(fDoMCTruth){\r
+ if(fV0Reader->HasSameMCMother() == kFALSE){\r
+ continue;\r
+ }\r
+ \r
+ TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();\r
+ TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();\r
+\r
+ if(negativeMC->GetPdgCode()!=11 || positiveMC->GetPdgCode()!=-11){\r
+ continue;\r
+ }\r
+ if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_XY", fV0Reader->GetX(),fV0Reader->GetY());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_OpeningAngle", fV0Reader->GetOpeningAngle());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Pt", fV0Reader->GetMotherCandidatePt());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Energy", fV0Reader->GetMotherCandidateEnergy());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Eta", fV0Reader->GetMotherCandidateEta());\r
+\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Phi", fV0Reader->GetMotherCandidatePhi());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Mass", fV0Reader->GetMotherCandidateMass());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Width", fV0Reader->GetMotherCandidateWidth());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_Chi2", fV0Reader->GetMotherCandidateChi2());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_NDF", fV0Reader->GetMotherCandidateNDF());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_R", fV0Reader->GetXYRadius());\r
+ fHistograms->FillHistogram("ESD_Match_Gamma_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());\r
+\r
+ //resolution\r
+ Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();\r
+ Double_t esdpt = fV0Reader->GetMotherCandidatePt();\r
+ Double_t resdPt = 0;\r
+ if(mcpt != 0){\r
+ resdPt = ((esdpt - mcpt)/mcpt)*100;\r
+ }\r
+\r
+ fHistograms->FillHistogram("Resolution_dPt", mcpt, resdPt);\r
+ fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);\r
+ fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);\r
+ \r
+ Double_t resdZ = 0;\r
+ if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){\r
+ resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100;\r
+ }\r
+ \r
+ fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);\r
+ fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());\r
+ fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());\r
+ \r
+ Double_t resdR = 0;\r
+ if(fV0Reader->GetNegativeMCParticle()->R() != 0){\r
+ resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100;\r
+ }\r
+ fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);\r
+ fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());\r
+ fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());\r
+ fHistograms->FillHistogram("Resolution_dR_dPt", resdR, resdPt);\r
+ }\r
+ }\r
+ }\r
+ fHistograms->FillHistogram("NumberOfSurvivingV0s", nSurvivingV0s);\r
+ fHistograms->FillHistogram("NumberOfV0s", fV0Reader->GetNumberOfV0s());\r
+}\r
+\r
+void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){\r
+ // see header file for documentation\r
+\r
+ for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){\r
+ for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){\r
+ AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];\r
+ AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];\r
+\r
+ \r
+ AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);\r
+\r
+ Double_t massTwoGammaCandidate =0.;\r
+ Double_t widthTwoGammaCandidate = 0.;\r
+ Double_t chi2TwoGammaCandidate =10000.; \r
+ twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);\r
+ if(twoGammaCandidate->GetNDF()>0){\r
+ chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();\r
+ if(chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()){\r
+\r
+ TVector3 vectorTwoGammaCandidate(twoGammaCandidate->Px(),twoGammaCandidate->Py(),twoGammaCandidate->Pz());\r
+\r
+ Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);\r
+\r
+ //Calculating by hand the radius\r
+ Double_t tmpX= twoGammaCandidate->GetX();\r
+ Double_t tmpY= twoGammaCandidate->GetY();\r
+ \r
+ Double_t radiusTwoGammaCandidate = TMath::Sqrt(tmpX*tmpX + tmpY*tmpY);\r
+\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_Energy", twoGammaCandidate->GetE());\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_Pt", sqrt(twoGammaCandidate->GetPx()*twoGammaCandidate->GetPx()+twoGammaCandidate->GetPy()*twoGammaCandidate->GetPy()));\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_Eta", vectorTwoGammaCandidate.Eta());\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_Phi", vectorTwoGammaCandidate.Phi());\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_Mass", massTwoGammaCandidate);\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_R", radiusTwoGammaCandidate);\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_ZR", tmpY, radiusTwoGammaCandidate);\r
+ fHistograms->FillHistogram("ESD_TwoGammaCombination_XY", tmpX, tmpY);\r
+ fHistograms->FillHistogram("InvMass_vs_Pt_Spectra",massTwoGammaCandidate ,sqrt(twoGammaCandidate->GetPx()*twoGammaCandidate->GetPx()+twoGammaCandidate->GetPy()*twoGammaCandidate->GetPy()));\r
+ }\r
+ }\r
+ delete twoGammaCandidate;\r
+ }\r
+ }\r
+}\r
+\r
+void AliAnalysisTaskGammaConversion::CalculateBackground(){\r
+ // see header file for documentation\r
+\r
+ vector<AliKFParticle> vectorCurrentEventGoodV0s = fV0Reader->GetCurrentEventGoodV0s();\r
+ vector<AliKFParticle> vectorPreviousEventGoodV0s = fV0Reader->GetPreviousEventGoodV0s();\r
+ for(UInt_t iCurrent=0;iCurrent<vectorCurrentEventGoodV0s.size();iCurrent++){\r
+ AliKFParticle * currentEventGoodV0 = &vectorCurrentEventGoodV0s.at(iCurrent);\r
+ for(UInt_t iPrevious=0;iPrevious<vectorPreviousEventGoodV0s.size();iPrevious++){\r
+ AliKFParticle * previousGoodV0 = &vectorPreviousEventGoodV0s.at(iPrevious);\r
+\r
+ AliKFParticle *backgroundCandidate = new AliKFParticle(*currentEventGoodV0,*previousGoodV0);\r
+\r
+ Double_t massBG =0.;\r
+ Double_t widthBG = 0.;\r
+ Double_t chi2BG =10000.; \r
+ backgroundCandidate->GetMass(massBG,widthBG);\r
+ if(backgroundCandidate->GetNDF()>0){\r
+ chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();\r
+ if(chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()){\r
+\r
+ TVector3 vectorBGCandidate(backgroundCandidate->Px(),backgroundCandidate->Py(),backgroundCandidate->Pz());\r
+\r
+ Double_t openingAngleBG = currentEventGoodV0->GetAngle(*previousGoodV0);\r
+\r
+ //Calculating by hand the radius (find a better way)\r
+ Double_t tmpX= backgroundCandidate->GetX();\r
+ Double_t tmpY= backgroundCandidate->GetY();\r
+ \r
+ Double_t radiusBG = TMath::Sqrt(tmpX*tmpX + tmpY*tmpY);\r
+\r
+ fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);\r
+ fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());\r
+ fHistograms->FillHistogram("ESD_Background_Pt", sqrt(backgroundCandidate->GetPx()*backgroundCandidate->GetPx()+backgroundCandidate->GetPy()*backgroundCandidate->GetPy()));\r
+ fHistograms->FillHistogram("ESD_Background_Eta", vectorBGCandidate.Eta());\r
+ fHistograms->FillHistogram("ESD_Background_Phi", vectorBGCandidate.Phi());\r
+ fHistograms->FillHistogram("ESD_Background_Mass", massBG);\r
+ fHistograms->FillHistogram("ESD_Background_R", radiusBG);\r
+ fHistograms->FillHistogram("ESD_Background_ZR", tmpY, radiusBG);\r
+ fHistograms->FillHistogram("ESD_Background_XY", tmpX, tmpY);\r
+ fHistograms->FillHistogram("Background_InvMass_vs_Pt_Spectra",massBG,sqrt(backgroundCandidate->GetPx()*backgroundCandidate->GetPx()+backgroundCandidate->GetPy()*backgroundCandidate->GetPy()));\r
+ }\r
+ }\r
+ delete backgroundCandidate; \r
+ }\r
+ }\r
+}\r
+\r
+void AliAnalysisTaskGammaConversion::Terminate(Option_t */*option*/)\r
+{\r
+ // Terminate analysis\r
+ //\r
+ AliDebug(1,"Do nothing in Terminate");\r
+}\r
+\r
+void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()\r
+{\r
+ // Create the output container\r
+ if(fOutputContainer != NULL){\r
+ delete fOutputContainer;\r
+ fOutputContainer = NULL;\r
+ }\r
+ if(fOutputContainer == NULL){\r
+ fOutputContainer = new TList();\r
+ }\r
+ fHistograms->GetOutputContainer(fOutputContainer);\r
+ fOutputContainer->SetName(GetName());\r
+}\r
+\r
+Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(TParticle* daughter0, TParticle* daughter1) const{\r
+ //helper function\r
+ TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());\r
+ TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());\r
+ return v3D0.Angle(v3D1);\r
+}\r
--- /dev/null
+#ifndef ALIANALYSISTASKGAMMACONVERSION_H\r
+#define ALIANALYSISTASKGAMMACONVERSION_H\r
+ \r
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice */\r
+\r
+////////////////////////////////////////////////\r
+//--------------------------------------------- \r
+// Class used to do analysis on conversion pairs\r
+//---------------------------------------------\r
+////////////////////////////////////////////////\r
+ \r
+#include "AliAnalysisTaskSE.h"\r
+#include <vector>\r
+#include "AliV0Reader.h"\r
+\r
+class AliGammaConversionHistograms;\r
+class AliESDv0;\r
+class AliKFParticle;\r
+class AliESDInputHandler;\r
+class AliESDEvent;\r
+class AliAODEvent;\r
+class TList;\r
+class AliStack;\r
+\r
+class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE\r
+{\r
+ public:\r
+ AliAnalysisTaskGammaConversion();\r
+ AliAnalysisTaskGammaConversion(const char* name);\r
+ virtual ~AliAnalysisTaskGammaConversion() ;// virtual destructor\r
+ \r
+ // Implementation of interface methods\r
+ virtual void UserCreateOutputObjects();\r
+ virtual void Init();\r
+ virtual void LocalInit() {Init();}\r
+ virtual void Exec(Option_t *option);\r
+ virtual void Terminate(Option_t *option);\r
+ virtual void ConnectInputData(Option_t *);\r
+ \r
+ void ProcessMCData();\r
+ void ProcessV0s();\r
+ void ProcessGammasForNeutralMesonAnalysis();\r
+ void SetHistograms(AliGammaConversionHistograms *histograms){fHistograms=histograms;}\r
+ void SetDoMCTruth(Bool_t flag){fDoMCTruth=flag;}\r
+ void SetElectronMass(Double_t electronMass){fElectronMass = electronMass;}\r
+ void SetGammaMass(Double_t gammaMass){fGammaMass = gammaMass;}\r
+ void SetGammaWidth(Double_t gammaWidth){fGammaWidth = gammaWidth;}\r
+ void SetPi0Mass(Double_t pi0Mass){fPi0Mass = pi0Mass;}\r
+ void SetPi0Width(Double_t pi0Width){fPi0Width = pi0Width;}\r
+ void SetEtaMass(Double_t etaMass){fEtaMass = etaMass;}\r
+ void SetEtaWidth(Double_t etaWidth){fEtaWidth = etaWidth;}\r
+ void SetV0Reader(AliV0Reader* reader){fV0Reader=reader;}\r
+ void SetCalculateBackground(Bool_t bg){fCalculateBackground=bg;}\r
+ void CalculateBackground();\r
+ Double_t GetMCOpeningAngle(TParticle* daughter0, TParticle* daughter1) const;\r
+\r
+ private:\r
+ AliAnalysisTaskGammaConversion(const AliAnalysisTaskGammaConversion&); // Not implemented\r
+ AliAnalysisTaskGammaConversion& operator=(const AliAnalysisTaskGammaConversion&); // Not implemented\r
+\r
+ AliV0Reader* fV0Reader;\r
+\r
+ AliStack * fStack;\r
+\r
+ TList * fOutputContainer ; // Histogram container\r
+\r
+ AliGammaConversionHistograms *fHistograms;\r
+\r
+ Bool_t fDoMCTruth;\r
+ \r
+ vector<TParticle*> fMCAllGammas;\r
+ vector<TParticle*> fMCPi0s;\r
+ vector<TParticle*> fMCEtas;\r
+ vector<TParticle*> fMCGammaChic;\r
+\r
+ vector<AliKFParticle> fKFReconstructedGammas;\r
+\r
+ //mass defines\r
+ Double_t fElectronMass;\r
+ Double_t fGammaMass;\r
+ Double_t fPi0Mass;\r
+ Double_t fEtaMass;\r
+\r
+ // width defines\r
+ Double_t fGammaWidth;\r
+ Double_t fPi0Width;\r
+ Double_t fEtaWidth;\r
+ Bool_t fCalculateBackground;\r
+\r
+\r
+ ClassDef(AliAnalysisTaskGammaConversion, 0); // Analysis task for gamma conversions\r
+};\r
+ \r
+#endif //ALIANALYSISTASKGAMMA_H\r
--- /dev/null
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *\r
+ * Version 1.0 *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+\r
+////////////////////////////////////////////////\r
+//--------------------------------------------- \r
+// Class used to do analysis on conversion pairs\r
+//---------------------------------------------\r
+////////////////////////////////////////////////\r
+\r
+#include "AliGammaConversionHistograms.h"\r
+#include "TMath.h"\r
+#include "TObjString.h"\r
+#include "TMap.h"\r
+#include "TList.h"\r
+#include "TH1F.h"\r
+#include "TH2F.h"\r
+\r
+\r
+using namespace std;\r
+\r
+ClassImp(AliGammaConversionHistograms)\r
+\r
+\r
+AliGammaConversionHistograms::AliGammaConversionHistograms() :\r
+ fHistogramMap(new TMap()),\r
+ fNPhiIndex(0),\r
+ fNRIndex(0),\r
+ fMinRadius(0.),\r
+ fMaxRadius(0.),\r
+ fDeltaR(0.),\r
+ fMinPhi(0.),\r
+ fMaxPhi(0.),\r
+ fDeltaPhi(0.),\r
+ fMappingContainer(NULL)\r
+{\r
+ // see header file for documenation\r
+}\r
+\r
+\r
+AliGammaConversionHistograms::AliGammaConversionHistograms(const AliGammaConversionHistograms & original) :\r
+ fHistogramMap(original.fHistogramMap),\r
+ fNPhiIndex(original.fNPhiIndex),\r
+ fNRIndex(original.fNRIndex),\r
+ fMinRadius(original.fMinRadius),\r
+ fMaxRadius(original.fMaxRadius),\r
+ fDeltaR(original.fDeltaR),\r
+ fMinPhi(original.fMinPhi),\r
+ fMaxPhi(original.fMaxPhi),\r
+ fDeltaPhi(original.fDeltaPhi),\r
+ fMappingContainer(original.fMappingContainer)\r
+{ \r
+ //see header file for documentation\r
+}\r
+\r
+\r
+AliGammaConversionHistograms & AliGammaConversionHistograms::operator = (const AliGammaConversionHistograms & /*original*/)\r
+{\r
+ // assignment operator\r
+ return *this;\r
+}\r
+\r
+\r
+AliGammaConversionHistograms::~AliGammaConversionHistograms() {\r
+ //destructor\r
+\r
+\r
+}\r
+\r
+void AliGammaConversionHistograms::AddHistogram(TString histogramName, TString histogramTitle, Int_t nXBins, Double_t firstX,Double_t lastX,TString xAxisTitle, TString yAxisTitle){\r
+ // see header file for documentation\r
+ TH1F *tmp = new TH1F(histogramName, histogramTitle,nXBins,firstX,lastX);\r
+ tmp->GetXaxis()->SetTitle(xAxisTitle);\r
+ tmp->GetYaxis()->SetTitle(yAxisTitle);\r
+ TObjString* tobjstring = new TObjString(histogramName.Data());\r
+ fHistogramMap->Add((TObject*)tobjstring,(TObject*)tmp);\r
+}\r
+\r
+void AliGammaConversionHistograms::AddHistogram(TString histogramName, TString histogramTitle, Int_t nXBins, Double_t firstX, Double_t lastX, Int_t nYBins, Double_t firstY, Double_t lastY, TString xAxisTitle, TString yAxisTitle){\r
+ // see header file for documentation\r
+ TH2F *tmp = new TH2F(histogramName, histogramTitle,nXBins,firstX,lastX,nYBins,firstY,lastY);\r
+ tmp->GetXaxis()->SetTitle(xAxisTitle);\r
+ tmp->GetYaxis()->SetTitle(yAxisTitle);\r
+ TObjString *tobjstring = new TObjString(histogramName.Data());\r
+ fHistogramMap->Add((TObject*)tobjstring,(TObject*)tmp);\r
+}\r
+\r
+void AliGammaConversionHistograms::FillHistogram(TString histogramName, Double_t xValue) const{\r
+ //see header file for documentation\r
+ TH1 *tmp = (TH1*)fHistogramMap->GetValue(histogramName.Data());\r
+ if(tmp){\r
+ tmp->Fill(xValue);\r
+ }\r
+ else{\r
+ cout<<"Histogram does not exist: "<<histogramName.Data()<<endl;\r
+ }\r
+}\r
+\r
+void AliGammaConversionHistograms::FillHistogram(TString histogramName, Double_t xValue, Double_t yValue) const{\r
+ //see header file for documentation\r
+ TH1 *tmp = (TH1*)fHistogramMap->GetValue(histogramName.Data());\r
+ if(tmp){\r
+ tmp->Fill(xValue, yValue);\r
+ }\r
+ else{\r
+ cout<<"Histogram does not exist: "<<histogramName.Data()<<endl;\r
+ }\r
+}\r
+\r
+void AliGammaConversionHistograms::GetOutputContainer(TList *fOutputContainer){\r
+ //checking if the container is alrerady created\r
+\r
+ if(fOutputContainer == NULL){\r
+ //print warning\r
+ return;\r
+ }\r
+\r
+ if(fHistogramMap != NULL){\r
+ TIter iter(fHistogramMap);\r
+ TObjString *histogramName;\r
+ while ((histogramName = (TObjString*) iter.Next())) {\r
+ TString histogramString = histogramName->GetString();\r
+ if(histogramString.Contains("Mapping")){// means it should be put in the mapping folder\r
+ if(fMappingContainer == NULL){\r
+ fMappingContainer = new TList();\r
+ fMappingContainer->SetName("Mapping histograms");\r
+ }\r
+ if(fMappingContainer != NULL){\r
+ fMappingContainer->Add((TH1*)fHistogramMap->GetValue(histogramString.Data()));\r
+ }\r
+ }\r
+ fOutputContainer->Add((TH1*)fHistogramMap->GetValue(histogramString.Data()));\r
+ histogramName = NULL;\r
+ } // end while\r
+ // fOutputContainer->Add(fMappingContainer);\r
+ }\r
+}\r
+\r
+Int_t AliGammaConversionHistograms::GetRBin(Double_t radius) const{\r
+ // see header file for documentation\r
+ Int_t iResult=0;\r
+ if(fDeltaR>0){\r
+ iResult = (Int_t)((radius - fMinRadius)/fDeltaR);\r
+ }\r
+ return iResult;\r
+}\r
+\r
+Int_t AliGammaConversionHistograms::GetPhiBin(Double_t phi) const{\r
+ // see header file for documentation\r
+ Int_t iResult=0;\r
+ if(fDeltaPhi>0){\r
+ if(phi>TMath::Pi()){\r
+ phi-=2*TMath::Pi();\r
+ }\r
+ iResult = (Int_t)((phi - fMinPhi)/fDeltaPhi);\r
+ }\r
+ return iResult;\r
+}\r
+\r
+\r
+\r
+void AliGammaConversionHistograms::InitializeMappingValues(Int_t nPhiIndex, Int_t nRIndex, Int_t nBinsR, Double_t minRadius, Double_t maxRadius,Int_t nBinsPhi, Double_t minPhi, Double_t maxPhi){\r
+ // Initializing the valuse for the mapping\r
+\r
+ fNPhiIndex = nPhiIndex;\r
+ fNRIndex = nRIndex;\r
+ fMinRadius = minRadius;\r
+ fMaxRadius = maxRadius;\r
+ if(nBinsR>0 && nRIndex!=0){\r
+ fDeltaR = (fMaxRadius - fMinRadius)/nRIndex;\r
+ }\r
+ fMinPhi = minPhi;\r
+ fMaxPhi = maxPhi;\r
+ if(nBinsPhi>0 && nPhiIndex!=0){\r
+ fDeltaPhi = (fMaxPhi-fMinPhi)/nPhiIndex;\r
+ }\r
+}\r
+\r
+\r
+//mapping\r
+void AliGammaConversionHistograms::AddMappingHistograms(Int_t nPhiIndex, Int_t nRIndex,Int_t nXBins, Double_t firstX, Double_t lastX, Int_t nYBins, Double_t firstY, Double_t lastY, TString xAxisTitle, TString yAxisTitle){\r
+ // see header file for documentation\r
+ \r
+ for(Int_t phi =0; phi<=fNPhiIndex;phi++){\r
+\r
+ for(Int_t r =0; r<fNRIndex;r++){\r
+\r
+ // setting axis to "" changes below\r
+ xAxisTitle="";\r
+ yAxisTitle="";\r
+ //Creating the axis titles\r
+ if(xAxisTitle.Length() == 0){\r
+ xAxisTitle.Form("Phi %02d",phi);\r
+ }\r
+ \r
+ if(yAxisTitle.Length() == 0){\r
+ yAxisTitle.Form("R %02d",phi);\r
+ }\r
+\r
+ //MC\r
+ TString nameMC="";\r
+ nameMC.Form("MC_EP_Mapping-Phi%02d-R%02d",phi,r);\r
+ TString titleMC="";\r
+ titleMC.Form("Electron-Positron MC Mapping-Phi%02d-R%02d",phi,r);\r
+\r
+ AddHistogram(nameMC, titleMC, nXBins, firstX, lastX, nYBins, firstY, lastY, xAxisTitle, yAxisTitle);\r
+\r
+ //ESD\r
+ TString nameESD="";\r
+ nameESD.Form("ESD_EP_Mapping-Phi%02d-R%02d",phi,r);\r
+ TString titleESD="";\r
+ titleESD.Form("Electron-Positron ESD Mapping-Phi%02d-R%02d",phi,r);\r
+\r
+ AddHistogram(nameESD, titleESD, nXBins, firstX, lastX, nYBins, firstY, lastY, xAxisTitle, yAxisTitle);\r
+ }\r
+ }\r
+\r
+\r
+ for(Int_t phi =0; phi<=nPhiIndex;phi++){ \r
+\r
+ // setting axis to "" changes below\r
+ xAxisTitle="";\r
+ yAxisTitle="";\r
+ //Creating the axis titles\r
+ if(xAxisTitle.Length() == 0){\r
+ xAxisTitle.Form("Phi %02d",phi);\r
+ }\r
+ if(yAxisTitle.Length() == 0){\r
+ yAxisTitle = "Counts";\r
+ }\r
+ \r
+ //MC\r
+ TString nameMC="";\r
+ nameMC.Form("MC_EP_Mapping-Phi%02d",phi);\r
+ TString titleMC="";\r
+ titleMC.Form("Electron-Positron MC Mapping-Phi%02d",phi);\r
+ \r
+ AddHistogram(nameMC, titleMC, nXBins, firstX, lastX, xAxisTitle, yAxisTitle);\r
+\r
+ //MC\r
+ TString nameESD="";\r
+ nameESD.Form("ESD_EP_Mapping-Phi%02d",phi);\r
+ TString titleESD="";\r
+ titleESD.Form("Electron-Positron ESD Mapping-Phi%02d",phi);\r
+ \r
+ AddHistogram(nameESD, titleESD, nXBins, firstX, lastX, xAxisTitle, yAxisTitle);\r
+ }\r
+\r
+\r
+ for(Int_t r =0; r<=nRIndex;r++){\r
+\r
+ // setting axis to "" changes below\r
+ xAxisTitle="";\r
+ yAxisTitle="";\r
+ //Creating the axis titles\r
+ if(xAxisTitle.Length() == 0){\r
+ xAxisTitle.Form("R %02d",r);\r
+ }\r
+ if(yAxisTitle.Length() == 0){\r
+ yAxisTitle = "Counts";\r
+ }\r
+ \r
+ //MC\r
+ TString nameMC="";\r
+ nameMC.Form("MC_EP_Mapping-R%02d",r);\r
+ TString titleMC="";\r
+ titleMC.Form("Electron-Positron MC Mapping-R%02d",r);\r
+ \r
+ AddHistogram(nameMC, titleMC, nXBins, firstX, lastX, xAxisTitle, yAxisTitle);\r
+\r
+ //ESD\r
+ TString nameESD="";\r
+ nameESD.Form("ESD_EP_Mapping-R%02d",r);\r
+ TString titleESD="";\r
+ titleESD.Form("Electron-Positron ESD Mapping-R%02d",r);\r
+ \r
+ AddHistogram(nameESD, titleESD, nXBins, firstX, lastX, xAxisTitle, yAxisTitle);\r
+\r
+ //Mapping Phi in R\r
+ TString nameMCPhiInR="";\r
+ nameMCPhiInR.Form("MC_EP_Mapping_Phi_R-%02d",r);\r
+ TString titleMCPhiInR="";\r
+ titleMCPhiInR.Form("MC Mapping of Phi in R%02d",r);\r
+ AddHistogram(nameMCPhiInR, titleMCPhiInR, nXBins, firstX, lastX, xAxisTitle, yAxisTitle);\r
+\r
+ //Mapping Phi in R\r
+ TString nameESDPhiInR="";\r
+ nameESDPhiInR.Form("ESD_EP_Mapping_Phi_R-%02d",r);\r
+ TString titleESDPhiInR="";\r
+ titleESDPhiInR.Form("ESD Mapping of Phi in R%02d",r);\r
+ AddHistogram(nameESDPhiInR, titleESDPhiInR, nXBins, firstX, lastX, xAxisTitle, yAxisTitle); \r
+ }\r
+}\r
--- /dev/null
+#ifndef ALIGAMMACONVERSIONHISTOGRAMS_H\r
+#define ALIGAMMACONVERSIONHISTOGRAMS_H\r
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice */\r
+\r
+////////////////////////////////////////////////\r
+//--------------------------------------------- \r
+// Class used to do analysis on conversion pairs\r
+//---------------------------------------------\r
+////////////////////////////////////////////////\r
+\r
+#include "TString.h"\r
+#include "Riostream.h"\r
+#include <vector>\r
+\r
+class TMap;\r
+class TList;\r
+class TH1F;\r
+class TH2F;\r
+\r
+class AliGammaConversionHistograms{\r
+\r
+ public: \r
+ \r
+ AliGammaConversionHistograms(); //constructor\r
+ AliGammaConversionHistograms(const AliGammaConversionHistograms & original); //copy constructor\r
+ AliGammaConversionHistograms & operator = (const AliGammaConversionHistograms & original); //assignment operator\r
+ virtual ~AliGammaConversionHistograms(); //virtual destructor\r
+ \r
+\r
+ // TList * GetOutputContainer();\r
+ void GetOutputContainer(TList *fOutputContainer);\r
+ \r
+ Int_t GetRBin(Double_t radius) const;\r
+ Int_t GetPhiBin(Double_t phi) const;\r
+\r
+ void InitializeMappingValues(Int_t nPhiHistograms, Int_t nRHistograms, Int_t nBinsR, Double_t minRadius, Double_t maxRadius,Int_t nBinsPhi, Double_t minPhi, Double_t maxPhi);\r
+\r
+ void AddMappingHistograms(Int_t nPhiHistograms, Int_t nRHistograms,Int_t nXBins, Double_t firstX, Double_t lastX, Int_t nYBins, Double_t firstY, Double_t lastY, TString xAxisTitle="", TString yAxisTitle="");\r
+\r
+ /*\r
+ * Adds a TH1F histogram to the histogram map and create a key for it \r
+ */\r
+ void AddHistogram(TString histogramName, TString histogramTitle, Int_t nXBins, Double_t firstX,Double_t lastX,TString xAxisTitle="", TString yAxisTitle="");\r
+\r
+ /*\r
+ * Adds a TH2F histogram to the histogram map and create a key for it \r
+ */\r
+ void AddHistogram(TString histogramName, TString histogramTitle, Int_t nXBins, Double_t firstX, Double_t lastX, Int_t nYBins, Double_t firstY, Double_t lastY, TString xAxisTitle="", TString yAxisTitle="");\r
+\r
+ /*\r
+ * Fills a TH1F histogram with the given name with the given value \r
+ */\r
+ void FillHistogram(TString histogramName, Double_t xValue) const;\r
+\r
+ /*\r
+ * Fills a TH2F histogram with the given name with the given value \r
+ */\r
+ void FillHistogram(TString histogramName, Double_t xValue, Double_t yValue) const;\r
+\r
+ private:\r
+ TMap* fHistogramMap;\r
+\r
+ Int_t fNPhiIndex;\r
+ Int_t fNRIndex;\r
+ Double_t fMinRadius;\r
+ Double_t fMaxRadius;\r
+ Double_t fDeltaR;\r
+ Double_t fMinPhi;\r
+ Double_t fMaxPhi;\r
+ Double_t fDeltaPhi;\r
+\r
+ TList * fMappingContainer;\r
+\r
+ \r
+ ClassDef(AliGammaConversionHistograms,1)\r
+} ;\r
+\r
+\r
+#endif\r
+\r
+\r
+\r
--- /dev/null
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *\r
+ * Version 1.0 *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+\r
+////////////////////////////////////////////////\r
+//--------------------------------------------- \r
+// Class used to do analysis on conversion pairs\r
+//---------------------------------------------\r
+////////////////////////////////////////////////\r
+\r
+// --- ROOT system ---\r
+#include <TMath.h>\r
+\r
+//---- ANALYSIS system ----\r
+#include "AliV0Reader.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliESDInputHandler.h"\r
+#include "AliMCEvent.h"\r
+#include "AliKFVertex.h"\r
+\r
+#include "AliStack.h"\r
+#include "AliMCEventHandler.h"\r
+\r
+\r
+class iostream;\r
+class AliESDv0;\r
+class TFormula;\r
+\r
+using namespace std;\r
+\r
+ClassImp(AliV0Reader)\r
+\r
+\r
+\r
+ AliV0Reader::AliV0Reader() :\r
+ TObject(),\r
+ fMCStack(NULL),\r
+ fMCTruth(NULL),\r
+ fChain(NULL),\r
+ fESDHandler(NULL),\r
+ fESDEvent(NULL),\r
+ fHistograms(NULL),\r
+ fCurrentV0IndexNumber(0),\r
+ fCurrentV0(NULL),\r
+ fCurrentNegativeKFParticle(NULL),\r
+ fCurrentPositiveKFParticle(NULL),\r
+ fCurrentMotherKFCandidate(NULL),\r
+ fCurrentNegativeESDTrack(NULL),\r
+ fCurrentPositiveESDTrack(NULL),\r
+ fNegativeTrackLorentzVector(NULL),\r
+ fPositiveTrackLorentzVector(NULL),\r
+ fMotherCandidateLorentzVector(NULL),\r
+ fCurrentXValue(0),\r
+ fCurrentYValue(0),\r
+ fCurrentZValue(0),\r
+ fPositiveTrackPID(0),\r
+ fNegativeTrackPID(0),\r
+ fNegativeMCParticle(NULL),\r
+ fPositiveMCParticle(NULL),\r
+ fMotherMCParticle(NULL),\r
+ fMotherCandidateKFMass(0),\r
+ fMotherCandidateKFWidth(0),\r
+ fUseKFParticle(kTRUE),\r
+ fUseESDTrack(kFALSE),\r
+ fDoMC(kFALSE),\r
+ fMaxR(10000),// 100 meter(outside of ALICE)\r
+ fEtaCut(0.),\r
+ fPtCut(0.),\r
+ fChi2CutConversion(0.),\r
+ fChi2CutMeson(0.),\r
+ fPIDProbabilityCutNegativeParticle(0),\r
+ fPIDProbabilityCutPositiveParticle(0),\r
+ fXVertexCut(0.),\r
+ fYVertexCut(0.),\r
+ fZVertexCut(0.),\r
+ fNSigmaMass(0.),\r
+ fUseImprovedVertex(kFALSE),\r
+ fCurrentEventGoodV0s(),\r
+ fPreviousEventGoodV0s()\r
+{\r
+\r
+}\r
+\r
+\r
+AliV0Reader::AliV0Reader(const AliV0Reader & original) :\r
+ TObject(original),\r
+ fMCStack(original.fMCStack),\r
+ fMCTruth(original.fMCTruth),\r
+ fChain(original.fChain),\r
+ fESDHandler(original.fESDHandler),\r
+ fESDEvent(original.fESDEvent),\r
+ fHistograms(original.fHistograms),\r
+ fCurrentV0IndexNumber(original.fCurrentV0IndexNumber),\r
+ fCurrentV0(original.fCurrentV0),\r
+ fCurrentNegativeKFParticle(original.fCurrentNegativeKFParticle),\r
+ fCurrentPositiveKFParticle(original.fCurrentPositiveKFParticle),\r
+ fCurrentMotherKFCandidate(original.fCurrentMotherKFCandidate),\r
+ fCurrentNegativeESDTrack(original.fCurrentNegativeESDTrack),\r
+ fCurrentPositiveESDTrack(original.fCurrentPositiveESDTrack),\r
+ fNegativeTrackLorentzVector(original.fNegativeTrackLorentzVector),\r
+ fPositiveTrackLorentzVector(original.fPositiveTrackLorentzVector),\r
+ fMotherCandidateLorentzVector(original.fMotherCandidateLorentzVector),\r
+ fCurrentXValue(original.fCurrentXValue),\r
+ fCurrentYValue(original.fCurrentYValue),\r
+ fCurrentZValue(original.fCurrentZValue),\r
+ fPositiveTrackPID(original.fPositiveTrackPID),\r
+ fNegativeTrackPID(original.fNegativeTrackPID),\r
+ fNegativeMCParticle(original.fNegativeMCParticle),\r
+ fPositiveMCParticle(original.fPositiveMCParticle),\r
+ fMotherMCParticle(original.fMotherMCParticle),\r
+ fMotherCandidateKFMass(original.fMotherCandidateKFMass),\r
+ fMotherCandidateKFWidth(original.fMotherCandidateKFWidth),\r
+ fUseKFParticle(kTRUE),\r
+ fUseESDTrack(kFALSE),\r
+ fDoMC(kFALSE),\r
+ fMaxR(original.fMaxR),\r
+ fEtaCut(original.fEtaCut),\r
+ fPtCut(original.fPtCut),\r
+ fChi2CutConversion(original.fChi2CutConversion),\r
+ fChi2CutMeson(original.fChi2CutMeson),\r
+ fPIDProbabilityCutNegativeParticle(original.fPIDProbabilityCutNegativeParticle),\r
+ fPIDProbabilityCutPositiveParticle(original.fPIDProbabilityCutPositiveParticle),\r
+ fXVertexCut(original.fXVertexCut),\r
+ fYVertexCut(original.fYVertexCut),\r
+ fZVertexCut(original.fZVertexCut),\r
+ fNSigmaMass(original.fNSigmaMass),\r
+ fUseImprovedVertex(original.fUseImprovedVertex),\r
+ fCurrentEventGoodV0s(original.fCurrentEventGoodV0s),\r
+ fPreviousEventGoodV0s(original.fPreviousEventGoodV0s)\r
+{\r
+\r
+}\r
+\r
+\r
+AliV0Reader & AliV0Reader::operator = (const AliV0Reader & /*source*/)\r
+{\r
+ // assignment operator\r
+ return *this;\r
+}\r
+\r
+void AliV0Reader::Initialize(){\r
+ //see header file for documentation\r
+\r
+ // Get the input handler from the manager\r
+ fESDHandler = (AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
+ if(fESDHandler == NULL){\r
+ //print warning here\r
+ }\r
+ \r
+ // Get pointer to esd event from input handler\r
+ fESDEvent = fESDHandler->GetEvent();\r
+ if(fESDEvent == NULL){\r
+ //print warning here\r
+ }\r
+\r
+ //Get pointer to MCTruth\r
+ fMCTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());\r
+ if(fMCTruth == NULL){\r
+ //print warning here\r
+ }\r
+\r
+ //Get pointer to the mc stack\r
+ fMCStack = fMCTruth->MCEvent()->Stack();\r
+ if(fMCStack == NULL){\r
+ //print warning here\r
+ }\r
+\r
+ AliKFParticle::SetField(fESDEvent->GetMagneticField());\r
+\r
+}\r
+\r
+AliESDv0* AliV0Reader::GetV0(Int_t index){\r
+ //see header file for documentation\r
+\r
+ fCurrentV0 = fESDEvent->GetV0(index);\r
+ UpdateV0Information();\r
+ return fCurrentV0;\r
+}\r
+\r
+\r
+Bool_t AliV0Reader::NextV0(){\r
+ //see header file for documentation\r
+\r
+ Bool_t iResult=kFALSE;\r
+ while(fCurrentV0IndexNumber<fESDEvent->GetNumberOfV0s()){\r
+ fCurrentV0 = fESDEvent->GetV0(fCurrentV0IndexNumber);\r
+ \r
+ //checks if on the fly mode is set\r
+ if ( !fCurrentV0->GetOnFlyStatus() ){\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut1",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+\r
+ if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) {//checks if we have a vertex\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut2",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+\r
+ if(CheckPIDProbability(fPIDProbabilityCutNegativeParticle,fPIDProbabilityCutPositiveParticle)==kFALSE){\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut3",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+\r
+\r
+ fCurrentV0->GetXYZ(fCurrentXValue,fCurrentYValue,fCurrentZValue);\r
+ \r
+ if(GetXYRadius()>fMaxR){ // cuts on distance from collision point\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut4",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+\r
+ UpdateV0Information();\r
+ \r
+ if(fUseKFParticle){\r
+ if(fCurrentMotherKFCandidate->GetNDF()<=0){\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut5",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+ Double_t chi2V0 = fCurrentMotherKFCandidate->GetChi2()/fCurrentMotherKFCandidate->GetNDF();\r
+ if(chi2V0 > fChi2CutConversion || chi2V0 <=0){\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut6",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+ \r
+ if(TMath::Abs(fMotherCandidateLorentzVector->Eta())> fEtaCut){\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut7",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+ \r
+ if(fMotherCandidateLorentzVector->Pt()<fPtCut){\r
+ fCurrentV0IndexNumber++;\r
+ fHistograms->FillHistogram("V0MassDebugCut8",GetMotherCandidateMass());\r
+ continue;\r
+ }\r
+ \r
+ }\r
+ else if(fUseESDTrack){\r
+ //TODO\r
+ }\r
+\r
+ iResult=kTRUE;//means we have a v0 who survived all the cuts applied\r
+\r
+ fCurrentV0IndexNumber++;\r
+ \r
+ break;\r
+ }\r
+ return iResult; \r
+}\r
+\r
+void AliV0Reader::UpdateV0Information(){\r
+ //see header file for documentation\r
+ \r
+ if(fCurrentNegativeKFParticle != NULL){\r
+ delete fCurrentNegativeKFParticle;\r
+ }\r
+ fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fNegativeTrackPID);\r
+ \r
+ if(fCurrentPositiveKFParticle != NULL){\r
+ delete fCurrentPositiveKFParticle;\r
+ }\r
+ fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fPositiveTrackPID);\r
+ \r
+ if(fCurrentMotherKFCandidate != NULL){\r
+ delete fCurrentMotherKFCandidate;\r
+ }\r
+ fCurrentMotherKFCandidate = new AliKFParticle(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle);\r
+\r
+ fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
+\r
+ fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
+\r
+ if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){\r
+ fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass);\r
+ }\r
+\r
+ if(fUseImprovedVertex == kTRUE){\r
+ AliKFVertex primaryVertexImproved(*GetPrimaryVertex());\r
+ primaryVertexImproved+=*fCurrentMotherKFCandidate;\r
+ fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved);\r
+ }\r
+\r
+ fCurrentMotherKFCandidate->GetMass(fMotherCandidateKFMass,fMotherCandidateKFWidth);\r
+\r
+\r
+ if(fNegativeTrackLorentzVector != NULL){\r
+ delete fNegativeTrackLorentzVector;\r
+ }\r
+ if(fUseKFParticle){\r
+ fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeKFParticle->Px(),fCurrentNegativeKFParticle->Py(),fCurrentNegativeKFParticle->Pz());\r
+ }\r
+ else if(fUseESDTrack){\r
+ fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeESDTrack->Px(),fCurrentNegativeESDTrack->Py(),fCurrentNegativeESDTrack->Pz());\r
+ }\r
+\r
+ if(fPositiveTrackLorentzVector != NULL){\r
+ delete fPositiveTrackLorentzVector;\r
+ }\r
+ if(fUseKFParticle){\r
+ fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveKFParticle->Px(),fCurrentPositiveKFParticle->Py(),fCurrentPositiveKFParticle->Pz());\r
+ }\r
+ else if(fUseESDTrack){\r
+ fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveESDTrack->Px(),fCurrentPositiveESDTrack->Py(),fCurrentPositiveESDTrack->Pz());\r
+ }\r
+\r
+ if(fMotherCandidateLorentzVector != NULL){\r
+ delete fMotherCandidateLorentzVector;\r
+ }\r
+ if(fUseKFParticle){\r
+ fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);\r
+ }\r
+ else if(fUseESDTrack){\r
+ fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);\r
+ }\r
+\r
+ if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){\r
+ fMotherCandidateLorentzVector->SetXYZM(fMotherCandidateLorentzVector->Px() ,fMotherCandidateLorentzVector->Py(),fMotherCandidateLorentzVector->Pz(),0.); \r
+ }\r
+ \r
+ if(fDoMC == kTRUE){\r
+ fNegativeMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetNindex())->GetLabel()));\r
+ fPositiveMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetPindex())->GetLabel()));\r
+ }\r
+ fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate);\r
+}\r
+\r
+Bool_t AliV0Reader::HasSameMCMother(){\r
+ //see header file for documentation\r
+\r
+ Bool_t iResult = kFALSE;\r
+ if(fDoMC == kTRUE){\r
+ if(fNegativeMCParticle != NULL && fPositiveMCParticle != NULL){\r
+ if(fNegativeMCParticle->GetMother(0) == fPositiveMCParticle->GetMother(0))\r
+ fMotherMCParticle = fMCStack->Particle(fPositiveMCParticle->GetMother(0));\r
+ iResult = kTRUE;\r
+ }\r
+ }\r
+ return iResult;\r
+}\r
+\r
+Bool_t AliV0Reader::CheckPIDProbability(Double_t negProbCut, Double_t posProbCut){\r
+ //see header file for documentation\r
+\r
+ Bool_t iResult=kFALSE;\r
+\r
+ Double_t *posProbArray = new Double_t[10];\r
+ Double_t *negProbArray = new Double_t[10];\r
+ AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
+ AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
+ \r
+ negTrack->GetTPCpid(negProbArray);\r
+ posTrack->GetTPCpid(posProbArray);\r
+\r
+ if(negProbArray!=NULL && posProbArray!=NULL){\r
+ if(negProbArray[GetSpeciesIndex(-1)]>=negProbCut && posProbArray[GetSpeciesIndex(1)]>=posProbCut){\r
+ iResult=kTRUE;\r
+ }\r
+ }\r
+ delete [] posProbArray;\r
+ delete [] negProbArray;\r
+ return iResult;\r
+}\r
+\r
+void AliV0Reader::UpdateEventByEventData(){\r
+ //see header file for documentation\r
+\r
+ if(fCurrentEventGoodV0s.size() >0 ){\r
+ fPreviousEventGoodV0s.clear();\r
+ fPreviousEventGoodV0s = fCurrentEventGoodV0s;\r
+ }\r
+ fCurrentEventGoodV0s.clear();\r
+ \r
+ fCurrentV0IndexNumber=0;\r
+}\r
+\r
+Double_t AliV0Reader::GetNegativeTrackPhi() const{\r
+ //see header file for documentation\r
+\r
+ Double_t offset=0;\r
+ if(fNegativeTrackLorentzVector->Phi()> TMath::Pi()){\r
+ offset = -2*TMath::Pi();\r
+ }\r
+ return fNegativeTrackLorentzVector->Phi()+offset;\r
+}\r
+\r
+Double_t AliV0Reader::GetPositiveTrackPhi() const{\r
+ //see header file for documentation\r
+\r
+ Double_t offset=0;\r
+ if(fPositiveTrackLorentzVector->Phi()> TMath::Pi()){\r
+ offset = -2*TMath::Pi();\r
+ }\r
+ return fPositiveTrackLorentzVector->Phi()+offset;\r
+}\r
+\r
+Double_t AliV0Reader::GetMotherCandidatePhi() const{\r
+ //see header file for documentation\r
+\r
+ Double_t offset=0;\r
+ if(fMotherCandidateLorentzVector->Phi()> TMath::Pi()){\r
+ offset = -2*TMath::Pi();\r
+ }\r
+ return fMotherCandidateLorentzVector->Phi()+offset;\r
+}\r
+\r
+Int_t AliV0Reader::GetSpeciesIndex(Int_t chargeOfTrack){\r
+ //see header file for documentation\r
+\r
+ Int_t iResult = 10; // Unknown particle\r
+\r
+ if(chargeOfTrack==-1){ //negative track\r
+ switch(abs(fNegativeTrackPID)){\r
+ case 11: //electron\r
+ iResult = 0;\r
+ break;\r
+ case 13: //muon\r
+ iResult = 1;\r
+ break;\r
+ case 211: //pion\r
+ iResult = 2;\r
+ break;\r
+ case 321: //kaon\r
+ iResult = 3;\r
+ break;\r
+ case 2212: //proton\r
+ iResult = 4;\r
+ break;\r
+ case 22: //photon\r
+ iResult = 5;\r
+ break;\r
+ case 111: //pi0\r
+ iResult = 6;\r
+ break;\r
+ case 2112: //neutron\r
+ iResult = 7;\r
+ break;\r
+ case 311: //K0\r
+ iResult = 8;\r
+ break;\r
+ \r
+ //Put in here for kSPECIES::kEleCon ????\r
+ }\r
+ }\r
+ else if(chargeOfTrack==1){ //positive track\r
+ switch(abs(fPositiveTrackPID)){\r
+ case 11: //electron\r
+ iResult = 0;\r
+ break;\r
+ case 13: //muon\r
+ iResult = 1;\r
+ break;\r
+ case 211: //pion\r
+ iResult = 2;\r
+ break;\r
+ case 321: //kaon\r
+ iResult = 3;\r
+ break;\r
+ case 2212: //proton\r
+ iResult = 4;\r
+ break;\r
+ case 22: //photon\r
+ iResult = 5;\r
+ break;\r
+ case 111: //pi0\r
+ iResult = 6;\r
+ break;\r
+ case 2112: //neutron\r
+ iResult = 7;\r
+ break;\r
+ case 311: //K0\r
+ iResult = 8;\r
+ break;\r
+\r
+ //Put in here for kSPECIES::kEleCon ????\r
+ }\r
+ }\r
+ else{\r
+ //Wrong parameter.. Print warning\r
+ }\r
+ return iResult;\r
+}\r
--- /dev/null
+#ifndef ALIV0READER_H\r
+#define ALIV0READER_H\r
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice */\r
+\r
+////////////////////////////////////////////////\r
+//--------------------------------------------- \r
+// Class used to do analysis on conversion pairs\r
+//---------------------------------------------\r
+////////////////////////////////////////////////\r
+\r
+// --- ROOT system ---\r
+#include "TObject.h" \r
+#include "AliESDv0.h"\r
+#include "AliESDEvent.h"\r
+#include "AliKFParticle.h"\r
+#include "TParticle.h"\r
+#include "AliGammaConversionHistograms.h"\r
+#include <vector>\r
+\r
+class TClonesArray; \r
+class TFormula;\r
+class Riostream;\r
+class TChain;\r
+//--- AliRoot system ---\r
+\r
+class AliStack;\r
+class AliESDEvent; \r
+class AliMCEventHandler;\r
+class AliESDInputHandler;\r
+class AliESDVertex;\r
+class AliLog;\r
+class TChain;\r
+class TChain;\r
+\r
+\r
+\r
+class AliV0Reader : public TObject {\r
+\r
+ public: \r
+\r
+ AliV0Reader(); //constructor\r
+ AliV0Reader(const AliV0Reader & g); //copy constructor\r
+ AliV0Reader & operator = (const AliV0Reader & g); //assignment operator\r
+ virtual ~AliV0Reader() {;} //virtual destructor\r
+ /*\r
+ *Initialize the reader\r
+ */\r
+ void Initialize();\r
+\r
+ /*\r
+ *Returns the number of v0s in the event, no cuts applied.\r
+ */\r
+ Int_t GetNumberOfV0s() const{return fESDEvent->GetNumberOfV0s();}\r
+\r
+ /*\r
+ * Check if there are any more good v0s left in the v0 stack\r
+ * if so, fCurrent v0 is set to this v0 and can be retrieved\r
+ * by GetCurrentV0 function.\r
+ * returns kFALSE if there is no more good v0s in the v0 stack\r
+ */\r
+ Bool_t NextV0();\r
+ \r
+ /*\r
+ * Returns the v0 at the given index, no checks are done on the v0. \r
+ */\r
+ AliESDv0* GetV0(Int_t index);\r
+\r
+ /*\r
+ * Returns the current v0\r
+ */\r
+ AliESDv0* GetCurrentV0() const{return fCurrentV0;}\r
+\r
+ /*\r
+ * Returns the negative ESD track which belongs to fCurrentV0\r
+ */\r
+ AliESDtrack* GetNegativeESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetNindex());}\r
+\r
+ /*\r
+ * Returns the positive ESD track which belongs to fCurrentV0\r
+ */\r
+ AliESDtrack* GetPositiveESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetPindex());}\r
+\r
+ /*\r
+ * Returns the negative KF particle which belongs to fCurrentV0\r
+ */\r
+ AliKFParticle* GetNegativeKFParticle() const{return fCurrentNegativeKFParticle;}\r
+\r
+ /*\r
+ * Returns the positive KF particle which belongs to fCurrentV0\r
+ */\r
+ AliKFParticle* GetPositiveKFParticle() const{return fCurrentPositiveKFParticle;}\r
+ /*\r
+ * Returns the KFParticle object of the 2 tracks.\r
+ */\r
+ AliKFParticle* GetMotherCandidateKFCombination() const{return fCurrentMotherKFCandidate;}\r
+ /*\r
+ * Checks the probablity that the PID of the particle is what we want it to be.\r
+ */\r
+ Bool_t CheckPIDProbability(Double_t negProbCut, Double_t posProbCut);\r
+\r
+ /*\r
+ *Get the negative MC TParticle from the stack \r
+ */\r
+ TParticle * GetNegativeMCParticle() const{return fNegativeMCParticle;}\r
+\r
+ /*\r
+ *Get the positive MC TParticle from the stack \r
+ */\r
+ TParticle * GetPositiveMCParticle() const{return fPositiveMCParticle;}\r
+\r
+ /*\r
+ *Get the mother MC TParticle from the stack \r
+ */\r
+ TParticle * GetMotherMCParticle() const{return fMotherMCParticle;}\r
+\r
+ Bool_t HasSameMCMother();\r
+\r
+ /*\r
+ *Get the MC stack \r
+ */\r
+ AliStack* GetMCStack() const{return fMCStack;}\r
+\r
+ /*\r
+ *Get the magnetic field from the ESD event \r
+ */\r
+ Double_t GetMagneticField() const{return fESDEvent->GetMagneticField();}\r
+\r
+ /*\r
+ *Get the primary vertex from the esd event\r
+ */\r
+ const AliESDVertex *GetPrimaryVertex() const {return fESDEvent->GetPrimaryVertex();}\r
+\r
+ /*\r
+ * Set the PID of the negative track\r
+ */\r
+ void SetNegativeTrackPID(Int_t negTrackPID){fNegativeTrackPID=negTrackPID;}\r
+\r
+ /*\r
+ * Set the PID of the positive track\r
+ */\r
+ void SetPositiveTrackPID(Int_t posTrackPID){fPositiveTrackPID=posTrackPID;}\r
+\r
+ /*\r
+ * Set the flag to use the kfparticle class. Will also disable the use of esd tracks\r
+ */\r
+ void UseKFParticle(){fUseKFParticle = kTRUE; fUseESDTrack = kFALSE;}\r
+\r
+ /*\r
+ * Set the flag to use the esd track class. Will also disable the use of kf particles\r
+ */\r
+ void UseESDTrack(){fUseESDTrack = kTRUE; fUseKFParticle = kFALSE;}\r
+\r
+ /*\r
+ * Set the flag to use improved vertex or not\r
+ */\r
+ void SetUseImprovedVertex(Bool_t useImprovedVertex){fUseImprovedVertex=useImprovedVertex;}\r
+\r
+ /*\r
+ * Return the number in the species array belonging to the negative or positive track pid.\r
+ */\r
+ Int_t GetSpeciesIndex(Int_t chargeOfTrack);\r
+\r
+ /*\r
+ * Return the x coordinate of the v0\r
+ */\r
+ Double_t GetX() const{return fCurrentXValue;}\r
+\r
+ /*\r
+ * Return the y coordinate of the v0\r
+ */\r
+ Double_t GetY() const{return fCurrentYValue;}\r
+\r
+ /*\r
+ * Return the Z coordinate of the v0\r
+ */\r
+ Double_t GetZ() const{return fCurrentZValue;}\r
+\r
+ /*\r
+ * Return the radius of the v0\r
+ */\r
+ Double_t GetXYRadius() const{return sqrt((Double_t)(fCurrentXValue*fCurrentXValue + fCurrentYValue*fCurrentYValue));}\r
+\r
+ /*\r
+ * Get the opening angle between the two tracks\r
+ */\r
+ Double_t GetOpeningAngle(){return fNegativeTrackLorentzVector->Angle(fPositiveTrackLorentzVector->Vect());}\r
+\r
+ Double_t GetNegativeTrackEnergy() const{return fCurrentNegativeKFParticle->E();}\r
+ Double_t GetPositiveTrackEnergy() const{return fCurrentPositiveKFParticle->E();}\r
+ Double_t GetMotherCandidateEnergy() const{return fCurrentMotherKFCandidate->E();}\r
+\r
+ Double_t GetNegativeTrackPt() const{return fNegativeTrackLorentzVector->Pt();}\r
+ Double_t GetPositiveTrackPt() const{return fPositiveTrackLorentzVector->Pt();}\r
+ Double_t GetMotherCandidatePt() const{return fMotherCandidateLorentzVector->Pt();}\r
+\r
+ Double_t GetNegativeTrackEta() const{return fNegativeTrackLorentzVector->Eta();}\r
+ Double_t GetPositiveTrackEta() const{return fPositiveTrackLorentzVector->Eta();}\r
+ Double_t GetMotherCandidateEta() const{return fMotherCandidateLorentzVector->Eta();}\r
+\r
+ Double_t GetMotherCandidateNDF() const{return fCurrentMotherKFCandidate->GetNDF();}\r
+ Double_t GetMotherCandidateChi2() const{return fCurrentMotherKFCandidate->GetChi2();}\r
+ Double_t GetMotherCandidateMass() const{return fMotherCandidateKFMass;}\r
+ Double_t GetMotherCandidateWidth() const{return fMotherCandidateKFWidth;}\r
+\r
+ Double_t GetNegativeTrackPhi() const;\r
+ Double_t GetPositiveTrackPhi() const;\r
+ Double_t GetMotherCandidatePhi() const;\r
+\r
+ void UpdateEventByEventData();\r
+ \r
+ Double_t GetMaxRCut() const{return fMaxR;}\r
+ Double_t GetEtaCut() const{return fEtaCut;}\r
+ Double_t GetPtCut() const{return fPtCut;}\r
+ Double_t GetChi2CutConversion() const{return fChi2CutConversion;}\r
+ Double_t GetChi2CutMeson() const{return fChi2CutMeson;}\r
+\r
+ void SetMaxRCut(Double_t maxR){fMaxR=maxR;}\r
+ void SetEtaCut(Double_t etaCut){fEtaCut=etaCut;}\r
+ void SetPtCut(Double_t ptCut){fPtCut=ptCut;}\r
+ void SetChi2CutConversion(Double_t chi2){fChi2CutConversion=chi2;}\r
+ void SetChi2CutMeson(Double_t chi2){fChi2CutMeson=chi2;}\r
+ \r
+ void SetXVertexCut(Double_t xVtx){fCurrentXValue=xVtx;}\r
+ void SetYVertexCut(Double_t yVtx){fCurrentYValue=yVtx;}\r
+ void SetZVertexCut(Double_t zVtx){fCurrentZValue=zVtx;}\r
+ void SetPIDProbability(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb; fPIDProbabilityCutNegativeParticle=pidProb;}\r
+ void SetPIDProbabilityNegativeParticle(Double_t pidProb){fPIDProbabilityCutNegativeParticle=pidProb;}\r
+ void SetPIDProbabilityPositiveParticle(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb;}\r
+ void SetSigmaMass(Double_t sigmaMass){fNSigmaMass=sigmaMass;}\r
+\r
+ void SetDoMCTruth(Bool_t doMC){fDoMC = doMC;}\r
+\r
+ void UpdateV0Information();\r
+\r
+ void SetHistograms(AliGammaConversionHistograms *histograms){fHistograms=histograms;}\r
+\r
+ vector<AliKFParticle> GetCurrentEventGoodV0s() const{return fCurrentEventGoodV0s;}\r
+ vector<AliKFParticle> GetPreviousEventGoodV0s() const{return fPreviousEventGoodV0s;}\r
+\r
+ private:\r
+ AliStack * fMCStack; // pointer to MonteCarlo particle stack \r
+ AliMCEventHandler* fMCTruth; // pointer to the MC event handler\r
+ TChain * fChain; // pointer to the TChain\r
+ \r
+ AliESDInputHandler* fESDHandler; //! pointer to esd object\r
+ AliESDEvent *fESDEvent; //! pointer to esd object\r
+\r
+ AliGammaConversionHistograms *fHistograms;\r
+ \r
+ Int_t fCurrentV0IndexNumber;\r
+ AliESDv0 * fCurrentV0; //! pointer to the current v0\r
+ AliKFParticle * fCurrentNegativeKFParticle; //! pointer to the negative KF particle\r
+ AliKFParticle * fCurrentPositiveKFParticle; //! pointer to the positive KF particle\r
+ AliKFParticle * fCurrentMotherKFCandidate; //! pointer to the positive KF particle\r
+\r
+ AliESDtrack * fCurrentNegativeESDTrack; //! pointer to the negative ESD track\r
+ AliESDtrack * fCurrentPositiveESDTrack; //! pointer to the positive ESD track\r
+ \r
+ TLorentzVector * fNegativeTrackLorentzVector; //! pointer to the negative Track Lorentz Vector\r
+ TLorentzVector * fPositiveTrackLorentzVector; //! pointer to the positive Track Lorentz Vector\r
+ TLorentzVector * fMotherCandidateLorentzVector; //! pointer to the mother candidate Track Lorentz Vector\r
+\r
+ Double_t fCurrentXValue;\r
+ Double_t fCurrentYValue;\r
+ Double_t fCurrentZValue;\r
+\r
+ Int_t fPositiveTrackPID;\r
+ Int_t fNegativeTrackPID;\r
+\r
+ TParticle *fNegativeMCParticle; //!\r
+ TParticle *fPositiveMCParticle; //!\r
+ TParticle *fMotherMCParticle; //!\r
+\r
+ Double_t fMotherCandidateKFMass;\r
+ Double_t fMotherCandidateKFWidth;\r
+\r
+ Bool_t fUseKFParticle;\r
+ Bool_t fUseESDTrack;\r
+ Bool_t fDoMC;\r
+\r
+ //cuts\r
+ Double_t fMaxR;\r
+ Double_t fEtaCut;\r
+ Double_t fPtCut;\r
+ Double_t fChi2CutConversion;\r
+ Double_t fChi2CutMeson;\r
+ Double_t fPIDProbabilityCutNegativeParticle;\r
+ Double_t fPIDProbabilityCutPositiveParticle;\r
+ Double_t fXVertexCut;\r
+ Double_t fYVertexCut;\r
+ Double_t fZVertexCut;\r
+ \r
+ Double_t fNSigmaMass;\r
+ \r
+ Bool_t fUseImprovedVertex;\r
+ \r
+ vector<AliKFParticle> fCurrentEventGoodV0s;\r
+ vector<AliKFParticle> fPreviousEventGoodV0s;\r
+\r
+ ClassDef(AliV0Reader,0)\r
+};\r
+\r
+\r
+#endif\r
+\r
+\r
+\r
--- /dev/null
+#ifdef __CINT__
+
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+
+#pragma link C++ class AliAnalysisTaskGammaConversion+;
+#pragma link C++ class AliV0Reader+;
+#pragma link C++ class AliGammaConversionHistograms+;
+
+#endif
--- /dev/null
+#-*- Mode: Makefile -*-
+
+SRCS = GammaConv/AliV0Reader.cxx GammaConv/AliAnalysisTaskGammaConversion.cxx GammaConv/AliGammaConversionHistograms.cxx
+
+HDRS:= $(SRCS:.cxx=.h)
+
+DHDR:= PWG4GammaConvLinkDef.h
+
+EXPORT:=$(SRCS:.cxx=.h)
+
+ifeq (win32gcc,$(ALICE_TARGET))
+PACKSOFLAGS:= $(SOFLAGS) -L$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET) -lSTEERBase \
+ -lESD -lAOD -lANALYSIS -lANALYSISalice \
+ -L$(shell root-config --libdir) -lEG
+endif