#include "TNtuple.h"
//#include "AliCFManager.h" // for CF
//#include "AliCFContainer.h" // for CF
+#include "AliESDInputHandler.h"
+#include "AliAnalysisManager.h"
#include "AliGammaConversionAODObject.h"
+#include "AliGammaConversionBGHandler.h"
+#include "AliESDCaloCluster.h" // for combining PHOS and GammaConv
+#include "AliKFVertex.h"
+#include "AliGenPythiaEventHeader.h"
+#include "AliGenDPMjetEventHeader.h"
+#include "AliGenEventHeader.h"
+#include <AliMCEventHandler.h>
+#include "TRandom3.h"
+#include "AliTriggerAnalysis.h"
+#include "AliESDCentrality.h"
+class AliESDTrackCuts;
class AliCFContainer;
class AliCFManager;
class AliKFVertex;
fOutputContainer(NULL),
fCFManager(0x0), // for CF
fHistograms(NULL),
+ fTriggerCINT1B(kFALSE),
fDoMCTruth(kFALSE),
fDoNeutralMeson(kFALSE),
+ fDoOmegaMeson(kFALSE),
fDoJet(kFALSE),
fDoChic(kFALSE),
+ fRecalculateV0ForGamma(kFALSE),
fKFReconstructedGammasTClone(NULL),
+ fKFReconstructedPi0sTClone(NULL),
+ fKFRecalculatedGammasTClone(NULL),
fCurrentEventPosElectronTClone(NULL),
fCurrentEventNegElectronTClone(NULL),
fKFReconstructedGammasCutTClone(NULL),
fPreviousEventTLVNegElectronTClone(NULL),
fPreviousEventTLVPosElectronTClone(NULL),
-// fKFReconstructedGammas(),
fElectronv1(),
fElectronv2(),
-// fCurrentEventPosElectron(),
-// fCurrentEventNegElectron(),
-// fKFReconstructedGammasCut(),
-// fPreviousEventTLVNegElectron(),
-// fPreviousEventTLVPosElectron(),
+ fGammav1(),
+ fGammav2(),
+ fElectronRecalculatedv1(),
+ fElectronRecalculatedv2(),
fElectronMass(-1),
fGammaMass(-1),
fPi0Mass(-1),
fLowPtMapping(1.),
fHighPtMapping(3.),
fDoCF(kFALSE),
- fAODBranch(NULL),
- fAODBranchName("GammaConv")//,
- // fAODObjects(NULL)
+ fAODGamma(NULL),
+ fAODPi0(NULL),
+ fAODOmega(NULL),
+ fAODBranchName("GammaConv"),
+ fKFForceAOD(kFALSE),
+ fKFDeltaAODFileName(""),
+ fDoNeutralMesonV0MCCheck(kFALSE),
+ fUseTrackMultiplicityForBG(kTRUE),
+ fMoveParticleAccordingToVertex(kFALSE),
+ fApplyChi2Cut(kFALSE),
+ fNRandomEventsForBG(15),
+ fNDegreesPMBackground(15),
+ fDoRotation(kTRUE),
+ fCheckBGProbability(kTRUE),
+ fKFReconstructedGammasV0Index(),
+ fRemovePileUp(kFALSE),
+ fSelectV0AND(kFALSE),
+ fTriggerAnalysis(NULL),
+ fMultiplicity(0),
+ fUseMultiplicity(0),
+ fUseMultiplicityBin(0),
+ fUseCentrality(0),
+ fUseCentralityBin(0)
{
// Default constructor
fOutputContainer(0x0),
fCFManager(0x0), // for CF
fHistograms(NULL),
+ fTriggerCINT1B(kFALSE),
fDoMCTruth(kFALSE),
fDoNeutralMeson(kFALSE),
+ fDoOmegaMeson(kFALSE),
fDoJet(kFALSE),
fDoChic(kFALSE),
+ fRecalculateV0ForGamma(kFALSE),
fKFReconstructedGammasTClone(NULL),
+ fKFReconstructedPi0sTClone(NULL),
+ fKFRecalculatedGammasTClone(NULL),
fCurrentEventPosElectronTClone(NULL),
fCurrentEventNegElectronTClone(NULL),
fKFReconstructedGammasCutTClone(NULL),
fPreviousEventTLVNegElectronTClone(NULL),
fPreviousEventTLVPosElectronTClone(NULL),
- // fKFReconstructedGammas(),
fElectronv1(),
fElectronv2(),
- // fCurrentEventPosElectron(),
- // fCurrentEventNegElectron(),
- // fKFReconstructedGammasCut(),
- // fPreviousEventTLVNegElectron(),
- // fPreviousEventTLVPosElectron(),
+ fGammav1(),
+ fGammav2(),
+ fElectronRecalculatedv1(),
+ fElectronRecalculatedv2(),
fElectronMass(-1),
fGammaMass(-1),
fPi0Mass(-1),
fLowPtMapping(1.),
fHighPtMapping(3.),
fDoCF(kFALSE),
- fAODBranch(NULL),
- fAODBranchName("GammaConv")//,
- // fAODObjects(NULL)
+ fAODGamma(NULL),
+ fAODPi0(NULL),
+ fAODOmega(NULL),
+ fAODBranchName("GammaConv"),
+ fKFForceAOD(kFALSE),
+ fKFDeltaAODFileName(""),
+ fDoNeutralMesonV0MCCheck(kFALSE),
+ fUseTrackMultiplicityForBG(kTRUE),
+ fMoveParticleAccordingToVertex(kFALSE),
+ fApplyChi2Cut(kFALSE),
+ fNRandomEventsForBG(15),
+ fNDegreesPMBackground(15),
+ fDoRotation(kTRUE),
+ fCheckBGProbability(kTRUE),
+ fKFReconstructedGammasV0Index(),
+ fRemovePileUp(kFALSE),
+ fSelectV0AND(kFALSE),
+ fTriggerAnalysis(NULL),
+ fMultiplicity(0),
+ fUseMultiplicity(0),
+ fUseMultiplicityBin(0),
+ fUseCentrality(0),
+ fUseCentralityBin(0)
{
// Common I/O in slot 0
DefineInput (0, TChain::Class());
// Define standard ESD track cuts for Gamma-hadron correlation
SetESDtrackCuts();
+
}
AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion()
delete fEsdTrackCuts;
}
- if (fAODBranch) {
- fAODBranch->Clear();
- delete fAODBranch ;
+ //Delete AODs
+ if (fAODGamma) {
+ fAODGamma->Clear();
+ delete fAODGamma;
+ }
+ fAODGamma = NULL;
+
+ if (fAODPi0) {
+ fAODPi0->Clear();
+ delete fAODPi0;
+ }
+ fAODPi0 = NULL;
+
+ if (fAODOmega) {
+ fAODOmega->Clear();
+ delete fAODOmega;
+ }
+ fAODOmega = NULL;
+
+ if(fTriggerAnalysis) {
+ delete fTriggerAnalysis;
}
+
+
}
fEsdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
//standard cuts from:
//http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup
- //fEsdTrackCuts->SetMinNClustersTPC(50);
- //fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
- //fEsdTrackCuts->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);
- fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
- fEsdTrackCuts->SetRequireITSRefit(kTRUE);
- fEsdTrackCuts->SetMaxNsigmaToVertex(3);
- fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE);
+
+ // Cuts used up to 3rd of March
+
+ // fEsdTrackCuts->SetMinNClustersTPC(50);
+// fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
+// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
+// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
+// fEsdTrackCuts->SetMaxNsigmaToVertex(3);
+// fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE);
+
+ //------- To be tested-----------
+ // Cuts used up to 26th of Agost
+// Int_t minNClustersTPC = 70;
+// Double_t maxChi2PerClusterTPC = 4.0;
+// Double_t maxDCAtoVertexXY = 2.4; // cm
+// Double_t maxDCAtoVertexZ = 3.2; // cm
+// fEsdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
+// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
+// // fEsdTrackCuts->SetRequireTPCStandAlone(kTRUE);
+// fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
+// fEsdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
+// fEsdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
+// fEsdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
+// fEsdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
+// fEsdTrackCuts->SetDCAToVertex2D(kTRUE);
+// fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
+// fEsdTrackCuts->SetPtRange(0.15);
+
+// fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
+
+
+// Using standard function for setting Cuts
+ Bool_t selectPrimaries=kTRUE;
+ fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(selectPrimaries);
+ fEsdTrackCuts->SetMaxDCAToVertexZ(2);
+ fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
+ fEsdTrackCuts->SetPtRange(0.15);
+
+ //----- From Jacek 10.03.03 ------------------/
+// minNClustersTPC = 70;
+// maxChi2PerClusterTPC = 4.0;
+// maxDCAtoVertexXY = 2.4; // cm
+// maxDCAtoVertexZ = 3.2; // cm
+
+// esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+// esdTrackCuts->SetRequireTPCRefit(kFALSE);
+// esdTrackCuts->SetRequireTPCStandAlone(kTRUE);
+// esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
+// esdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
+// esdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
+// esdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
+// esdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
+// esdTrackCuts->SetDCAToVertex2D(kTRUE);
+
+
+
// fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
-
+ // fV0Reader->SetESDtrackCuts(fEsdTrackCuts);
}
-void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/)
+void AliAnalysisTaskGammaConversion::UserExec(Option_t */*option*/)
{
// Execute analysis for current event
-
- ConnectInputData("");
-
- //Each event needs an empty branch
- // fAODBranch->Clear();
- fAODBranch->Delete();
+
+ // Load the esdpid from the esdhandler if exists (tender was applied) otherwise set the Bethe Bloch parameters
+ Int_t eventQuality=-1;
+
+ AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+ AliESDInputHandler *esdHandler=0x0;
+ if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){
+ AliV0Reader::SetESDpid(esdHandler->GetESDpid());
+ } else {
+ //load esd pid bethe bloch parameters depending on the existance of the MC handler
+ // yes: MC parameters
+ // no: data parameters
+ if (!AliV0Reader::GetESDpid()){
+ if (fMCEvent ) {
+ AliV0Reader::InitESDpid();
+ } else {
+ AliV0Reader::InitESDpid(1);
+ }
+ }
+ }
+
+ //Must set fForceAOD to true for the AOD to get filled. Should only be done when running independent chain / train.
+ if(fKFForceAOD) {
+ if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
+ AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
+ }
+
+ if(fV0Reader == NULL){
+ // Write warning here cuts and so on are default if this ever happens
+ }
+
+ if (fMCEvent ) {
+ // To avoid crashes due to unzip errors. Sometimes the trees are not there.
+
+ AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
+ if (!mcHandler){
+ AliError("Could not retrive MC event handler!");
+ return;
+
+ eventQuality=0;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ }
+ if (!mcHandler->InitOk() ){
+ eventQuality=0;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+ if (!mcHandler->TreeK() ){
+ eventQuality=0;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+ if (!mcHandler->TreeTR() ) {
+ eventQuality=0;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+ }
+
+ fV0Reader->SetInputAndMCEvent(InputEvent(), MCEvent());
+
+ fV0Reader->Initialize();
+ fDoMCTruth = fV0Reader->GetDoMCTruth();
+
+ if(fAODGamma) fAODGamma->Delete();
+ if(fAODPi0) fAODPi0->Delete();
+ if(fAODOmega) fAODOmega->Delete();
if(fKFReconstructedGammasTClone == NULL){
fKFReconstructedGammasTClone = new TClonesArray("AliKFParticle",0);
if(fChargedParticles == NULL){
fChargedParticles = new TClonesArray("AliESDtrack",0);
}
-
+
+ if(fKFReconstructedPi0sTClone == NULL){
+ fKFReconstructedPi0sTClone = new TClonesArray("AliKFParticle",0);
+ }
+
+ if(fKFRecalculatedGammasTClone == NULL){
+ fKFRecalculatedGammasTClone = new TClonesArray("AliKFParticle",0);
+ }
+
+ if(fTriggerAnalysis== NULL){
+ fTriggerAnalysis = new AliTriggerAnalysis;
+ }
+
//clear TClones
fKFReconstructedGammasTClone->Delete();
fCurrentEventPosElectronTClone->Delete();
fKFReconstructedGammasCutTClone->Delete();
fPreviousEventTLVNegElectronTClone->Delete();
fPreviousEventTLVPosElectronTClone->Delete();
- //fKFReconstructedGammasTClone->Clear();
- //fCurrentEventPosElectronTClone->Clear();
- //fCurrentEventNegElectronTClone->Clear();
- //fKFReconstructedGammasCutTClone->Clear();
- //fPreviousEventTLVNegElectronTClone->Clear();
- //fPreviousEventTLVPosElectronTClone->Clear();
-
+ fKFReconstructedPi0sTClone->Delete();
+ fKFRecalculatedGammasTClone->Delete();
+
//clear vectors
- // fKFReconstructedGammas.clear();
fElectronv1.clear();
fElectronv2.clear();
- // fCurrentEventPosElectron.clear();
- // fCurrentEventNegElectron.clear();
- // fKFReconstructedGammasCut.clear();
+
+ fGammav1.clear();
+ fGammav2.clear();
+
+ fElectronRecalculatedv1.clear();
+ fElectronRecalculatedv2.clear();
+
fChargedParticles->Delete();
- //fChargedParticles->Clear();
+
fChargedParticlesId.clear();
+
+ fKFReconstructedGammasV0Index.clear();
//Clear the data in the v0Reader
- fV0Reader->UpdateEventByEventData();
-
-
+ // fV0Reader->UpdateEventByEventData();
+
+ //Take Only events with proper trigger
+ /*
+ if(fTriggerCINT1B){
+ if(!fV0Reader->GetESDEvent()->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL")) return;
+ }
+ */
+
+ if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
+ // cout<< "Event not taken"<< endl;
+ eventQuality=1;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ if(fDoMCTruth){
+ CheckMesonProcessTypeEventQuality(eventQuality);
+ }
+ return; // aborts if the primary vertex does not have contributors.
+ }
+
+
+ if(!fV0Reader->CheckForPrimaryVertexZ() ){
+ eventQuality=2;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ if(fDoMCTruth){
+ CheckMesonProcessTypeEventQuality(eventQuality);
+ }
+ return;
+ }
+
+ if(fRemovePileUp && fV0Reader->GetESDEvent()->IsPileupFromSPD()) {
+ eventQuality=4;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+
+ Bool_t v0A = fTriggerAnalysis->IsOfflineTriggerFired(fV0Reader->GetESDEvent(), AliTriggerAnalysis::kV0A);
+ Bool_t v0C = fTriggerAnalysis->IsOfflineTriggerFired(fV0Reader->GetESDEvent(), AliTriggerAnalysis::kV0C);
+ Bool_t v0AND = v0A && v0C;
+
+ if(fSelectV0AND && !v0AND){
+ eventQuality=5;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+ fMultiplicity = fEsdTrackCuts->CountAcceptedTracks(fV0Reader->GetESDEvent());
+
+ if( CalculateMultiplicityBin() != fUseMultiplicityBin){
+ eventQuality=6;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+
+ if(fV0Reader->GetIsHeavyIon()){
+ if(fUseCentrality>0){
+ AliESDCentrality *esdCentrality = fV0Reader->GetESDEvent()->GetCentrality();
+ Int_t centralityC = -1;
+
+ if(fUseCentrality==1){
+ centralityC = esdCentrality->GetCentralityClass10("V0M");
+ if( centralityC != fUseCentralityBin ){
+ eventQuality=7;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+ }
+
+ if(fUseCentrality==2){
+ centralityC = esdCentrality->GetCentralityClass10("CL1");
+ if( centralityC != fUseCentralityBin ){
+ eventQuality=7;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+ return;
+ }
+ }
+ }
+ }
+ eventQuality=3;
+ fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
+
+
+
+ fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",fMultiplicity);
+ if (fV0Reader->GetNumberOfContributorsVtx()>=1){
+ fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",fMultiplicity);
+ }
+
+
+
// Process the MC information
if(fDoMCTruth){
ProcessMCData();
//Process the v0 information with no cuts
ProcessV0sNoCut();
-
+
// Process the v0 information
ProcessV0s();
-
+
+
//Fill Gamma AOD
FillAODWithConversionGammas() ;
- //calculate background if flag is set
- if(fCalculateBackground){
- CalculateBackground();
- }
-
// Process reconstructed gammas
if(fDoNeutralMeson == kTRUE){
ProcessGammasForNeutralMesonAnalysis();
+
+ }
+
+ if(fDoMCTruth == kTRUE){
+ CheckV0Efficiency();
}
-
- CheckV0Efficiency();
-
//Process reconstructed gammas electrons for Chi_c Analysis
if(fDoChic == kTRUE){
ProcessGammaElectronsForChicAnalysis();
ProcessGammasForGammaJetAnalysis();
}
+ //calculate background if flag is set
+ if(fCalculateBackground){
+ CalculateBackground();
+ }
+
+ if(fDoNeutralMeson == kTRUE){
+// ProcessConvPHOSGammasForNeutralMesonAnalysis();
+ if(fDoOmegaMeson == kTRUE){
+ ProcessGammasForOmegaMesonAnalysis();
+ }
+ }
+
+ //Clear the data in the v0Reader
+ fV0Reader->UpdateEventByEventData();
+ if(fRecalculateV0ForGamma==kTRUE){
+ RecalculateV0ForGamma();
+ }
PostData(1, fOutputContainer);
PostData(2, fCFManager->GetParticleContainer()); // for CF
}
-void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t */*option*/){
- // see header file for documentation
-
- if(fV0Reader == NULL){
- // Write warning here cuts and so on are default if this ever happens
- }
- fV0Reader->Initialize();
-}
+// void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t *option){
+// // see header file for documentation
+// // printf(" ConnectInputData %s\n", GetName());
+
+// AliAnalysisTaskSE::ConnectInputData(option);
+
+// if(fV0Reader == NULL){
+// // Write warning here cuts and so on are default if this ever happens
+// }
+// fV0Reader->Initialize();
+// fDoMCTruth = fV0Reader->GetDoMCTruth();
+// }
+
+void AliAnalysisTaskGammaConversion::CheckMesonProcessTypeEventQuality(Int_t evtQ){
+ // Check meson process type event quality
+ fStack= MCEvent()->Stack();
+ fGCMCEvent=MCEvent();
+
+ for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
+ TParticle* particle = (TParticle *)fStack->Particle(iTracks);
+ if (!particle) {
+ //print warning here
+ continue;
+ }
+ if(particle->GetPdgCode()!=111){ //Pi0
+ continue;
+ }
+ if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
+ if(evtQ==1){
+ switch(GetProcessType(fGCMCEvent)){
+ case kProcSD:
+ fHistograms->FillHistogram("MC_SD_EvtQ1_Pi0_Pt", particle->Pt());
+ break;
+ case kProcDD:
+ fHistograms->FillHistogram("MC_DD_EvtQ1_Pi0_Pt", particle->Pt());
+ break;
+ case kProcND:
+ fHistograms->FillHistogram("MC_ND_EvtQ1_Pi0_Pt", particle->Pt());
+ break;
+ default:
+ AliError("Unknown Process");
+ }
+ }
+ if(evtQ==2){
+ switch(GetProcessType(fGCMCEvent)){
+ case kProcSD:
+ fHistograms->FillHistogram("MC_SD_EvtQ2_Pi0_Pt", particle->Pt());
+ break;
+ case kProcDD:
+ fHistograms->FillHistogram("MC_DD_EvtQ2_Pi0_Pt", particle->Pt());
+ break;
+ case kProcND:
+ fHistograms->FillHistogram("MC_ND_EvtQ2_Pi0_Pt", particle->Pt());
+ break;
+ default:
+ AliError("Unknown Process");
+ }
+ }
+ }
+
+}
void AliAnalysisTaskGammaConversion::ProcessMCData(){
// see header file for documentation
-
+ //InputEvent(), MCEvent());
+ /* TestAnaMarin
fStack = fV0Reader->GetMCStack();
fMCTruth = fV0Reader->GetMCTruth(); // for CF
fGCMCEvent = fV0Reader->GetMCEvent(); // for CF
-
+ */
+ fStack= MCEvent()->Stack();
+ fGCMCEvent=MCEvent();
// for CF
Double_t containerInput[3];
if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
return; // aborts if the primary vertex does not have contributors.
}
-
- for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
+ for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
+ // for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
TParticle* particle = (TParticle *)fStack->Particle(iTracks);
-
+
+
+
if (!particle) {
//print warning here
continue;
}
///////////////////////Begin Chic Analysis/////////////////////////////
-
- if(particle->GetPdgCode() == 443){//Is JPsi
- if(particle->GetNDaughters()==2){
- if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 &&
- TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){
- TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
- TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
- if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9)
- fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance
+ if(fDoChic) {
+ if(particle->GetPdgCode() == 443){//Is JPsi
+ if(particle->GetNDaughters()==2){
+ if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 &&
+ TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){
+
+ TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
+ TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
+ if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9)
+ fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance
- if( TMath::Abs(daug0->Eta()) < 0.9){
- if(daug0->GetPdgCode() == -11)
- fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
- else
- fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
+ if( TMath::Abs(daug0->Eta()) < 0.9){
+ if(daug0->GetPdgCode() == -11)
+ fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
+ else
+ fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
- }
- if(TMath::Abs(daug1->Eta()) < 0.9){
- if(daug1->GetPdgCode() == -11)
- fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
- else
- fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
+ }
+ if(TMath::Abs(daug1->Eta()) < 0.9){
+ if(daug1->GetPdgCode() == -11)
+ fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
+ else
+ fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
+ }
}
}
}
- }
- // const int CHI_C0 = 10441;
- // const int CHI_C1 = 20443;
- // const int CHI_C2 = 445
- if(particle->GetPdgCode() == 22){//gamma from JPsi
- if(particle->GetMother(0) > -1){
- if(fStack->Particle(particle->GetMother(0))->GetPdgCode() == 10441 ||
- fStack->Particle(particle->GetMother(0))->GetPdgCode() == 20443 ||
- fStack->Particle(particle->GetMother(0))->GetPdgCode() == 445){
- if(TMath::Abs(particle->Eta()) < 1.2)
- fHistograms->FillTable("Table_Electrons",17);// gamma from chic inside accptance
+ // const int CHI_C0 = 10441;
+ // const int CHI_C1 = 20443;
+ // const int CHI_C2 = 445
+ if(particle->GetPdgCode() == 22){//gamma from JPsi
+ if(particle->GetMother(0) > -1){
+ if(fStack->Particle(particle->GetMother(0))->GetPdgCode() == 10441 ||
+ fStack->Particle(particle->GetMother(0))->GetPdgCode() == 20443 ||
+ fStack->Particle(particle->GetMother(0))->GetPdgCode() == 445){
+ if(TMath::Abs(particle->Eta()) < 1.2)
+ fHistograms->FillTable("Table_Electrons",17);// gamma from chic inside accptance
+ }
}
}
- }
- if(particle->GetPdgCode() == 10441 || particle->GetPdgCode() == 20443 || particle->GetPdgCode() == 445){
- if( particle->GetNDaughters() == 2){
- TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
- TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
+ if(particle->GetPdgCode() == 10441 || particle->GetPdgCode() == 20443 || particle->GetPdgCode() == 445){
+ if( particle->GetNDaughters() == 2){
+ TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
+ TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
- if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){
- if( daug0->GetPdgCode() == 443){
- TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter());
- TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter());
- if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
- fHistograms->FillTable("Table_Electrons",18);
+ if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){
+ if( daug0->GetPdgCode() == 443){
+ TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter());
+ TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter());
+ if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
+ fHistograms->FillTable("Table_Electrons",18);
- }//if
- else if (daug1->GetPdgCode() == 443){
- TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter());
- TParticle* daugE1 = fStack->Particle(daug1->GetLastDaughter());
- if( TMath::Abs(daug0->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
- fHistograms->FillTable("Table_Electrons",18);
- }//else if
- }//gamma o Jpsi
- }//GetNDaughters
+ }//if
+ else if (daug1->GetPdgCode() == 443){
+ TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter());
+ TParticle* daugE1 = fStack->Particle(daug1->GetLastDaughter());
+ if( TMath::Abs(daug0->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
+ fHistograms->FillTable("Table_Electrons",18);
+ }//else if
+ }//gamma o Jpsi
+ }//GetNDaughters
+ }
}
-
/////////////////////End Chic Analysis////////////////////////////
- if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
+ // if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
if(particle->R()>fV0Reader->GetMaxRCut()) continue; // cuts on distance from collision point
rapidity = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())));
}
+
+
+ if(iTracks<=fStack->GetNprimary() ){
+ if ( particle->GetPdgCode()== -211 || particle->GetPdgCode()== 211 ||
+ particle->GetPdgCode()== 2212 || particle->GetPdgCode()==-2212 ||
+ particle->GetPdgCode()== 321 || particle->GetPdgCode()==-321 ){
+ if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
+ fHistograms->FillHistogram("MC_PhysicalPrimaryCharged_Pt", particle->Pt());
+ }
+ }
+
+
+
//process the gammas
if (particle->GetPdgCode() == 22){
-
+ if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
if(particle->GetMother(0) >-1 && fStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
continue; // no photon as mothers!
if(particle->GetMother(0) >= fStack->GetNprimary()){
continue; // the gamma has a mother, and it is not a primary particle
}
-
+
+ if(particle->GetMother(0) >-1){
+ fHistograms->FillHistogram("MC_DecayAllGamma_Pt", particle->Pt()); // All
+ switch(fStack->Particle(particle->GetMother(0))->GetPdgCode()){
+ case 111: // Pi0
+ fHistograms->FillHistogram("MC_DecayPi0Gamma_Pt", particle->Pt());
+ break;
+ case 113: // Rho0
+ fHistograms->FillHistogram("MC_DecayRho0Gamma_Pt", particle->Pt());
+ break;
+ case 221: // Eta
+ fHistograms->FillHistogram("MC_DecayEtaGamma_Pt", particle->Pt());
+ break;
+ case 223: // Omega
+ fHistograms->FillHistogram("MC_DecayOmegaGamma_Pt", particle->Pt());
+ break;
+ case 310: // K_s0
+ fHistograms->FillHistogram("MC_DecayK0sGamma_Pt", particle->Pt());
+ break;
+ case 331: // Eta'
+ fHistograms->FillHistogram("MC_DecayEtapGamma_Pt", particle->Pt());
+ break;
+ }
+ }
+
fHistograms->FillHistogram("MC_allGamma_Energy", particle->Energy());
fHistograms->FillHistogram("MC_allGamma_Pt", particle->Pt());
fHistograms->FillHistogram("MC_allGamma_Eta", particle->Eta());
Int_t rBin = fHistograms->GetRBin(ePos->R());
Int_t zBin = fHistograms->GetZBin(ePos->Vz());
Int_t phiBin = fHistograms->GetPhiBin(particle->Phi());
+ Double_t rFMD=30;
+ Double_t rITSTPCMin=50;
+ Double_t rITSTPCMax=80;
TVector3 vtxPos(ePos->Vx(),ePos->Vy(),ePos->Vz());
TString nameMCMappingPhiR="";
nameMCMappingPhiR.Form("MC_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
- fHistograms->FillHistogram(nameMCMappingPhiR, ePos->Vz(), particle->Eta());
+ // fHistograms->FillHistogram(nameMCMappingPhiR, ePos->Vz(), particle->Eta());
TString nameMCMappingPhi="";
nameMCMappingPhi.Form("MC_Conversion_Mapping_Phi%02d",phiBin);
// fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());
- fHistograms->FillHistogram(nameMCMappingPhi, ePos->Vz(), particle->Eta());
+ //fHistograms->FillHistogram(nameMCMappingPhi, ePos->Vz(), particle->Eta());
TString nameMCMappingR="";
nameMCMappingR.Form("MC_Conversion_Mapping_R%02d",rBin);
// fHistograms->FillHistogram(nameMCMappingR, particle->Eta());
- fHistograms->FillHistogram(nameMCMappingR,ePos->Vz(), particle->Eta());
+ //fHistograms->FillHistogram(nameMCMappingR,ePos->Vz(), particle->Eta());
TString nameMCMappingPhiInR="";
nameMCMappingPhiInR.Form("MC_Conversion_Mapping_Phi_in_R_%02d",rBin);
// fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
fHistograms->FillHistogram(nameMCMappingPhiInZ, vtxPos.Phi());
+
+ if(ePos->R()<rFMD){
+ TString nameMCMappingFMDPhiInZ="";
+ nameMCMappingFMDPhiInZ.Form("MC_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
+ fHistograms->FillHistogram(nameMCMappingFMDPhiInZ, vtxPos.Phi());
+ }
+
+ if(ePos->R()>rITSTPCMin && ePos->R()<rITSTPCMax){
+ TString nameMCMappingITSTPCPhiInZ="";
+ nameMCMappingITSTPCPhiInZ.Form("MC_Conversion_Mapping_ITSTPC_Phi_in_Z_%02d",zBin);
+ fHistograms->FillHistogram(nameMCMappingITSTPCPhiInZ, vtxPos.Phi());
+ }
+
TString nameMCMappingRInZ="";
nameMCMappingRInZ.Form("MC_Conversion_Mapping_R_in_Z_%02d",zBin);
fHistograms->FillHistogram(nameMCMappingRInZ,ePos->R() );
TString nameMCMappingMidPtPhiInZ="";
nameMCMappingMidPtPhiInZ.Form("MC_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
fHistograms->FillHistogram(nameMCMappingMidPtPhiInZ, vtxPos.Phi());
+
+
+ if(ePos->R()<rFMD){
+ TString nameMCMappingMidPtFMDPhiInZ="";
+ nameMCMappingMidPtFMDPhiInZ.Form("MC_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
+ fHistograms->FillHistogram(nameMCMappingMidPtFMDPhiInZ, vtxPos.Phi());
+ }
+
TString nameMCMappingMidPtRInZ="";
nameMCMappingMidPtRInZ.Form("MC_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());
if(daughter0->GetPdgCode() != 22 || daughter1->GetPdgCode() != 22) continue; //check for gamma gamma daughters
-
+
+ if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ) continue;
+
// Check the acceptance for both gammas
Bool_t gammaEtaCut = kTRUE;
if(TMath::Abs(daughter0->Eta()) > fV0Reader->GetEtaCut() || TMath::Abs(daughter1->Eta()) > fV0Reader->GetEtaCut() ) gammaEtaCut = kFALSE;
}
-
-
if(particle->GetPdgCode()==111){ //Pi0
if( iTracks >= fStack->GetNprimary()){
fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());
fHistograms->FillHistogram("MC_Pi0_Rapid", rapidity);
fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);
fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());
+ fHistograms->FillHistogram("MC_Pi0_Pt_vs_Rapid", particle->Pt(),rapidity);
fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());
fHistograms->FillHistogram("MC_Pi0_R", particle->R());
fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());
fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
-
+ if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_Fiducial", particle->Pt());
+
+ switch(GetProcessType(fGCMCEvent)){
+ case kProcSD:
+ fHistograms->FillHistogram("MC_SD_EvtQ3_Pi0_Pt", particle->Pt());
+ break;
+ case kProcDD:
+ fHistograms->FillHistogram("MC_DD_EvtQ3_Pi0_Pt", particle->Pt());
+ break;
+ case kProcND:
+ fHistograms->FillHistogram("MC_ND_EvtQ3_Pi0_Pt", particle->Pt());
+ break;
+ default:
+ AliError("Unknown Process");
+ }
+
if(gammaEtaCut && gammaRCut){
// if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
fHistograms->FillHistogram("MC_Pi0_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
+ if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_withinAcceptance_Fiducial", particle->Pt());
+
if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
fHistograms->FillHistogram("MC_Pi0_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
fHistograms->FillHistogram("MC_Pi0_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
+ fHistograms->FillHistogram("MC_Pi0_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
+ fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
+ fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
+
+ Double_t alfa=0.;
+ if((daughter0->Energy()+daughter1->Energy()) > 0.){
+ alfa= TMath::Abs((daughter0->Energy()-daughter1->Energy())/(daughter0->Energy()+daughter1->Energy()));
+ }
+ fHistograms->FillHistogram("MC_Pi0_alpha",alfa);
+ if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_ConvGamma_withinAcceptance_Fiducial", particle->Pt());
+
}
}
}
fHistograms->FillHistogram("MC_Eta_Rapid", rapidity);
fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);
fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());
+ fHistograms->FillHistogram("MC_Eta_Pt_vs_Rapid", particle->Pt(),rapidity);
fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());
fHistograms->FillHistogram("MC_Eta_R", particle->R());
fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());
fHistograms->FillHistogram("MC_Eta_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
fHistograms->FillHistogram("MC_Eta_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
fHistograms->FillHistogram("MC_Eta_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
+ fHistograms->FillHistogram("MC_Eta_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
+ fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
+ fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
+
}
}
/*AliESDv0 * cV0 = */fV0Reader->GetV0(i);
if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
- return;
+ continue;
}
- if( !fV0Reader->GetV0(i)->GetOnFlyStatus()){
- return;
+ // if( !fV0Reader->GetV0(i)->GetOnFlyStatus()){
+ if( !fV0Reader->CheckV0FinderStatus(i)){
+ continue;
}
- if( !(fV0Reader->GetNegativeESDTrack()->GetStatus() & AliESDtrack::kTPCrefit) ||
- !(fV0Reader->GetPositiveESDTrack()->GetStatus() & AliESDtrack::kTPCrefit) ){
- return;
+
+ if( !((fV0Reader->GetNegativeESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ||
+ !((fV0Reader->GetPositiveESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ){
+ continue;
+ }
+
+ if( fV0Reader->GetNegativeESDTrack()->GetSign()== fV0Reader->GetPositiveESDTrack()->GetSign()){
+ continue;
+ }
+
+ if( fV0Reader->GetNegativeESDTrack()->GetKinkIndex(0) > 0 ||
+ fV0Reader->GetPositiveESDTrack()->GetKinkIndex(0) > 0) {
+ continue;
+ }
+ if(TMath::Abs(fV0Reader->GetMotherCandidateEta())> fV0Reader->GetEtaCut()){
+ continue;
+ }
+ if(TMath::Abs(fV0Reader->GetPositiveTrackEta())> fV0Reader->GetEtaCut()){
+ continue;
+ }
+ if(TMath::Abs(fV0Reader->GetNegativeTrackEta())> fV0Reader->GetEtaCut()){
+ continue;
+ }
+ if((TMath::Abs(fV0Reader->GetZ())*fV0Reader->GetLineCutZRSlope())-fV0Reader->GetLineCutZValue() > fV0Reader->GetXYRadius() ){ // cuts out regions where we do not reconstruct
+ continue;
}
-
if(fDoMCTruth){
if(fV0Reader->HasSameMCMother() == kFALSE){
if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
continue;
}
- if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
+ if(negativeMC->GetPdgCode() == positiveMC->GetPdgCode()){
continue;
}
- if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){
+ if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){ // id 5 is conversion
continue;
}
void AliAnalysisTaskGammaConversion::ProcessV0s(){
// see header file for documentation
-
+
+
if(fWriteNtuple == kTRUE){
FillNtuple();
}
Int_t nSurvivingV0s=0;
+ fV0Reader->ResetNGoodV0s();
while(fV0Reader->NextV0()){
nSurvivingV0s++;
-
+
+ TVector3 vtxConv(fV0Reader->GetX(),fV0Reader->GetY(), fV0Reader->GetZ());
+
//-------------------------- filling v0 information -------------------------------------
fHistograms->FillHistogram("ESD_Conversion_R", fV0Reader->GetXYRadius());
fHistograms->FillHistogram("ESD_Conversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
fHistograms->FillHistogram("ESD_Conversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
fHistograms->FillHistogram("ESD_Conversion_OpeningAngle", fV0Reader->GetOpeningAngle());
-
+
+ // Specific histograms for beam pipe studies
+ if( TMath::Abs(fV0Reader->GetZ()) < fV0Reader->GetLineCutZValue() ){
+ fHistograms->FillHistogram("ESD_Conversion_XY_BeamPipe", fV0Reader->GetX(),fV0Reader->GetY());
+ fHistograms->FillHistogram("ESD_Conversion_RPhi_BeamPipe", vtxConv.Phi(),fV0Reader->GetXYRadius());
+ }
+
+
fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());
fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());
fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());
fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());
-
+ fHistograms->FillHistogram("ESD_E_nTPCClusters", fV0Reader->GetNegativeTracknTPCClusters());
+ fHistograms->FillHistogram("ESD_E_nITSClusters", fV0Reader->GetNegativeTracknITSClusters());
+ if(fV0Reader->GetNegativeTracknTPCFClusters()!=0 && fV0Reader->GetNegativeTracknTPCClusters()!=0 ){
+ Double_t eClsToF= (Double_t)fV0Reader->GetNegativeTracknTPCClusters()/(Double_t)fV0Reader->GetNegativeTracknTPCFClusters();
+ fHistograms->FillHistogram("ESD_E_nTPCClustersToFP", fV0Reader->GetNegativeTrackP(),eClsToF );
+ fHistograms->FillHistogram("ESD_E_TPCchi2", fV0Reader->GetNegativeTrackTPCchi2()/(Double_t)fV0Reader->GetNegativeTracknTPCClusters());
+ }
+
+
+
fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());
fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());
fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());
fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());
-
+ fHistograms->FillHistogram("ESD_P_nTPCClusters", fV0Reader->GetPositiveTracknTPCClusters());
+ fHistograms->FillHistogram("ESD_P_nITSClusters", fV0Reader->GetPositiveTracknITSClusters());
+ if(fV0Reader->GetPositiveTracknTPCFClusters()!=0 && (Double_t)fV0Reader->GetPositiveTracknTPCClusters()!=0 ){
+ Double_t pClsToF= (Double_t)fV0Reader->GetPositiveTracknTPCClusters()/(Double_t)fV0Reader->GetPositiveTracknTPCFClusters();
+ fHistograms->FillHistogram("ESD_P_nTPCClustersToFP",fV0Reader->GetPositiveTrackP(), pClsToF);
+ fHistograms->FillHistogram("ESD_P_TPCchi2", fV0Reader->GetPositiveTrackTPCchi2()/(Double_t)fV0Reader->GetPositiveTracknTPCClusters());
+ }
+
+
fHistograms->FillHistogram("ESD_ConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
fHistograms->FillHistogram("ESD_ConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
fHistograms->FillHistogram("ESD_ConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
fHistograms->FillHistogram("ESD_ConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
fHistograms->FillHistogram("ESD_ConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
-
-
+ Double_t negPID=0;
+ Double_t posPID=0;
+ fV0Reader->GetPIDProbability(negPID,posPID);
+ fHistograms->FillHistogram("ESD_ConvGamma_E_EProbP",fV0Reader->GetNegativeTrackP(),negPID);
+ fHistograms->FillHistogram("ESD_ConvGamma_P_EProbP",fV0Reader->GetPositiveTrackP(),posPID);
+
+ Double_t negPIDmupi=0;
+ Double_t posPIDmupi=0;
+ fV0Reader->GetPIDProbabilityMuonPion(negPIDmupi,posPIDmupi);
+ fHistograms->FillHistogram("ESD_ConvGamma_E_mupiProbP",fV0Reader->GetNegativeTrackP(),negPIDmupi);
+ fHistograms->FillHistogram("ESD_ConvGamma_P_mupiProbP",fV0Reader->GetPositiveTrackP(),posPIDmupi);
+
+ Double_t armenterosQtAlfa[2];
+ fV0Reader->GetArmenterosQtAlfa(fV0Reader-> GetNegativeKFParticle(),
+ fV0Reader-> GetPositiveKFParticle(),
+ fV0Reader->GetMotherCandidateKFCombination(),
+ armenterosQtAlfa);
+
+ fHistograms->FillHistogram("ESD_ConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
+
+
// begin mapping
Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());
Int_t zBin = fHistograms->GetZBin(fV0Reader->GetZ());
Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());
- TVector3 vtxConv(fV0Reader->GetX(),fV0Reader->GetY(), fV0Reader->GetZ());
+ Double_t rFMD=30;
+ Double_t rITSTPCMin=50;
+ Double_t rITSTPCMax=80;
- Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();
+
+ // Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();
TString nameESDMappingPhiR="";
nameESDMappingPhiR.Form("ESD_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
- fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);
+ //fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);
TString nameESDMappingPhi="";
nameESDMappingPhi.Form("ESD_Conversion_Mapping_Phi%02d",phiBin);
- fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);
+ //fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);
TString nameESDMappingR="";
nameESDMappingR.Form("ESD_Conversion_Mapping_R%02d",rBin);
- fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta);
+ //fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta);
TString nameESDMappingPhiInR="";
nameESDMappingPhiInR.Form("ESD_Conversion_Mapping_Phi_in_R_%02d",rBin);
// fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
fHistograms->FillHistogram(nameESDMappingPhiInZ, vtxConv.Phi());
+ if(fV0Reader->GetXYRadius()<rFMD){
+ TString nameESDMappingFMDPhiInZ="";
+ nameESDMappingFMDPhiInZ.Form("ESD_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
+ fHistograms->FillHistogram(nameESDMappingFMDPhiInZ, vtxConv.Phi());
+ }
+
+ if(fV0Reader->GetXYRadius()>rITSTPCMin && fV0Reader->GetXYRadius()<rITSTPCMax){
+ TString nameESDMappingITSTPCPhiInZ="";
+ nameESDMappingITSTPCPhiInZ.Form("ESD_Conversion_Mapping_ITSTPC_Phi_in_Z_%02d",zBin);
+ fHistograms->FillHistogram(nameESDMappingITSTPCPhiInZ, vtxConv.Phi());
+ }
+
TString nameESDMappingRInZ="";
nameESDMappingRInZ.Form("ESD_Conversion_Mapping_R_in_Z_%02d",zBin);
fHistograms->FillHistogram(nameESDMappingRInZ, fV0Reader->GetXYRadius());
TString nameESDMappingMidPtPhiInZ="";
nameESDMappingMidPtPhiInZ.Form("ESD_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
fHistograms->FillHistogram(nameESDMappingMidPtPhiInZ, vtxConv.Phi());
+ if(fV0Reader->GetXYRadius()<rFMD){
+ TString nameESDMappingMidPtFMDPhiInZ="";
+ nameESDMappingMidPtFMDPhiInZ.Form("ESD_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
+ fHistograms->FillHistogram(nameESDMappingMidPtFMDPhiInZ, vtxConv.Phi());
+ }
+
TString nameESDMappingMidPtRInZ="";
nameESDMappingMidPtRInZ.Form("ESD_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
// end mapping
new((*fKFReconstructedGammasTClone)[fKFReconstructedGammasTClone->GetEntriesFast()]) AliKFParticle(*fV0Reader->GetMotherCandidateKFCombination());
-
+ fKFReconstructedGammasV0Index.push_back(fV0Reader->GetCurrentV0IndexNumber()-1);
// fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());
fElectronv1.push_back(fV0Reader->GetCurrentV0()->GetPindex());
fElectronv2.push_back(fV0Reader->GetCurrentV0()->GetNindex());
-
//----------------------------------- checking for "real" conversions (MC match) --------------------------------------
if(fDoMCTruth){
+ TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
+ TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
if(fV0Reader->HasSameMCMother() == kFALSE){
+ fHistograms->FillHistogram("ESD_TrueConvCombinatorial_R", fV0Reader->GetXYRadius());
+ if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
+ fHistograms->FillHistogram("ESD_TrueConvCombinatorialElec_R", fV0Reader->GetXYRadius());
+ }
continue;
}
- TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
- TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
+ // Moved up to check true electron background
+ // TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
+ // TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
continue;
if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
continue;
}
+ if( (negativeMC->GetUniqueID() == 4 && positiveMC->GetUniqueID() ==4) ||
+ (negativeMC->GetUniqueID() == 0 && positiveMC->GetUniqueID() ==0) ){// fill r distribution for Dalitz decays
+ if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 111){ //pi0
+ fHistograms->FillHistogram("ESD_TrueDalitzContamination_R", fV0Reader->GetXYRadius());
+ fHistograms->FillHistogram("ESD_TrueConvDalitzPi0_R", fV0Reader->GetXYRadius());
+ }
+ if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 221){ //eta
+ fHistograms->FillHistogram("ESD_TrueConvDalitzEta_R", fV0Reader->GetXYRadius());
+ }
+
+ }
if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){// check if the daughters come from a conversion
continue;
if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
+ if(fDoCF){
+ Double_t containerInput[3];
+ containerInput[0] = fV0Reader->GetMotherCandidatePt();
+ containerInput[1] = fV0Reader->GetMotherCandidateEta();
+ containerInput[2] = fV0Reader->GetMotherCandidateMass();
+ fCFManager->GetParticleContainer()->Fill(containerInput,kStepTrueGamma); // for CF
+ }
+
fHistograms->FillHistogram("ESD_TrueConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
fHistograms->FillHistogram("ESD_TrueConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
fHistograms->FillHistogram("ESD_TrueConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
fHistograms->FillHistogram("ESD_TrueConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
fHistograms->FillHistogram("ESD_TrueConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
fHistograms->FillHistogram("ESD_TrueConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
-
- fHistograms->FillHistogram("ESD_TrueConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
- fHistograms->FillHistogram("ESD_TrueConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
+ if (fV0Reader->GetMotherCandidateP() != 0) {
+ fHistograms->FillHistogram("ESD_TrueConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
+ fHistograms->FillHistogram("ESD_TrueConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
+ } else { cout << "Error::fV0Reader->GetNegativeTrackP() == 0 !!!" << endl; }
+
fHistograms->FillHistogram("ESD_TrueConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
fHistograms->FillHistogram("ESD_TrueConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
-
+ fHistograms->FillHistogram("ESD_TrueConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
+
//store MCTruth properties
fHistograms->FillHistogram("ESD_TrueConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
fHistograms->FillHistogram("ESD_TrueConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
fHistograms->FillHistogram("ESD_TrueConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
-
+
//resolution
Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();
Double_t esdpt = fV0Reader->GetMotherCandidatePt();
- Double_t resdPt = 0;
+ Double_t resdPt = 0.;
if(mcpt > 0){
- resdPt = ((esdpt - mcpt)/mcpt)*100;
+ resdPt = ((esdpt - mcpt)/mcpt)*100.;
}
else if(mcpt < 0){
cout<<"Pt of MC particle is negative, this will cause wrong calculation of resPt"<<endl;
}
- fHistograms->FillHistogram("Resolution_dPt", mcpt, resdPt);
+ fHistograms->FillHistogram("Resolution_Gamma_dPt_Pt", mcpt, resdPt);
fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);
fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);
+ fHistograms->FillHistogram("Resolution_Gamma_dPt_Phi", fV0Reader->GetMotherCandidatePhi(), resdPt);
- Double_t resdZ = 0;
+ Double_t resdZ = 0.;
if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){
- resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100;
+ resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100.;
}
-
+ Double_t resdZAbs = 0.;
+ resdZAbs = (fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz());
+
+ fHistograms->FillHistogram("Resolution_dZAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdZAbs);
fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);
fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());
fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());
-
- Double_t resdR = 0;
+
+ // new for dPt_Pt-histograms for Electron and Positron
+ Double_t mcEpt = fV0Reader->GetNegativeMCParticle()->Pt();
+ Double_t resEdPt = 0.;
+ if (mcEpt > 0){
+ resEdPt = ((fV0Reader->GetNegativeTrackPt()-mcEpt)/mcEpt)*100.;
+ }
+ UInt_t statusN = fV0Reader->GetNegativeESDTrack()->GetStatus();
+ // AliESDtrack * negTrk = fV0Reader->GetNegativeESDTrack();
+ UInt_t kTRDoutN = (statusN & AliESDtrack::kTRDout);
+
+ Int_t nITSclsE= fV0Reader->GetNegativeTracknITSClusters();
+ // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
+ switch(nITSclsE){
+ case 0: // 0 ITS clusters
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS0", mcEpt, resEdPt);
+ break;
+ case 1: // 1 ITS cluster
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS1", mcEpt, resEdPt);
+ break;
+ case 2: // 2 ITS clusters
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS2", mcEpt, resEdPt);
+ break;
+ case 3: // 3 ITS clusters
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS3", mcEpt, resEdPt);
+ break;
+ case 4: // 4 ITS clusters
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS4", mcEpt, resEdPt);
+ break;
+ case 5: // 5 ITS clusters
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS5", mcEpt, resEdPt);
+ break;
+ case 6: // 6 ITS clusters
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS6", mcEpt, resEdPt);
+ break;
+ }
+ //Filling histograms with respect to Electron resolution
+ fHistograms->FillHistogram("Resolution_E_dPt_Pt", mcEpt, resEdPt);
+ fHistograms->FillHistogram("Resolution_E_dPt_Phi", fV0Reader->GetNegativeTrackPhi(), resEdPt);
+ if(kTRDoutN){
+ fHistograms->FillHistogram("Resolution_E_nTRDtracklets_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
+ fHistograms->FillHistogram("Resolution_E_nTRDtracklets_MCPt", mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
+ fHistograms->FillHistogram("Resolution_E_nTRDclusters_ESDPt",fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDncls());
+ fHistograms->FillHistogram("Resolution_E_nTRDclusters_MCPt",mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDncls());
+ fHistograms->FillHistogram("Resolution_E_TRDsignal_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDsignal());
+ }
+
+ Double_t mcPpt = fV0Reader->GetPositiveMCParticle()->Pt();
+ Double_t resPdPt = 0;
+ if (mcPpt > 0){
+ resPdPt = ((fV0Reader->GetPositiveTrackPt()-mcPpt)/mcPpt)*100.;
+ }
+
+ UInt_t statusP = fV0Reader->GetPositiveESDTrack()->GetStatus();
+// AliESDtrack * posTr= fV0Reader->GetPositiveESDTrack();
+ UInt_t kTRDoutP = (statusP & AliESDtrack::kTRDout);
+
+ Int_t nITSclsP = fV0Reader->GetPositiveTracknITSClusters();
+ // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
+ switch(nITSclsP){
+ case 0: // 0 ITS clusters
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS0", mcPpt, resPdPt);
+ break;
+ case 1: // 1 ITS cluster
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS1", mcPpt, resPdPt);
+ break;
+ case 2: // 2 ITS clusters
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS2", mcPpt, resPdPt);
+ break;
+ case 3: // 3 ITS clusters
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS3", mcPpt, resPdPt);
+ break;
+ case 4: // 4 ITS clusters
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS4", mcPpt, resPdPt);
+ break;
+ case 5: // 5 ITS clusters
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS5", mcPpt, resPdPt);
+ break;
+ case 6: // 6 ITS clusters
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS6", mcPpt, resPdPt);
+ break;
+ }
+ //Filling histograms with respect to Positron resolution
+ fHistograms->FillHistogram("Resolution_P_dPt_Pt", mcPpt, resPdPt);
+ fHistograms->FillHistogram("Resolution_P_dPt_Phi", fV0Reader->GetPositiveTrackPhi(), resPdPt);
+ if(kTRDoutP){
+ fHistograms->FillHistogram("Resolution_P_nTRDtracklets_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
+ fHistograms->FillHistogram("Resolution_P_nTRDtracklets_MCPt", mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
+ fHistograms->FillHistogram("Resolution_P_nTRDclusters_ESDPt",fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDncls());
+ fHistograms->FillHistogram("Resolution_P_nTRDclusters_MCPt",mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDncls());
+ fHistograms->FillHistogram("Resolution_P_TRDsignal_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDsignal());
+ }
+
+
+ Double_t resdR = 0.;
if(fV0Reader->GetNegativeMCParticle()->R() != 0){
- resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100;
+ resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100.;
}
-
+ Double_t resdRAbs = 0.;
+ resdRAbs = (fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R());
+
+ fHistograms->FillHistogram("Resolution_dRAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdRAbs);
fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);
fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());
fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());
- fHistograms->FillHistogram("Resolution_dR_dPt", resdR, resdPt);
+ fHistograms->FillHistogram("Resolution_R_dPt", fV0Reader->GetNegativeMCParticle()->R(), resdPt);
+
+ Double_t resdPhiAbs=0.;
+ resdPhiAbs=0.;
+ resdPhiAbs= (fV0Reader->GetMotherCandidatePhi()-fV0Reader->GetNegativeMCParticle()->Phi());
+ fHistograms->FillHistogram("Resolution_dPhiAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdPhiAbs);
+
}//if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22)
}//if(fDoMCTruth)
}//while(fV0Reader->NextV0)
fHistograms->FillHistogram("ESD_NumberOfSurvivingV0s", nSurvivingV0s);
fHistograms->FillHistogram("ESD_NumberOfV0s", fV0Reader->GetNumberOfV0s());
+ fHistograms->FillHistogram("ESD_NumberOfContributorsVtx", fV0Reader->GetNumberOfContributorsVtx());
+
+ fV0Reader->ResetV0IndexNumber();
}
void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){
for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
// for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
//Create AOD particle object from AliKFParticle
-
- /* AliKFParticle * gammakf = &fKFReconstructedGammas[gammaIndex];
//You could add directly AliKFParticle objects to the AOD, avoiding dependences with PartCorr
//but this means that I have to work a little bit more in my side.
//AODPWG4Particle objects are simpler and lighter, I think
+ /*
+ AliKFParticle * gammakf = dynamic_cast<AliKFParticle*>(fKFReconstructedGammasTClone->At(gammaIndex));
AliAODPWG4Particle gamma = AliAODPWG4Particle(gammakf->Px(),gammakf->Py(),gammakf->Pz(), gammakf->E());
- gamma.SetLabel(-1);//How to get the MC label of the reconstructed gamma?
- gamma.SetCaloLabel(-1,-1); //How to get the MC label of the 2 electrons that form the gamma?
+ //gamma.SetLabel(-1);//How to get the MC label of the reconstructed gamma?
+ gamma.SetTrackLabel( fElectronv1[gammaIndex], fElectronv2[gammaIndex] ); //How to get the MC label of the 2 electrons that form the gamma?
gamma.SetDetector("CTS"); //tag the gamma as reconstructed in the central barrel
- gamma.SetPdg(AliCaloPID::kPhotonConv); //photon id
+ gamma.SetPdg(AliPID::kEleCon); //photon id
gamma.SetTag(-1); //Here I usually put a flag saying that montecarlo says it is prompt, decay fragmentation photon, or hadrons or whatever
-
- //Add it to the aod list
+ gamma.SetChi2(gammakf->Chi2());
Int_t i = fAODBranch->GetEntriesFast();
new((*fAODBranch)[i]) AliAODPWG4Particle(gamma);
*/
- // AliKFParticle * gammakf = &fKFReconstructedGammas[gammaIndex];
+
AliKFParticle * gammakf = (AliKFParticle *)fKFReconstructedGammasTClone->At(gammaIndex);
AliGammaConversionAODObject aodObject;
aodObject.SetPx(gammakf->GetPx());
aodObject.SetPz(gammakf->GetPz());
aodObject.SetLabel1(fElectronv1[gammaIndex]);
aodObject.SetLabel2(fElectronv2[gammaIndex]);
- Int_t i = fAODBranch->GetEntriesFast();
- new((*fAODBranch)[i]) AliGammaConversionAODObject(aodObject);
+ aodObject.SetChi2(gammakf->Chi2());
+ aodObject.SetE(gammakf->E());
+ Int_t i = fAODGamma->GetEntriesFast();
+ new((*fAODGamma)[i]) AliGammaConversionAODObject(aodObject);
}
}
+void AliAnalysisTaskGammaConversion::ProcessGammasForOmegaMesonAnalysis(){
+ // omega meson analysis pi0+gamma decay
+ for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
+ AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
+ for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
+
+ AliKFParticle * omegaCandidateGammaDaughter = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
+ if(fGammav1[firstPi0Index]==firstGammaIndex || fGammav2[firstPi0Index]==firstGammaIndex){
+ continue;
+ }
+
+ AliKFParticle omegaCandidate(*omegaCandidatePi0Daughter,*omegaCandidateGammaDaughter);
+ Double_t massOmegaCandidate = 0.;
+ Double_t widthOmegaCandidate = 0.;
+
+ omegaCandidate.GetMass(massOmegaCandidate,widthOmegaCandidate);
+
+ if ( massOmegaCandidate > 733 && massOmegaCandidate < 833 ) {
+ AddOmegaToAOD(&omegaCandidate, massOmegaCandidate, firstPi0Index, firstGammaIndex);
+ }
+
+ fHistograms->FillHistogram("ESD_Omega_InvMass_vs_Pt",massOmegaCandidate ,omegaCandidate.GetPt());
+ fHistograms->FillHistogram("ESD_Omega_InvMass",massOmegaCandidate);
+
+ //delete omegaCandidate;
+
+ }// end of omega reconstruction in pi0+gamma channel
+
+ if(fDoJet == kTRUE){
+ AliKFParticle* negPiKF=NULL;
+ AliKFParticle* posPiKF=NULL;
+
+ // look at the pi+pi+pi0 channel
+ for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
+ AliESDtrack* posTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
+ if (posTrack->GetSign()<0) continue;
+ if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion))>2.) continue;
+ if (posPiKF) delete posPiKF; posPiKF=NULL;
+ posPiKF = new AliKFParticle( *(posTrack) ,211);
+
+ for(Int_t jCh=0;jCh<fChargedParticles->GetEntriesFast();jCh++){
+ AliESDtrack* negTrack = (AliESDtrack*)(fChargedParticles->At(jCh));
+ if( negTrack->GetSign()>0) continue;
+ if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion))>2.) continue;
+ if (negPiKF) delete negPiKF; negPiKF=NULL;
+ negPiKF = new AliKFParticle( *(negTrack) ,-211);
+ AliKFParticle omegaCandidatePipPinPi0(*omegaCandidatePi0Daughter,*posPiKF,*negPiKF);
+ Double_t massOmegaCandidatePipPinPi0 = 0.;
+ Double_t widthOmegaCandidatePipPinPi0 = 0.;
+
+ omegaCandidatePipPinPi0.GetMass(massOmegaCandidatePipPinPi0,widthOmegaCandidatePipPinPi0);
+
+ if ( massOmegaCandidatePipPinPi0 > 733 && massOmegaCandidatePipPinPi0 < 833 ) {
+ AddOmegaToAOD(&omegaCandidatePipPinPi0, massOmegaCandidatePipPinPi0, -1, -1);
+ }
+
+ fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass_vs_Pt",massOmegaCandidatePipPinPi0 ,omegaCandidatePipPinPi0.GetPt());
+ fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass",massOmegaCandidatePipPinPi0);
+
+ // delete omegaCandidatePipPinPi0;
+ }
+ }
+
+ if (posPiKF) delete posPiKF; posPiKF=NULL; if (negPiKF) delete negPiKF; negPiKF=NULL;
+
+ } // checking ig gammajet because in that case the chargedparticle list is created
+
+ }
+
+ if(fCalculateBackground){
+
+ AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
+
+ Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
+ Int_t mbin = 0;
+ if(fUseTrackMultiplicityForBG == kTRUE){
+ mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
+ }
+ else{
+ mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
+ }
+
+ AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
+
+ // Background calculation for the omega
+ for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
+ AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);
+
+ if(fMoveParticleAccordingToVertex == kTRUE){
+ bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
+ }
+ for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
+ AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
+
+ if(fMoveParticleAccordingToVertex == kTRUE){
+ MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
+ }
+
+ for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
+ AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
+ AliKFParticle * omegaBckCandidate = new AliKFParticle(*omegaCandidatePi0Daughter,previousGoodV0);
+ Double_t massOmegaBckCandidate = 0.;
+ Double_t widthOmegaBckCandidate = 0.;
+
+ omegaBckCandidate->GetMass(massOmegaBckCandidate,widthOmegaBckCandidate);
+
+
+ fHistograms->FillHistogram("ESD_Omega_Bck_InvMass_vs_Pt",massOmegaBckCandidate ,omegaBckCandidate->GetPt());
+ fHistograms->FillHistogram("ESD_Omega_Bck_InvMass",massOmegaBckCandidate);
+
+ delete omegaBckCandidate;
+ }
+ }
+ }
+ } // end of checking if background calculation is available
+}
+
+
+void AliAnalysisTaskGammaConversion::AddOmegaToAOD(const AliKFParticle * const omegakf, Double_t mass, Int_t omegaDaughter, Int_t gammaDaughter) {
+ //See header file for documentation
+ AliGammaConversionAODObject omega;
+ omega.SetPx(omegakf->GetPx());
+ omega.SetPy(omegakf->GetPy());
+ omega.SetPz(omegakf->GetPz());
+ omega.SetChi2(omegakf->GetChi2());
+ omega.SetE(omegakf->GetE());
+ omega.SetIMass(mass);
+ omega.SetLabel1(omegaDaughter);
+ //dynamic_cast<AliGammaConversionAODObject*>(fAODBranch->At(daughter1))->SetTagged(kTRUE);
+ omega.SetLabel2(gammaDaughter);
+ new((*fAODOmega)[fAODOmega->GetEntriesFast()]) AliGammaConversionAODObject(omega);
+}
+
+
+
+void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
+ // see header file for documentation
+
+ // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
+ // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
+
+ fESDEvent = fV0Reader->GetESDEvent();
+
+ if(fKFReconstructedGammasTClone->GetEntriesFast()>fV0Reader->GetNumberOfV0s()){
+ cout<<"Warning, number of entries in the tclone is bigger than number of v0s"<<endl;
+ }
+
+ for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
+ for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
+
+ // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
+ // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
+
+ AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
+ AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
+
+ if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
+ continue;
+ }
+ if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
+ continue;
+ }
+
+ AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
+
+ Double_t massTwoGammaCandidate = 0.;
+ Double_t widthTwoGammaCandidate = 0.;
+ Double_t chi2TwoGammaCandidate =10000.;
+ twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
+ // if(twoGammaCandidate->GetNDF()>0){
+ // chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
+ chi2TwoGammaCandidate = twoGammaCandidate->GetChi2();
+
+ fHistograms->FillHistogram("ESD_Mother_Chi2",chi2TwoGammaCandidate);
+ if((chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
+
+ TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
+ TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
+
+ Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
+ Double_t rapidity;
+ if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() <= 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() <= 0){
+ cout << "Error: |Pz| > E !!!! " << endl;
+ rapidity=0;
+ }
+ else{
+ rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
+ }
+
+ if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut()){
+ delete twoGammaCandidate;
+ continue; // rapidity cut
+ }
+
+
+ Double_t alfa=0.0;
+ if( (twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()) != 0){
+ alfa=TMath::Abs((twoGammaDecayCandidateDaughter0->GetE()-twoGammaDecayCandidateDaughter1->GetE())
+ /(twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()));
+ }
+
+ if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
+ delete twoGammaCandidate;
+ continue; // minimum opening angle to avoid using ghosttracks
+ }
+
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
+ fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
+ fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
+ fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
+ fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_Mother_alfa", alfa);
+ if(massTwoGammaCandidate>0.1 && massTwoGammaCandidate<0.15){
+ fHistograms->FillHistogram("ESD_Mother_alfa_Pi0", alfa);
+ }
+ fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
+ fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ }
+ if(alfa<0.1){
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_E_alpha",massTwoGammaCandidate ,twoGammaCandidate->GetE());
+ }
+
+
+ if(fCalculateBackground){
+ /* Kenneth, just for testing*/
+ AliGammaConversionBGHandler * bgHandlerTest = fV0Reader->GetBGHandler();
+
+ Int_t zbin= bgHandlerTest->GetZBinIndex(fV0Reader->GetVertexZ());
+ Int_t mbin=0;
+ Int_t multKAA=0;
+ if(fUseTrackMultiplicityForBG == kTRUE){
+ multKAA=fV0Reader->CountESDTracks();
+ mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
+ }
+ else{// means we use #v0s for multiplicity
+ multKAA=fV0Reader->GetNGoodV0s();
+ mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
+ }
+ // cout<<"Filling bin number "<<zbin<<" and "<<mbin<<endl;
+ // cout<<"Corresponding to z = "<<fV0Reader->GetVertexZ()<<" and m = "<<multKAA<<endl;
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass",zbin,mbin),massTwoGammaCandidate);
+ fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass_vs_Pt",zbin,mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ /* end Kenneth, just for testing*/
+ fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ }
+ }
+ /* if(fCalculateBackground){
+ AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
+ Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
+ fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ }*/
+ // if(fDoNeutralMesonV0MCCheck){
+ if(fDoMCTruth){
+ //Kenneth: Checking the eta of the gamma to check the difference between 0.9 and 1.2
+ Int_t indexKF1 = fKFReconstructedGammasV0Index.at(firstGammaIndex);
+ if(indexKF1<fV0Reader->GetNumberOfV0s()){
+ fV0Reader->GetV0(indexKF1);//updates to the correct v0
+ Double_t eta1 = fV0Reader->GetMotherCandidateEta();
+ Bool_t isRealPi0=kFALSE;
+ Bool_t isRealEta=kFALSE;
+ Int_t gamma1MotherLabel=-1;
+ if(fV0Reader->HasSameMCMother() == kTRUE){
+ //cout<<"This v0 is a real v0!!!!"<<endl;
+ TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
+ TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
+ if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
+ if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
+ if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
+ gamma1MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
+ }
+ }
+ }
+ }
+ Int_t indexKF2 = fKFReconstructedGammasV0Index.at(secondGammaIndex);
+ if(indexKF1 == indexKF2){
+ cout<<"index of the two KF particles are the same.... should not happen"<<endl;
+ }
+ if(indexKF2<fV0Reader->GetNumberOfV0s()){
+ fV0Reader->GetV0(indexKF2);
+ Double_t eta2 = fV0Reader->GetMotherCandidateEta();
+ Int_t gamma2MotherLabel=-1;
+ if(fV0Reader->HasSameMCMother() == kTRUE){
+ TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
+ TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
+ if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
+ if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
+ if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
+ gamma2MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
+ }
+ }
+ }
+ }
+ if(gamma1MotherLabel>=0 && gamma1MotherLabel==gamma2MotherLabel){
+ if(fV0Reader->CheckIfPi0IsMother(gamma1MotherLabel)){
+ isRealPi0=kTRUE;
+ }
+ if(fV0Reader->CheckIfEtaIsMother(gamma1MotherLabel)){
+ isRealEta=kTRUE;
+ }
+
+ }
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ if(TMath::Abs(eta1)>0.9 && TMath::Abs(eta2)>0.9){
+ // fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
+ // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ if(isRealPi0 || isRealEta){
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_1212",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_1212",openingAngleTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ }
+
+ if(!isRealPi0 && !isRealEta){
+ if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
+ fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
+ }else{
+ fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
+ }
+ }
+
+ }
+ else if(TMath::Abs(eta1)>0.9 || TMath::Abs(eta2)>0.9){
+ // fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
+ // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+
+ if(isRealPi0 || isRealEta){
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_0912",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0912",openingAngleTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ }
+ if(!isRealPi0 && !isRealEta){
+ if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
+ fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
+ }else{
+ fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
+ }
+ }
+ }
+ else{
+ // fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
+ // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ if(isRealPi0 || isRealEta){
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_0909",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0909",openingAngleTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ if(gamma1MotherLabel > fV0Reader->GetMCStack()->GetNprimary()){
+ fHistograms->FillHistogram("ESD_TruePi0Sec_InvMass",massTwoGammaCandidate);
+ }
+ }
+ if(!isRealPi0 && !isRealEta){
+ if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
+ fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
+ }else{
+ fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
+ }
+
+ if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
+ fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ }
+ else if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 || TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
+ fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ }
+ else{
+ fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
+ }
+
+ Double_t lowMassPi0=0.1;
+ Double_t highMassPi0=0.15;
+ if (massTwoGammaCandidate > lowMassPi0 && massTwoGammaCandidate < highMassPi0 ){
+ new((*fKFReconstructedPi0sTClone)[fKFReconstructedPi0sTClone->GetEntriesFast()]) AliKFParticle(*twoGammaCandidate);
+ fGammav1.push_back(firstGammaIndex);
+ fGammav2.push_back(secondGammaIndex);
+ AddPionToAOD(twoGammaCandidate, massTwoGammaCandidate, firstGammaIndex, secondGammaIndex);
+ }
+ }
+
+ }
+ //}
+ delete twoGammaCandidate;
+ }
+ }
+}
+
+void AliAnalysisTaskGammaConversion::AddPionToAOD(const AliKFParticle * const pionkf, Double_t mass, Int_t daughter1, Int_t daughter2) {
+ //See header file for documentation
+ AliGammaConversionAODObject pion;
+ pion.SetPx(pionkf->GetPx());
+ pion.SetPy(pionkf->GetPy());
+ pion.SetPz(pionkf->GetPz());
+ pion.SetChi2(pionkf->GetChi2());
+ pion.SetE(pionkf->GetE());
+ pion.SetIMass(mass);
+ pion.SetLabel1(daughter1);
+ //dynamic_cast<AliGammaConversionAODObject*>(fAODBranch->At(daughter1))->SetTagged(kTRUE);
+ pion.SetLabel2(daughter2);
+ new((*fAODPi0)[fAODPi0->GetEntriesFast()]) AliGammaConversionAODObject(pion);
+
+}
+
+
+/*
+ void AliAnalysisTaskGammaConversion::ProcessConvPHOSGammasForNeutralMesonAnalysis(){
-void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
// see header file for documentation
+ // Analyse Pi0 with one photon from Phos and 1 photon from conversions
- // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
- // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
+
+
+ Double_t vtx[3];
+ vtx[0] = fV0Reader->GetPrimaryVertex()->GetX();
+ vtx[1] = fV0Reader->GetPrimaryVertex()->GetY();
+ vtx[2] = fV0Reader->GetPrimaryVertex()->GetZ();
+
+
+ // Loop over all CaloClusters and consider only the PHOS ones:
+ AliESDCaloCluster *clu;
+ TLorentzVector pPHOS;
+ TLorentzVector gammaPHOS;
+ TLorentzVector gammaGammaConv;
+ TLorentzVector pi0GammaConvPHOS;
+ TLorentzVector gammaGammaConvBck;
+ TLorentzVector pi0GammaConvPHOSBck;
+
+
+ for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
+ clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
+ if ( !clu->IsPHOS() || clu->E()<0.1 ) continue;
+ clu ->GetMomentum(pPHOS ,vtx);
for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
- for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
-
- // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
- // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
-
- AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
- AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
-
- if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
- continue;
+ AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
+ gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz(),0.);
+ gammaPHOS.SetXYZM(pPHOS.Px(),pPHOS.Py(),pPHOS.Pz(),0.);
+ pi0GammaConvPHOS=gammaGammaConv+gammaPHOS;
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS",pi0GammaConvPHOS.M());
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvPHOS",pi0GammaConvPHOS.M(),pi0GammaConvPHOS.Pt());
+
+ TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz());
+ TVector3 v3D1(gammaPHOS.Px(),gammaPHOS.Py(),gammaPHOS.Pz());
+ Double_t opanConvPHOS= v3D0.Angle(v3D1);
+ if ( opanConvPHOS < 0.35){
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanLow",pi0GammaConvPHOS.M());
+ }else{
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanHigh",pi0GammaConvPHOS.M());
+ }
+
+ }
+
+ // Now the LorentVector pPHOS is obtained and can be paired with the converted proton
+ }
+ //==== End of the PHOS cluster selection ============
+ TLorentzVector pEMCAL;
+ TLorentzVector gammaEMCAL;
+ TLorentzVector pi0GammaConvEMCAL;
+ TLorentzVector pi0GammaConvEMCALBck;
+
+ for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
+ clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
+ if ( !clu->IsEMCAL() || clu->E()<0.1 ) continue;
+ if (clu->GetNCells() <= 1) continue;
+ if ( clu->GetTOF()*1e9 < 550 || clu->GetTOF()*1e9 > 750) continue;
+
+ clu ->GetMomentum(pEMCAL ,vtx);
+ for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
+ AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
+ gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),
+ twoGammaDecayCandidateDaughter0->Py(),
+ twoGammaDecayCandidateDaughter0->Pz(),0.);
+ gammaEMCAL.SetXYZM(pEMCAL.Px(),pEMCAL.Py(),pEMCAL.Pz(),0.);
+ pi0GammaConvEMCAL=gammaGammaConv+gammaEMCAL;
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL",pi0GammaConvEMCAL.M());
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL",pi0GammaConvEMCAL.M(),pi0GammaConvEMCAL.Pt());
+ TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),
+ twoGammaDecayCandidateDaughter0->Py(),
+ twoGammaDecayCandidateDaughter0->Pz());
+ TVector3 v3D1(gammaEMCAL.Px(),gammaEMCAL.Py(),gammaEMCAL.Pz());
+
+
+ Double_t opanConvEMCAL= v3D0.Angle(v3D1);
+ if ( opanConvEMCAL < 0.35){
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanLow",pi0GammaConvEMCAL.M());
+ }else{
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanHigh",pi0GammaConvEMCAL.M());
+ }
+
+ }
+ if(fCalculateBackground){
+ for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
+ AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
+ for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
+ AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
+ gammaGammaConvBck.SetXYZM(previousGoodV0.Px(),
+ previousGoodV0.Py(),
+ previousGoodV0.Pz(),0.);
+ pi0GammaConvEMCALBck=gammaGammaConvBck+gammaEMCAL;
+ fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M());
+ fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M(),
+ pi0GammaConvEMCALBck.Pt());
+ }
+ }
+
+ // Now the LorentVector pEMCAL is obtained and can be paired with the converted proton
+ } // end of checking if background photons are available
+ }
+ //==== End of the PHOS cluster selection ============
+
+ }
+*/
+
+void AliAnalysisTaskGammaConversion::MoveParticleAccordingToVertex(AliKFParticle * particle,const AliGammaConversionBGHandler::GammaConversionVertex *vertex){
+ //see header file for documentation
+
+ Double_t dx = vertex->fX - fESDEvent->GetPrimaryVertex()->GetX();
+ Double_t dy = vertex->fY - fESDEvent->GetPrimaryVertex()->GetY();
+ Double_t dz = vertex->fZ - fESDEvent->GetPrimaryVertex()->GetZ();
+
+ // cout<<"dx, dy, dz: ["<<dx<<","<<dy<<","<<dz<<"]"<<endl;
+ particle->X() = particle->GetX() - dx;
+ particle->Y() = particle->GetY() - dy;
+ particle->Z() = particle->GetZ() - dz;
+}
+
+void AliAnalysisTaskGammaConversion::RotateKFParticle(AliKFParticle * kfParticle,Double_t angle){
+ // Rotate the kf particle
+ Double_t c = cos(angle);
+ Double_t s = sin(angle);
+
+ Double_t mA[7][ 7];
+ for( Int_t i=0; i<7; i++ ){
+ for( Int_t j=0; j<7; j++){
+ mA[i][j] = 0;
+ }
+ }
+ for( int i=0; i<7; i++ ){
+ mA[i][i] = 1;
+ }
+ mA[0][0] = c; mA[0][1] = s;
+ mA[1][0] = -s; mA[1][1] = c;
+ mA[3][3] = c; mA[3][4] = s;
+ mA[4][3] = -s; mA[4][4] = c;
+
+ Double_t mAC[7][7];
+ Double_t mAp[7];
+
+ for( Int_t i=0; i<7; i++ ){
+ mAp[i] = 0;
+ for( Int_t k=0; k<7; k++){
+ mAp[i]+=mA[i][k] * kfParticle->GetParameter(k);
+ }
+ }
+
+ for( Int_t i=0; i<7; i++){
+ kfParticle->Parameter(i) = mAp[i];
+ }
+
+ for( Int_t i=0; i<7; i++ ){
+ for( Int_t j=0; j<7; j++ ){
+ mAC[i][j] = 0;
+ for( Int_t k=0; k<7; k++ ){
+ mAC[i][j]+= mA[i][k] * kfParticle->GetCovariance(k,j);
}
- if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
- continue;
+ }
+ }
+
+ for( Int_t i=0; i<7; i++ ){
+ for( Int_t j=0; j<=i; j++ ){
+ Double_t xx = 0;
+ for( Int_t k=0; k<7; k++){
+ xx+= mAC[i][k]*mA[j][k];
}
-
- AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
-
- Double_t massTwoGammaCandidate = 0.;
- Double_t widthTwoGammaCandidate = 0.;
- Double_t chi2TwoGammaCandidate =10000.;
- twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
- if(twoGammaCandidate->GetNDF()>0){
- chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
-
- if(chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()){
-
- TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
- TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
-
- Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
- Double_t rapidity;
- if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() == 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() == 0){
- rapidity=0;
- }
- else{
- rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
+ kfParticle->Covariance(i,j) = xx;
+ }
+ }
+}
+
+
+void AliAnalysisTaskGammaConversion::CalculateBackground(){
+ // see header file for documentation
+
+
+ TClonesArray * currentEventV0s = fV0Reader->GetCurrentEventGoodV0s();
+
+ AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
+
+ Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
+ Int_t mbin = 0;
+ if(fUseTrackMultiplicityForBG == kTRUE){
+ mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
+ }
+ else{
+ mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
+ }
+
+ if(fDoRotation == kTRUE){
+ TRandom3 *random = new TRandom3();
+ for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
+ AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
+ for(Int_t iCurrent2=iCurrent+1;iCurrent2<currentEventV0s->GetEntriesFast();iCurrent2++){
+ for(Int_t nRandom=0;nRandom<fNRandomEventsForBG;nRandom++){
+
+ AliKFParticle currentEventGoodV02 = *(AliKFParticle *)(currentEventV0s->At(iCurrent2));
+
+ if(fCheckBGProbability == kTRUE){
+ Double_t massBGprob =0.;
+ Double_t widthBGprob = 0.;
+ AliKFParticle *backgroundCandidateProb = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
+ backgroundCandidateProb->GetMass(massBGprob,widthBGprob);
+ if(massBGprob>0.1 && massBGprob<0.14){
+ if(random->Rndm()>bgHandler->GetBGProb(zbin,mbin)){
+ delete backgroundCandidateProb;
+ continue;
+ }
+ }
+ delete backgroundCandidateProb;
}
+
+ Double_t nRadiansPM = fNDegreesPMBackground*TMath::Pi()/180;
+
+ Double_t rotationValue = random->Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
+
+ RotateKFParticle(¤tEventGoodV02,rotationValue);
+
+ AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
+
+ Double_t massBG =0.;
+ Double_t widthBG = 0.;
+ Double_t chi2BG =10000.;
+ backgroundCandidate->GetMass(massBG,widthBG);
+
+ // if(backgroundCandidate->GetNDF()>0){
+ chi2BG = backgroundCandidate->GetChi2();
+ if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
+
+ TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
+ TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
+
+ Double_t openingAngleBG = currentEventGoodV0.GetAngle(currentEventGoodV02);
+
+ Double_t rapidity;
+ if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
+ else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
+
+ if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
+ delete backgroundCandidate;
+ continue; // rapidity cut
+ }
- if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
- delete twoGammaCandidate;
- continue; // minimum opening angle to avoid using ghosttracks
- }
-
- fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
- fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
- fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
- fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
- fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
- fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
- fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
- fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
- fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
- fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
- fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
- fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
-
- if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
- fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
- fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
+
+ Double_t alfa=0.0;
+ if( (currentEventGoodV0.GetE()+currentEventGoodV02.GetE()) != 0){
+ alfa=TMath::Abs((currentEventGoodV0.GetE()-currentEventGoodV02.GetE())
+ /(currentEventGoodV0.GetE()+currentEventGoodV02.GetE()));
+ }
+
+
+ if(openingAngleBG < fMinOpeningAngleGhostCut ){
+ delete backgroundCandidate;
+ continue; // minimum opening angle to avoid using ghosttracks
+ }
+
+ // original
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
+ fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
+ fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
+ fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
+ fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
+ fHistograms->FillHistogram("ESD_Background_Mass", massBG);
+ fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
+ fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
+
+
+ if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
+ }
+
+ fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
+ fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
+
+ if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
+ }
+ }
+ if(alfa<0.1){
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
+ }
+
}
+ //}
+ delete backgroundCandidate;
}
}
- delete twoGammaCandidate;
}
+ delete random;
}
-}
+ else{ // means no rotation
+ AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
+
+ if(fUseTrackMultiplicityForBG){
+ // cout<<"Using charged track multiplicity for background calculation"<<endl;
+ for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
-void AliAnalysisTaskGammaConversion::CalculateBackground(){
- // see header file for documentation
-
- vector<AliKFParticle> vectorCurrentEventGoodV0s = fV0Reader->GetCurrentEventGoodV0s();
- vector<AliKFParticle> vectorPreviousEventGoodV0s = fV0Reader->GetPreviousEventGoodV0s();
-
- for(UInt_t iCurrent=0;iCurrent<vectorCurrentEventGoodV0s.size();iCurrent++){
- AliKFParticle * currentEventGoodV0 = &vectorCurrentEventGoodV0s.at(iCurrent);
- for(UInt_t iPrevious=0;iPrevious<vectorPreviousEventGoodV0s.size();iPrevious++){
- AliKFParticle * previousGoodV0 = &vectorPreviousEventGoodV0s.at(iPrevious);
-
- AliKFParticle *backgroundCandidate = new AliKFParticle(*currentEventGoodV0,*previousGoodV0);
-
- Double_t massBG =0.;
- Double_t widthBG = 0.;
- Double_t chi2BG =10000.;
- backgroundCandidate->GetMass(massBG,widthBG);
- if(backgroundCandidate->GetNDF()>0){
- chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
- if(chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()){
+ AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);//fV0Reader->GetBGGoodV0s(nEventsInBG);
+
+ if(fMoveParticleAccordingToVertex == kTRUE){
+ bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
+ }
+
+ for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
+ AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
+ for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
+ AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
+ AliKFParticle previousGoodV0test = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
+
+ //cout<<"Primary Vertex event: ["<<fESDEvent->GetPrimaryVertex()->GetX()<<","<<fESDEvent->GetPrimaryVertex()->GetY()<<","<<fESDEvent->GetPrimaryVertex()->GetZ()<<"]"<<endl;
+ //cout<<"BG prim Vertex event: ["<<bgEventVertex->fX<<","<<bgEventVertex->fY<<","<<bgEventVertex->fZ<<"]"<<endl;
+
+ //cout<<"XYZ of particle before transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
+ if(fMoveParticleAccordingToVertex == kTRUE){
+ MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
+ }
+ //cout<<"XYZ of particle after transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
+
+ AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
+
+ Double_t massBG =0.;
+ Double_t widthBG = 0.;
+ Double_t chi2BG =10000.;
+ backgroundCandidate->GetMass(massBG,widthBG);
+ // if(backgroundCandidate->GetNDF()>0){
+ // chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
+ chi2BG = backgroundCandidate->GetChi2();
+ if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
- TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
- TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
+ TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
+ TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
- Double_t openingAngleBG = currentEventGoodV0->GetAngle(*previousGoodV0);
+ Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
- Double_t rapidity;
- if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
- else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
+ Double_t rapidity;
+
+ if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() <= 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() <= 0){
+ cout << "Error: |Pz| > E !!!! " << endl;
+ rapidity=0;
+ } else {
+ rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
+ }
+ if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
+ delete backgroundCandidate;
+ continue; // rapidity cut
+ }
+
+
+ Double_t alfa=0.0;
+ if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
+ alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
+ /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
+ }
+
+ if(openingAngleBG < fMinOpeningAngleGhostCut ){
+ delete backgroundCandidate;
+ continue; // minimum opening angle to avoid using ghosttracks
+ }
+
+ // original
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
+ fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
+ fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
+ fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
+ fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
+ fHistograms->FillHistogram("ESD_Background_Mass", massBG);
+ fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
+ fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
+
+
+ if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
+ }
+
+ // test
+ fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
+ fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
+
+ if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
+ }
+ // }
+ }
+ if(alfa<0.1){
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
+ }
+
+ }
+ delete backgroundCandidate;
+ }
+ }
+ }
+ }
+ else{ // means using #V0s for multiplicity
+
+ // cout<<"Using the v0 multiplicity to calculate background"<<endl;
+
+ fHistograms->FillHistogram("ESD_Background_z_m",zbin,mbin);
+ fHistograms->FillHistogram("ESD_Mother_multpilicityVSv0s",fV0Reader->CountESDTracks(),fV0Reader->GetNumberOfV0s());
+
+ for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
+ AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);// fV0Reader->GetBGGoodV0s(nEventsInBG);
+ if(previousEventV0s){
+
+ if(fMoveParticleAccordingToVertex == kTRUE){
+ bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
+ }
+
+ for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
+ AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
+ for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
+ AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
+
+ if(fMoveParticleAccordingToVertex == kTRUE){
+ MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
+ }
+
+ AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
+ Double_t massBG =0.;
+ Double_t widthBG = 0.;
+ Double_t chi2BG =10000.;
+ backgroundCandidate->GetMass(massBG,widthBG);
+ /* if(backgroundCandidate->GetNDF()>0){
+ chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
+ {//remember to remove
+ TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
+ TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
+
+ Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
+ fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle_nochi2", openingAngleBG);
+ }
+ */
+ chi2BG = backgroundCandidate->GetChi2();
+ if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
+ TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
+ TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
+ Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
+ Double_t rapidity;
+ if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
+ else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
- if(openingAngleBG < fMinOpeningAngleGhostCut ){
- delete backgroundCandidate;
- continue; // minimum opening angle to avoid using ghosttracks
- }
+ if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
+ delete backgroundCandidate;
+ continue; // rapidity cut
+ }
+
+
+ Double_t alfa=0.0;
+ if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
+ alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
+ /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
+ }
+
- fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
- fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
- fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
- fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
- fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
- fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
- fHistograms->FillHistogram("ESD_Background_Mass", massBG);
- fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
- fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
- fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
- fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
- fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
-
- if ( TMath::Abs(currentEventGoodV0->GetEta())<0.9 && TMath::Abs(previousGoodV0->GetEta())<0.9 ){
- fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
- fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
- }
+ if(openingAngleBG < fMinOpeningAngleGhostCut ){
+ delete backgroundCandidate;
+ continue; // minimum opening angle to avoid using ghosttracks
+ }
+ if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
+ fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
+ fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
+ fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
+ fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
+ fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
+ fHistograms->FillHistogram("ESD_Background_Mass", massBG);
+ fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
+ fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
+
+
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
+
+
+ if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
+ }
+
+ if(massBG>0.5 && massBG<0.6){
+ fHistograms->FillHistogram("ESD_Background_alfa_pt0506",momentumVectorbackgroundCandidate.Pt(),alfa);
+ }
+ if(massBG>0.3 && massBG<0.4){
+ fHistograms->FillHistogram("ESD_Background_alfa_pt0304",momentumVectorbackgroundCandidate.Pt(),alfa);
+ }
+
+ // test
+ fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
+ fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
+ fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
+
+ if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
+ fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
+ }
+ }
+
+ if(alfa<0.1){
+ fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
+ }
+ // }
+ }
+ delete backgroundCandidate;
+ }
+ }
}
}
- delete backgroundCandidate;
- }
- }
+ } // end else (means use #v0s as multiplicity)
+ } // end no rotation
}
-
void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){
//ProcessGammasForGammaJetAnalysis
}
}
+//____________________________________________________________________
+Bool_t AliAnalysisTaskGammaConversion::IsGoodImpPar(const AliESDtrack *const track)
+{
+//
+// check whether particle has good DCAr(Pt) impact
+// parameter. Only for TPC+ITS tracks (7*sigma cut)
+// Origin: Andrea Dainese
+//
+
+Float_t d0z0[2],covd0z0[3];
+track->GetImpactParameters(d0z0,covd0z0);
+Float_t sigma= 0.0050+0.0060/TMath::Power(track->Pt(),0.9);
+Float_t d0max = 7.*sigma;
+if(TMath::Abs(d0z0[0]) < d0max) return kTRUE;
+
+return kFALSE;
+}
+
+
void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){
// CreateListOfChargedParticles
fESDEvent = fV0Reader->GetESDEvent();
+ Int_t numberOfESDTracks=0;
for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
if(!curTrack){
continue;
}
+ // Not needed if Standard function used.
+// if(!IsGoodImpPar(curTrack)){
+// continue;
+// }
if(fEsdTrackCuts->AcceptTrack(curTrack) ){
new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack);
// fChargedParticles.push_back(curTrack);
fChargedParticlesId.push_back(iTracks);
+ numberOfESDTracks++;
}
}
+// Moved to UserExec using CountAcceptedTracks function. runjet is not needed by default
+// fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",numberOfESDTracks);
+// cout<<"esdtracks::"<< numberOfESDTracks<<endl;
+// if (fV0Reader->GetNumberOfContributorsVtx()>=1){
+// fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",numberOfESDTracks);
+// }
+}
+void AliAnalysisTaskGammaConversion::RecalculateV0ForGamma(){
+ //recalculates v0 for gamma
+
+ Double_t massE=0.00051099892;
+ TLorentzVector curElecPos;
+ TLorentzVector curElecNeg;
+ TLorentzVector curGamma;
+
+ TLorentzVector curGammaAt;
+ TLorentzVector curElecPosAt;
+ TLorentzVector curElecNegAt;
+ AliKFVertex primVtxGamma(*(fESDEvent->GetPrimaryVertex()));
+ AliKFVertex primVtxImprovedGamma = primVtxGamma;
+
+ const AliESDVertex *vtxT3D=fESDEvent->GetPrimaryVertex();
+
+ Double_t xPrimaryVertex=vtxT3D->GetXv();
+ Double_t yPrimaryVertex=vtxT3D->GetYv();
+ Double_t zPrimaryVertex=vtxT3D->GetZv();
+ // Float_t primvertex[3]={xPrimaryVertex,yPrimaryVertex,zPrimaryVertex};
+
+ Float_t nsigmaTPCtrackPos;
+ Float_t nsigmaTPCtrackNeg;
+ Float_t nsigmaTPCtrackPosToPion;
+ Float_t nsigmaTPCtrackNegToPion;
+ AliKFParticle* negKF=NULL;
+ AliKFParticle* posKF=NULL;
+
+ for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
+ AliESDtrack* posTrack = fESDEvent->GetTrack(iTracks);
+ if(!posTrack){
+ continue;
+ }
+ if (posKF) delete posKF; posKF=NULL;
+ if(posTrack->GetSign()<0) continue;
+ if(!(posTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
+ if(posTrack->GetKinkIndex(0)>0 ) continue;
+ if(posTrack->GetNcls(1)<50)continue;
+ Double_t momPos[3];
+ // posTrack->GetConstrainedPxPyPz(momPos);
+ posTrack->GetPxPyPz(momPos);
+ AliESDtrack *ptrk=fESDEvent->GetTrack(iTracks);
+ curElecPos.SetXYZM(momPos[0],momPos[1],momPos[2],massE);
+ if(TMath::Abs(curElecPos.Eta())<0.9) continue;
+ posKF = new AliKFParticle( *(posTrack),-11);
+
+ nsigmaTPCtrackPos = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kElectron);
+ nsigmaTPCtrackPosToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion);
+
+ if ( nsigmaTPCtrackPos>5.|| nsigmaTPCtrackPos<-2.){
+ continue;
+ }
+
+ if(pow((momPos[0]*momPos[0]+momPos[1]*momPos[1]+momPos[2]*momPos[2]),0.5)>0.5 && nsigmaTPCtrackPosToPion<1){
+ continue;
+ }
+
+
+
+ for(Int_t jTracks = 0; jTracks < fESDEvent->GetNumberOfTracks(); jTracks++){
+ AliESDtrack* negTrack = fESDEvent->GetTrack(jTracks);
+ if(!negTrack){
+ continue;
+ }
+ if (negKF) delete negKF; negKF=NULL;
+ if(negTrack->GetSign()>0) continue;
+ if(!(negTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
+ if(negTrack->GetKinkIndex(0)>0 ) continue;
+ if(negTrack->GetNcls(1)<50)continue;
+ Double_t momNeg[3];
+ // negTrack->GetConstrainedPxPyPz(momNeg);
+ negTrack->GetPxPyPz(momNeg);
+
+ nsigmaTPCtrackNeg = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kElectron);
+ nsigmaTPCtrackNegToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion);
+ if ( nsigmaTPCtrackNeg>5. || nsigmaTPCtrackNeg<-2.){
+ continue;
+ }
+ if(pow((momNeg[0]*momNeg[0]+momNeg[1]*momNeg[1]+momNeg[2]*momNeg[2]),0.5)>0.5 && nsigmaTPCtrackNegToPion<1){
+ continue;
+ }
+ AliESDtrack *ntrk=fESDEvent->GetTrack(jTracks);
+ curElecNeg.SetXYZM(momNeg[0],momNeg[1],momNeg[2],massE);
+ if(TMath::Abs(curElecNeg.Eta())<0.9) continue;
+ negKF = new AliKFParticle( *(negTrack) ,11);
+
+ Double_t b=fESDEvent->GetMagneticField();
+ Double_t xn, xp, dca=ntrk->GetDCA(ptrk,b,xn,xp);
+ AliExternalTrackParam nt(*ntrk), pt(*ptrk);
+ nt.PropagateTo(xn,b); pt.PropagateTo(xp,b);
+
+
+ //--- Like in ITSV0Finder
+ AliExternalTrackParam ntAt0(*ntrk), ptAt0(*ptrk);
+ Double_t xxP,yyP,alphaP;
+ Double_t rP[3];
+
+ // if (!ptAt0.GetGlobalXYZat(ptAt0->GetX(),xxP,yyP,zzP)) continue;
+ if (!ptAt0.GetXYZAt(ptAt0.GetX(),b,rP)) continue;
+ xxP=rP[0];
+ yyP=rP[1];
+ alphaP = TMath::ATan2(yyP,xxP);
+
+
+ ptAt0.Propagate(alphaP,0,b);
+ Float_t ptfacP = (1.+100.*TMath::Abs(ptAt0.GetC(b)));
+
+ // Double_t distP = ptAt0.GetY();
+ Double_t normP = ptfacP*TMath::Sqrt(ptAt0.GetSigmaY2());
+ Double_t normdist0P = TMath::Abs(ptAt0.GetY()/normP);
+ Double_t normdist1P = TMath::Abs((ptAt0.GetZ()-zPrimaryVertex)/(ptfacP*TMath::Sqrt(ptAt0.GetSigmaZ2())));
+ Double_t normdistP = TMath::Sqrt(normdist0P*normdist0P+normdist1P*normdist1P);
+
+
+ Double_t xxN,yyN,alphaN;
+ Double_t rN[3];
+ // if (!ntAt0.GetGlobalXYZat(ntAt0->GetX(),xxN,yyN,zzN)) continue;
+ if (!ntAt0.GetXYZAt(ntAt0.GetX(),b,rN)) continue;
+ xxN=rN[0];
+ yyN=rN[1];
+
+ alphaN = TMath::ATan2(yyN,xxN);
+
+ ntAt0.Propagate(alphaN,0,b);
+
+ Float_t ptfacN = (1.+100.*TMath::Abs(ntAt0.GetC(b)));
+ // Double_t distN = ntAt0.GetY();
+ Double_t normN = ptfacN*TMath::Sqrt(ntAt0.GetSigmaY2());
+ Double_t normdist0N = TMath::Abs(ntAt0.GetY()/normN);
+ Double_t normdist1N = TMath::Abs((ntAt0.GetZ()-zPrimaryVertex)/(ptfacN*TMath::Sqrt(ntAt0.GetSigmaZ2())));
+ Double_t normdistN = TMath::Sqrt(normdist0N*normdist0N+normdist1N*normdist1N);
+
+ //-----------------------------
+
+ Double_t momNegAt[3];
+ nt.GetPxPyPz(momNegAt);
+ curElecNegAt.SetXYZM(momNegAt[0],momNegAt[1],momNegAt[2],massE);
+
+ Double_t momPosAt[3];
+ pt.GetPxPyPz(momPosAt);
+ curElecPosAt.SetXYZM(momPosAt[0],momPosAt[1],momPosAt[2],massE);
+ if(dca>1){
+ continue;
+ }
+
+ // Double_t dneg= negTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
+ // Double_t dpos= posTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
+ AliESDv0 vertex(nt,jTracks,pt,iTracks);
+
+
+ Float_t cpa=vertex.GetV0CosineOfPointingAngle(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex);
+
+
+
+ // cout<< "v0Rr::"<< v0Rr<<endl;
+ // if (pvertex.GetRr()<0.5){
+ // continue;
+ //}
+ // cout<<"vertex.GetChi2V0()"<<vertex.GetChi2V0()<<endl;
+ if(cpa<0.9)continue;
+ // if (vertex.GetChi2V0() > 30) continue;
+ // cout<<"xp+xn::"<<xp<<" "<<xn<<endl;
+ if ((xn+xp) < 0.4) continue;
+ if (TMath::Abs(ntrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05)
+ if (TMath::Abs(ptrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05) continue;
+
+ //cout<<"pass"<<endl;
+
+ AliKFParticle v0GammaC;
+ v0GammaC+=(*negKF);
+ v0GammaC+=(*posKF);
+ v0GammaC.SetMassConstraint(0,0.001);
+ primVtxImprovedGamma+=v0GammaC;
+ v0GammaC.SetProductionVertex(primVtxImprovedGamma);
+
+
+ curGamma=curElecNeg+curElecPos;
+ curGammaAt=curElecNegAt+curElecPosAt;
+
+ // invariant mass versus pt of K0short
+
+ Double_t chi2V0GammaC=100000.;
+ if( v0GammaC.GetNDF() != 0) {
+ chi2V0GammaC = v0GammaC.GetChi2()/v0GammaC.GetNDF();
+ }else{
+ cout<< "ERROR::v0K0C.GetNDF()" << endl;
+ }
+
+ if(chi2V0GammaC<200 &&chi2V0GammaC>0 ){
+ if(fHistograms != NULL){
+ fHistograms->FillHistogram("ESD_RecalculateV0_InvMass",v0GammaC.GetMass());
+ fHistograms->FillHistogram("ESD_RecalculateV0_Pt",v0GammaC.GetPt());
+ fHistograms->FillHistogram("ESD_RecalculateV0_E_dEdxP",curElecNegAt.P(),negTrack->GetTPCsignal());
+ fHistograms->FillHistogram("ESD_RecalculateV0_P_dEdxP",curElecPosAt.P(),posTrack->GetTPCsignal());
+ fHistograms->FillHistogram("ESD_RecalculateV0_cpa",cpa);
+ fHistograms->FillHistogram("ESD_RecalculateV0_dca",dca);
+ fHistograms->FillHistogram("ESD_RecalculateV0_normdistP",normdistP);
+ fHistograms->FillHistogram("ESD_RecalculateV0_normdistN",normdistN);
+
+ new((*fKFRecalculatedGammasTClone)[fKFRecalculatedGammasTClone->GetEntriesFast()]) AliKFParticle(v0GammaC);
+ fElectronRecalculatedv1.push_back(iTracks);
+ fElectronRecalculatedv2.push_back(jTracks);
+ }
+ }
+ }
+ }
+
+ for(Int_t firstGammaIndex=0;firstGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();firstGammaIndex++){
+ for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();secondGammaIndex++){
+ AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(firstGammaIndex);
+ AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(secondGammaIndex);
+
+ if(fElectronRecalculatedv1[firstGammaIndex]==fElectronRecalculatedv1[secondGammaIndex]){
+ continue;
+ }
+ if( fElectronRecalculatedv2[firstGammaIndex]==fElectronRecalculatedv2[secondGammaIndex]){
+ continue;
+ }
+
+ AliKFParticle twoGammaCandidate(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
+ fHistograms->FillHistogram("ESD_RecalculateGG_InvMass",twoGammaCandidate.GetMass());
+ fHistograms->FillHistogram("ESD_RecalculateGG_InvMass_vs_Pt",twoGammaCandidate.GetMass(),twoGammaCandidate.GetPt());
+ }
+ }
}
void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){
// CaculateJetCone
void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()
{
//AOD
- if(fAODBranch==NULL){
- fAODBranch = new TClonesArray("AliGammaConversionAODObject", 0);
+ if(!fAODGamma) fAODGamma = new TClonesArray("AliGammaConversionAODObject", 0);
+ else fAODGamma->Delete();
+ fAODGamma->SetName(Form("%s_gamma", fAODBranchName.Data()));
+
+ if(!fAODPi0) fAODPi0 = new TClonesArray("AliGammaConversionAODObject", 0);
+ else fAODPi0->Delete();
+ fAODPi0->SetName(Form("%s_Pi0", fAODBranchName.Data()));
+
+ if(!fAODOmega) fAODOmega = new TClonesArray("AliGammaConversionAODObject", 0);
+ else fAODOmega->Delete();
+ fAODOmega->SetName(Form("%s_Omega", fAODBranchName.Data()));
+
+ //If delta AOD file name set, add in separate file. Else add in standard aod file.
+ if(fKFDeltaAODFileName.Length() > 0) {
+ AddAODBranch("TClonesArray", &fAODGamma, fKFDeltaAODFileName.Data());
+ AddAODBranch("TClonesArray", &fAODPi0, fKFDeltaAODFileName.Data());
+ AddAODBranch("TClonesArray", &fAODOmega, fKFDeltaAODFileName.Data());
+ AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(fKFDeltaAODFileName.Data());
+ } else {
+ AddAODBranch("TClonesArray", &fAODGamma);
+ AddAODBranch("TClonesArray", &fAODPi0);
+ AddAODBranch("TClonesArray", &fAODOmega);
}
- fAODBranch->Delete();
- fAODBranch->SetName(fAODBranchName);
- AddAODBranch("TClonesArray", &fAODBranch);
-
+
// Create the output container
if(fOutputContainer != NULL){
delete fOutputContainer;
}
if(fOutputContainer == NULL){
fOutputContainer = new TList();
+ fOutputContainer->SetOwner(kTRUE);
}
//Adding the histograms to the output container
fNeutralMesonNtuple = new TNtuple("NeutralMesonNtuple","NeutralMesonNtuple","test");
}
TList * ntupleTList = new TList();
+ ntupleTList->SetOwner(kTRUE);
ntupleTList->SetName("Ntuple");
ntupleTList->Add((TNtuple*)fGammaNtuple);
fOutputContainer->Add(ntupleTList);
fOutputContainer->SetName(GetName());
}
-Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(TParticle* const daughter0, TParticle* const daughter1) const{
+Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(const TParticle* const daughter0, const TParticle* const daughter1) const{
//helper function
TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());
TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());
void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){
// see header file for documentation
-
+
vector<Int_t> indexOfGammaParticle;
fStack = fV0Reader->GetMCStack();
- Int_t labelMC = TMath::Abs(curTrack->GetLabel());
- TParticle* curParticle = fStack->Particle(labelMC);
-
-
TLorentzVector curElec;
curElec.SetXYZM(p[0],p[1],p[2],fElectronMass);
+ if(fDoMCTruth){
+ Int_t labelMC = TMath::Abs(curTrack->GetLabel());
+ TParticle* curParticle = fStack->Particle(labelMC);
+ if(curTrack->GetSign() > 0){
+ if( pid == 0){
+ fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
+ fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
+ }
+ else{
+ fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
+ fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
+ }
+ }
+ }
if(curTrack->GetSign() > 0){
fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
fHistograms->FillHistogram("ESD_ElectronPosPt",curElec.Pt());
- fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
+ // fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
- fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
+ // fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
// vESDePosTemp.push_back(curTrack);
new((*vESDePosTemp)[vESDePosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
}
else {
- // vESDxNegTemp.push_back(curTrack);
- /* if(vESDxNegTemp == NULL){
- cout<<"TCloes is zero"<<endl;
- }
- if(curTrack == NULL){
- cout<<"curTrack is zero"<<endl;
- }
- */
+
new((*vESDxNegTemp)[vESDxNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
if( pid == 0){
fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
fHistograms->FillHistogram("ESD_ElectronNegPt",curElec.Pt());
- fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
- fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
- //vESDeNegTemp.push_back(curTrack);
new((*vESDeNegTemp)[vESDeNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
}
}
-void AliAnalysisTaskGammaConversion::GetPID(AliESDtrack *track, Stat_t &pid, Stat_t &weight)
+void AliAnalysisTaskGammaConversion::GetPID(const AliESDtrack *track, Stat_t &pid, Stat_t &weight)
{
// see header file for documentation
pid = -1;
pid = ipid;
weight = max;
}
-double AliAnalysisTaskGammaConversion::GetSigmaToVertex(AliESDtrack* t)
+double AliAnalysisTaskGammaConversion::GetSigmaToVertex(const AliESDtrack* t)
{
// Calculates the number of sigma to the vertex.
// return tlVtrack;
}
+Int_t AliAnalysisTaskGammaConversion::GetProcessType(const AliMCEvent * mcEvt) {
+
+ // Determine if the event was generated with pythia or phojet and return the process type
+
+ // Check if mcEvt is fine
+ if (!mcEvt) {
+ AliFatal("NULL mc event");
+ }
+
+ // Determine if it was a pythia or phojet header, and return the correct process type
+ AliGenPythiaEventHeader * headPy = 0;
+ AliGenDPMjetEventHeader * headPho = 0;
+ AliGenEventHeader * htmp = mcEvt->GenEventHeader();
+ if(!htmp) {
+ AliFatal("Cannot Get MC Header!!");
+ }
+ if( TString(htmp->IsA()->GetName()) == "AliGenPythiaEventHeader") {
+ headPy = (AliGenPythiaEventHeader*) htmp;
+ } else if (TString(htmp->IsA()->GetName()) == "AliGenDPMjetEventHeader") {
+ headPho = (AliGenDPMjetEventHeader*) htmp;
+ } else {
+ AliError("Unknown header");
+ }
+
+ // Determine process type
+ if(headPy) {
+ if(headPy->ProcessType() == 92 || headPy->ProcessType() == 93) {
+ // single difractive
+ return kProcSD;
+ } else if (headPy->ProcessType() == 94) {
+ // double diffractive
+ return kProcDD;
+ }
+ else if(headPy->ProcessType() != 92 && headPy->ProcessType() != 93 && headPy->ProcessType() != 94) {
+ // non difractive
+ return kProcND;
+ }
+ } else if (headPho) {
+ if(headPho->ProcessType() == 5 || headPho->ProcessType() == 6 ) {
+ // single difractive
+ return kProcSD;
+ } else if (headPho->ProcessType() == 7) {
+ // double diffractive
+ return kProcDD;
+ } else if(headPho->ProcessType() != 5 && headPho->ProcessType() != 6 && headPho->ProcessType() != 7 ) {
+ // non difractive
+ return kProcND;
+ }
+ }
+
+
+ // no process type found?
+ AliError(Form("Unknown header: %s", htmp->IsA()->GetName()));
+ return kProcUnknown;
+}
+
+
+Int_t AliAnalysisTaskGammaConversion::CalculateMultiplicityBin(){
+ // Get Centrality bin
+
+ Int_t multiplicity = 0;
+ if ( fUseMultiplicity == 1 ) {
+
+ if (fMultiplicity>= 0 && fMultiplicity<= 5) multiplicity=1;
+ if (fMultiplicity>= 6 && fMultiplicity<= 9) multiplicity=2;
+ if (fMultiplicity>=10 && fMultiplicity<=14) multiplicity=3;
+ if (fMultiplicity>=15 && fMultiplicity<=22) multiplicity=4;
+ if (fMultiplicity>=23 ) multiplicity=5;
+
+ }
+ return multiplicity;
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff