DIFFRACTIVE/example/AliAnalysisTaskCDex.cxx
DIFFRACTIVE/example/AliCDMesonUtilsStripped.cxx
DIFFRACTIVE/example/AliCDMesonBaseStripped.cxx
+ DIFFRACTIVE/xsAndTwoProng/AliAnalysisTaskCDMeson.cxx
+ DIFFRACTIVE/xsAndTwoProng/AliAnalysisTaskCDskimESD.cxx
+ DIFFRACTIVE/xsAndTwoProng/AliCDMesonUtils.cxx
+ DIFFRACTIVE/xsAndTwoProng/AliCDMesonBase.cxx
DIFFRACTIVE/xsAndTwoProng/AliCDMesonTracks.cxx)
string ( REPLACE ".cxx" ".h" HDRS "${SRCS}" )
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+//
+// Select events according to gap conditions, analyze two track events in pp
+// collisions
+//
+// Author:
+// Xianguo Lu <lu@physi.uni-heidelberg.de>
+// continued by
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#include <TH1.h>
+#include <TH2.h>
+#include <TList.h>
+#include <THnSparse.h>
+
+#include <TRandom3.h>
+
+#include "AliAODInputHandler.h"
+#include "AliAODHandler.h"
+#include "AliESDInputHandler.h"
+#include "AliPIDResponse.h"
+#include "AliPhysicsSelection.h"
+#include "AliGenEventHeader.h"
+#include "AliGenPythiaEventHeader.h"
+#include "AliGenDPMjetEventHeader.h"
+#include "AliMCEvent.h"
+#include "AliESDtrack.h"
+#include "AliStack.h"
+
+#include "AliCDMesonBase.h"
+#include "AliCDMesonUtils.h"
+#include "AliCDMesonTracks.h"
+#include "AliAnalysisTaskCDMeson.h"
+
+
+//------------------------------------------------------------------------------
+AliAnalysisTaskCDMeson::AliAnalysisTaskCDMeson(const char* name, Long_t state):
+ AliAnalysisTaskSE(name)
+ , fDoAOD(kFALSE)
+ , fAnalysisStatus(state)
+ , fAnalysisStatusString(new TObjString(Form("0x%lx", fAnalysisStatus)))
+ , fNoGapFraction(0.01)
+ , fReducedGapFraction(0.02)
+ , fMaxVtxDst(0.5) // value to be checked with the vertex study histograms
+
+ , fESDEvent(0x0)
+ , fAODEvent(0x0)
+ , fPIDResponse(0x0)
+ , fPhysicsSelection(0x0)
+ , fTracks(new AliCDMesonTracks())
+ , fVtxDst(-1.)
+ , fVtxZ(-20)
+ , fResidualTracks(0)
+ , fResidualTracklets(0)
+ , fMCprocessType(0)
+ , fMCprocess(-1)
+
+ , fRun(-999)
+ , fPIDmode(0)
+
+ , fTheta(0.)
+ , fPhi(0.)
+ , fMass(0.)
+ , fCurrentGapCondition(0)
+
+ , fGapRun(0x0)
+ , fHist(0x0)
+ , fThnMother(0x0)
+ , fThnMotherSoft(0x0)
+ , fThnMultiplicity(0x0)
+ , fThnMC(0x0)
+ , fThnEmptyEvents(0x0)
+ , fPWAtree(0x0)
+
+ , fv0ntrk(0x0)
+ , fv0fmdntrk(0x0)
+ , fv0fmdspdntrk(0x0)
+ , fv0fmdspdtpcntrk(0x0)
+
+ , fv0Rmntrk(0x0)
+ , fv0fmdRmntrk(0x0)
+ , fv0fmdspdRmntrk(0x0)
+ , fv0fmdspdtpcRmntrk(0x0)
+
+ , fMultStudy(0x0)
+ , fMultStudyV0dg(0x0)
+ , fMultStudyV0FMDdg(0x0)
+ , fMultStudyV0FMDSPDdg(0x0)
+ , fMultStudyV0FMDSPDTPCdg(0x0)
+
+ , fMultResponseMC(0x0)
+ , fMultResponseV0dgMC(0x0)
+ , fMultResponseV0FMDdgMC(0x0)
+ , fMultResponseV0FMDSPDdgMC(0x0)
+ , fMultResponseV0FMDSPDTPCdgMC(0x0)
+
+ , fMultRegionsMC(0x0)
+
+ , fhspdV0dg(0x0)
+ , fhspdV0FMDdg(0x0)
+ , fhspdV0FMDSPDdg(0x0)
+ , fhspdV0FMDSPDTPCdg(0x0)
+ , fhspdAfterCuts(0x0)
+
+ , fGapResponseMCv0Dg(0x0)
+ , fGapResponseMCv0fmdDg(0x0)
+ , fGapResponseMCv0fmdspdDg(0x0)
+ , fGapResponseMCv0fmdspdtpcDg(0x0)
+
+ , fhspd(0x0)
+ , fhfo(0x0)
+ , fhpriVtxDist(0x0)
+ , fhpriVtxPos(0x0)
+ , fhpriv(0x0)
+ , fhntrk(0x0)
+ , fhStatsFlow(0x0)
+ , fhFOchans(0x0)
+
+ , fVZEROhists(0x0)
+
+ , fv0Map(0x0)
+ , fv0fmdMap(0x0)
+ , fv0fmdspdMap(0x0)
+ , fv0fmdspdtpcMap(0x0)
+
+ , fv0MapCutted(0x0)
+ , fv0fmdMapCutted(0x0)
+ , fv0fmdspdMapCutted(0x0)
+ , fv0fmdspdtpcMapCutted(0x0)
+
+ , fHitMapSPDinner(0x0)
+ , fHitMapSPDouter(0x0)
+ , fHitMapSPDtrklt(0x0)
+ , fHitMapFMDa(0x0)
+ , fHitMapFMDc(0x0)
+
+ , fFMDsum1I(0x0)
+ , fFMDsum2I(0x0)
+ , fFMDsum2O(0x0)
+ , fFMDsum3I(0x0)
+ , fFMDsum3O(0x0)
+
+ , fPriVtxX(0x0)
+ , fPriVtxY(0x0)
+ , fPriVtxZ(0x0)
+ , fPriVtxDst(0x0)
+ , fPriVtxDstV0dg(0x0)
+ , fPriVtxDstV0FMDdg(0x0)
+ , fPriVtxDstV0FMDSPDdg(0x0)
+ , fPriVtxDstV0FMDSPDTPCdg(0x0)
+
+ , fTPCGapDCAaSide(0x0)
+ , fTPCGapDCAcSide(0x0)
+
+ , fComb2trkPIDuls(0x0)
+ , fComb2trkPIDls(0x0)
+ , fComb2trkPIDulsDG(0x0)
+ , fComb2trkPIDlsDG(0x0)
+ , fComb2trkPIDulsMC(0x0)
+ , fComb2trkPIDlsMC(0x0)
+ , fComb2trkPIDulsDGmc(0x0)
+ , fComb2trkPIDlsDGmc(0x0)
+ , fMCProcessUls(0x0)
+ , fMCProcessLs(0x0)
+
+ , fMAllTrackMass(0x0)
+{
+ //
+ // standard constructor (the one which should be used)
+ //
+ // slot in TaskSE must start from 1
+ Int_t iOutputSlot = 1;
+ DefineOutput(iOutputSlot++, TList::Class());
+ if (!(fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) {
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitTHnMother)) {
+ DefineOutput(iOutputSlot++, THnSparse::Class());
+ }
+ if (!fAnalysisStatus
+ || ((fAnalysisStatus & AliCDMesonBase::kBitSoftTracks)
+ && (fAnalysisStatus & AliCDMesonBase::kBitTHnMother))) {
+ DefineOutput(iOutputSlot++, THnSparse::Class());
+ }
+ if (fAnalysisStatus & AliCDMesonBase::kBitMultStudy) {
+ DefineOutput(iOutputSlot++, THnSparse::Class());
+ }
+ if (fAnalysisStatus & AliCDMesonBase::kBitTHnMC) {
+ DefineOutput(iOutputSlot++, THnSparse::Class());
+ }
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitPWAtree)) {
+ DefineOutput(iOutputSlot++, TTree::Class());
+ }
+ }
+ else { // create empty event THnSparse
+ DefineOutput(iOutputSlot++, THnSparse::Class());
+ }
+
+ if (fAnalysisStatus & AliCDMesonBase::kBitEEStudy) { // empty event study
+ fPhysicsSelection = new AliPhysicsSelection;
+ fPhysicsSelection->SetSkipTriggerClassSelection();
+ // accept all (A,C,E and I) triggers of a certain class
+ }
+
+ for (Int_t iGap = 0; iGap < kMax; ++iGap) {
+ fGapInformation[iGap] = 0;
+ }
+
+ // reset track pair pointers
+ fTrkPair[0] = 0x0;
+ fTrkPair[1] = 0x0;
+}
+
+
+//------------------------------------------------------------------------------
+AliAnalysisTaskCDMeson::AliAnalysisTaskCDMeson():
+ AliAnalysisTaskSE()
+ , fDoAOD(kFALSE)
+ , fAnalysisStatus(0x0)
+ , fAnalysisStatusString(0x0)
+ , fNoGapFraction(0.01)
+ , fReducedGapFraction(0.02)
+ , fMaxVtxDst(0.5)
+
+ , fESDEvent(0x0)
+ , fAODEvent(0x0)
+ , fPIDResponse(0x0)
+ , fPhysicsSelection(0x0)
+ , fTracks(0x0)
+ , fVtxDst(-1.)
+ , fVtxZ(-20)
+ , fResidualTracks(0)
+ , fResidualTracklets(0)
+ , fMCprocessType(0)
+ , fMCprocess(-1)
+
+ , fRun(-999)
+ , fPIDmode(0)
+
+ , fTheta(0.)
+ , fPhi(0.)
+ , fMass(0.)
+ , fCurrentGapCondition(0)
+
+ , fGapRun(0x0)
+ , fHist(0x0)
+ , fThnMother(0x0)
+ , fThnMotherSoft(0x0)
+ , fThnMultiplicity(0x0)
+ , fThnMC(0x0)
+ , fThnEmptyEvents(0x0)
+ , fPWAtree(0x0)
+
+ , fv0ntrk(0x0)
+ , fv0fmdntrk(0x0)
+ , fv0fmdspdntrk(0x0)
+ , fv0fmdspdtpcntrk(0x0)
+
+ , fv0Rmntrk(0x0)
+ , fv0fmdRmntrk(0x0)
+ , fv0fmdspdRmntrk(0x0)
+ , fv0fmdspdtpcRmntrk(0x0)
+
+ , fMultStudy(0x0)
+ , fMultStudyV0dg(0x0)
+ , fMultStudyV0FMDdg(0x0)
+ , fMultStudyV0FMDSPDdg(0x0)
+ , fMultStudyV0FMDSPDTPCdg(0x0)
+
+ , fMultResponseMC(0x0)
+ , fMultResponseV0dgMC(0x0)
+ , fMultResponseV0FMDdgMC(0x0)
+ , fMultResponseV0FMDSPDdgMC(0x0)
+ , fMultResponseV0FMDSPDTPCdgMC(0x0)
+
+ , fMultRegionsMC(0x0)
+
+ , fhspdV0dg(0x0)
+ , fhspdV0FMDdg(0x0)
+ , fhspdV0FMDSPDdg(0x0)
+ , fhspdV0FMDSPDTPCdg(0x0)
+ , fhspdAfterCuts(0x0)
+
+ , fGapResponseMCv0Dg(0x0)
+ , fGapResponseMCv0fmdDg(0x0)
+ , fGapResponseMCv0fmdspdDg(0x0)
+ , fGapResponseMCv0fmdspdtpcDg(0x0)
+
+ , fhspd(0x0)
+ , fhfo(0x0)
+ , fhpriVtxDist(0x0)
+ , fhpriVtxPos(0x0)
+ , fhpriv(0x0)
+ , fhntrk(0x0)
+ , fhStatsFlow(0x0)
+ , fhFOchans(0x0)
+
+ , fVZEROhists(0x0)
+
+ , fv0Map(0x0)
+ , fv0fmdMap(0x0)
+ , fv0fmdspdMap(0x0)
+ , fv0fmdspdtpcMap(0x0)
+
+ , fv0MapCutted(0x0)
+ , fv0fmdMapCutted(0x0)
+ , fv0fmdspdMapCutted(0x0)
+ , fv0fmdspdtpcMapCutted(0x0)
+
+ , fHitMapSPDinner(0x0)
+ , fHitMapSPDouter(0x0)
+ , fHitMapSPDtrklt(0x0)
+ , fHitMapFMDa(0x0)
+ , fHitMapFMDc(0x0)
+
+ , fFMDsum1I(0x0)
+ , fFMDsum2I(0x0)
+ , fFMDsum2O(0x0)
+ , fFMDsum3I(0x0)
+ , fFMDsum3O(0x0)
+
+ , fPriVtxX(0x0)
+ , fPriVtxY(0x0)
+ , fPriVtxZ(0x0)
+ , fPriVtxDst(0x0)
+ , fPriVtxDstV0dg(0x0)
+ , fPriVtxDstV0FMDdg(0x0)
+ , fPriVtxDstV0FMDSPDdg(0x0)
+ , fPriVtxDstV0FMDSPDTPCdg(0x0)
+
+ , fTPCGapDCAaSide(0x0)
+ , fTPCGapDCAcSide(0x0)
+
+ , fComb2trkPIDuls(0x0)
+ , fComb2trkPIDls(0x0)
+ , fComb2trkPIDulsDG(0x0)
+ , fComb2trkPIDlsDG(0x0)
+ , fComb2trkPIDulsMC(0x0)
+ , fComb2trkPIDlsMC(0x0)
+ , fComb2trkPIDulsDGmc(0x0)
+ , fComb2trkPIDlsDGmc(0x0)
+ , fMCProcessUls(0x0)
+ , fMCProcessLs(0x0)
+
+ , fMAllTrackMass(0x0)
+{
+ //
+ // default constructor (should not be used)
+ //
+
+ fTrkPair[0] = 0x0;
+ fTrkPair[1] = 0x0;
+}
+
+
+//------------------------------------------------------------------------------
+AliAnalysisTaskCDMeson::~AliAnalysisTaskCDMeson()
+{
+ //
+ //Destructor
+ //
+ //delete fESDEvent;
+ //delete fESDpid;
+
+ if (!(fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) { // normal operation
+ if ((!fAnalysisStatus || (fAnalysisStatus && AliCDMesonBase::kBitTHnMother))
+ && fThnMother
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fThnMother;
+ fThnMother = 0x0;
+ }
+ if ((!fAnalysisStatus
+ || ((fAnalysisStatus && AliCDMesonBase::kBitSoftTracks)
+ && (fAnalysisStatus && AliCDMesonBase::kBitTHnMother)))
+ && fThnMotherSoft
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fThnMotherSoft;
+ fThnMotherSoft = 0x0;
+ }
+ if ((!fAnalysisStatus
+ || (fAnalysisStatus && AliCDMesonBase::kBitMultStudy))
+ && fThnMultiplicity
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fThnMultiplicity;
+ fThnMultiplicity = 0x0;
+ }
+ if ((!fAnalysisStatus || (fAnalysisStatus && AliCDMesonBase::kBitTHnMC))
+ && fThnMC
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fThnMC;
+ fThnMC = 0x0;
+ }
+ }
+ else { // empty event study
+ if (fPhysicsSelection
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fPhysicsSelection;
+ fPhysicsSelection = 0x0;
+ }
+ if (fThnEmptyEvents
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fThnEmptyEvents;
+ fThnEmptyEvents = 0x0;
+ }
+ }
+
+ /*
+ if (fHist
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ fHist->Clear();
+ delete fHist;
+ fHist = 0x0;
+ }
+ */
+
+ if (fTracks
+ && (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() !=
+ AliAnalysisManager::kProofAnalysis)) {
+ delete fTracks;
+ fTracks = 0x0;
+ }
+
+ if (fVZEROhists) {
+ delete fVZEROhists;
+ fVZEROhists = 0x0;
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::UserCreateOutputObjects()
+{
+ //
+ //CreateOutputObjects
+ //
+
+ //= THnSparse ================================================================
+ if (!(fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) { // normal operation
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitTHnMother)) {
+ fThnMother = AliCDMesonBase::GetThnMother(/*"CDMeson_Mother"*/);
+ }
+ if (!fAnalysisStatus
+ || ((fAnalysisStatus & AliCDMesonBase::kBitSoftTracks)
+ && (fAnalysisStatus & AliCDMesonBase::kBitTHnMother))) {
+ fThnMotherSoft = AliCDMesonBase::GetThnMother("CDMeson_MotherSoft");
+ }
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitMultStudy)){
+ fThnMultiplicity = AliCDMesonBase::GetThnMultiplicity();
+ }
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitTHnMC)) {
+ fThnMC = AliCDMesonBase::GetThnMother("CDMeson_MotherMC");
+ }
+ }
+ else { // empty event studies
+ fThnEmptyEvents = AliCDMesonBase::GetThnEmptyEvents();
+ }
+
+ //= TList for Histograms =====================================================
+ fHist = new TList;
+ fHist->SetOwner(); // ensures that the histograms are all deleted on exit!
+
+ //fHist->Add(fAnalysisStatusString); // WARNING: CRASHES MERGING
+
+ //= GAP RUN =
+ if (!(fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) { // normal operation
+ // LHC10 data
+ // runa1=114649, runb1=117630, runc0=117631, runc1=121526, rund0=121527,
+ // rund1=126460, rune0=126461, rune1 = 130930; runf0=130931
+ // all checked on elog
+ // LHC11 data
+ //const Int_t rune0 = 160676;
+ //const Int_t runf1 = 165462;
+
+ // LHC10 + LHC11 data
+ const Int_t runb0 = 114650; // first run of LHC10b
+ const Int_t runf1 = 165462; // last run of LHC11f
+
+ const Int_t run0 = runb0-1, run1 = runf1 + 1;
+
+ const Int_t nrun = (Int_t)(run1-run0);
+
+ const Int_t bins[2] = {nrun, (Double_t)AliCDMesonBase::kBitGapMax - 1.};
+ const Double_t xmin[2] = {(Double_t)run0, 1.};
+ const Double_t xmax[2] = {
+ (Double_t)run1, (Double_t)AliCDMesonBase::kBitGapMax
+ };
+
+ fGapRun = new THnSparseI("GapRun", ";run number;gap condition", 2, bins,
+ xmin, xmax);
+ fHist->Add(fGapRun);
+ }
+
+ //= MULTIPLICITY PER GAP CONDITION =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitMultPerGapHists)) {
+ fv0ntrk = new TH2D("b00_v0ntrk", ";number of tracks;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0
+ fHist->Add(fv0ntrk);
+
+ fv0fmdntrk = new TH2D("b01_v0fmdntrk", ";number of tracks;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0FMD
+ fHist->Add(fv0fmdntrk);
+
+ fv0fmdspdntrk =
+ new TH2D("b02_v0fmdspdntrk", ";number of tracks;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0FMDSPD
+ fHist->Add(fv0fmdspdntrk);
+
+ fv0fmdspdtpcntrk =
+ new TH2D("b03_v0fmdspdtpcntrk", ";number of tracks;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0FMDSPDTPC
+ fHist->Add(fv0fmdspdtpcntrk);
+ }
+
+ //= MULTIPLICITY REMOVED PER GAP CONDITION =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitRmMultPerGapHists)) {
+ fv0Rmntrk = new TH2F("b05_v0Rmntrk",
+ ";number of removed tracks due to cuts;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0
+ fHist->Add(fv0Rmntrk);
+
+ fv0fmdRmntrk =
+ new TH2F("b06_v0fmdRmntrk",
+ ";number of removed tracks due to cuts;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0FMD
+ fHist->Add(fv0fmdRmntrk);
+
+ fv0fmdspdRmntrk =
+ new TH2F("b07_v0fmdspdRmntrk",
+ ";number of removed tracks due to cuts;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0FMDSPD
+ fHist->Add(fv0fmdspdRmntrk);
+
+ fv0fmdspdtpcRmntrk =
+ new TH2F("b08_v0fmdspdtpcRmntrk",
+ ";number of removed tracks due to cuts;gap condition",
+ 80, 0., 80., 4, 1., 5.);
+ //x: ntrk; y: V0FMDSPDTPC
+ fHist->Add(fv0fmdspdtpcRmntrk);
+ }
+
+ //= SOFT TRACK INFLUENCE =
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitSoftTracks)) {
+ fMultStudy = new TH2F("b10_MultStudySoftInfluence",
+ ";multiplicity;multiplicity (including ITSsa)",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultStudy);
+
+ fMultStudyV0dg =
+ new TH2F("b11_MultStudySoftInfluence_V0_DoubleGap",
+ ";multiplicity;multiplicity (including ITSsa)",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultStudyV0dg);
+
+ fMultStudyV0FMDdg =
+ new TH2F("b12_MultStudySoftInfluence_V0FMD_DoubleGap",
+ ";multiplicity;multiplicity (including ITSsa)",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultStudyV0FMDdg);
+
+ fMultStudyV0FMDSPDdg
+ = new TH2F("b13_MultStudySoftInfluence_V0FMDSPD_DoubleGap",
+ ";multiplicity;multiplicity (including ITSsa)",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultStudyV0FMDSPDdg);
+
+ fMultStudyV0FMDSPDTPCdg
+ = new TH2F("b14_MultStudySoftInfluence_V0FMDSPDTPC_DoubleGap",
+ ";multiplicity;multiplicity (including ITSsa)",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultStudyV0FMDSPDTPCdg);
+ }
+
+ //= MULTIPLICITY RESPONSE DATA VS MC =
+ if (fAnalysisStatus & AliCDMesonBase::kBitMultResponseMC) {
+ fMultResponseMC
+ = new TH2F("b15_MultResponseMC",
+ ";multiplicity from MC truth; reconstructed multiplicity",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultResponseMC);
+
+ fMultResponseV0dgMC
+ = new TH2F("b16_MultResponseV0dgMC",
+ ";multiplicity from MC truth; reconstructed multiplicity",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultResponseV0dgMC);
+
+ fMultResponseV0FMDdgMC
+ = new TH2F("b17_MultResponseV0FMDdgMC",
+ ";multiplicity from MC truth; reconstructed multiplicity",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultResponseV0FMDdgMC);
+
+ fMultResponseV0FMDSPDdgMC
+ = new TH2F("b18_MultResponseV0FMDSPDdgMC",
+ ";multiplicity from MC truth; reconstructed multiplicity",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultResponseV0FMDSPDdgMC);
+
+ fMultResponseV0FMDSPDTPCdgMC
+ = new TH2F("b19_MultResponseV0FMDSPDTPCdgMC",
+ ";multiplicity from MC truth; reconstructed multiplicity",
+ 20, 0., 20., 20, 0., 20.);
+ fHist->Add(fMultResponseV0FMDSPDTPCdgMC);
+
+ fMultRegionsMC
+ = new TH2F("b20_fMultRegionsMC", ";multiplicity;",
+ 50, 0., 50., 4, 0., 4.);
+ TAxis* a = fMultRegionsMC->GetYaxis();
+ a->SetBinLabel(1, "A-side");
+ a->SetBinLabel(2, "Center");
+ a->SetBinLabel(3, "C-side");
+ a->SetBinLabel(4, "A and C side");
+ fHist->Add(fMultRegionsMC);
+
+ fGapResponseMCv0Dg =
+ new TH1I("b26_fGapResponseMCv0Dg", ";;counts", 3, 0., 3.);
+ a = fGapResponseMCv0Dg->GetXaxis();
+ a->SetBinLabel(1, "DG in MC, no DG in data");
+ a->SetBinLabel(2, "DG in both MC and data");
+ a->SetBinLabel(3, "DG in data, no DG in MC");
+ fHist->Add(fGapResponseMCv0Dg);
+ fGapResponseMCv0fmdDg =
+ new TH1I("b27_fGapResponseMCv0fmdDg", ";;counts", 3, 0., 3.);
+ a = fGapResponseMCv0fmdDg->GetXaxis();
+ a->SetBinLabel(1, "DG in MC, no DG in data");
+ a->SetBinLabel(2, "DG in both MC and data");
+ a->SetBinLabel(3, "DG in data, no DG in MC");
+ fHist->Add(fGapResponseMCv0fmdDg);
+ fGapResponseMCv0fmdspdDg =
+ new TH1I("b28_fGapResponseMCv0fmdspdDg", ";;counts", 3, 0., 3.);
+ a = fGapResponseMCv0fmdspdDg->GetXaxis();
+ a->SetBinLabel(1, "DG in MC, no DG in data");
+ a->SetBinLabel(2, "DG in both MC and data");
+ a->SetBinLabel(3, "DG in data, no DG in MC");
+ fHist->Add(fGapResponseMCv0fmdspdDg);
+ fGapResponseMCv0fmdspdtpcDg
+ = new TH1I("b29_fGapResponseMCv0fmdspdtpcDg", ";;counts", 3, 0., 3.);
+ a = fGapResponseMCv0fmdspdtpcDg->GetXaxis();
+ a->SetBinLabel(1, "DG in MC, no DG in data");
+ a->SetBinLabel(2, "DG in both MC and data");
+ a->SetBinLabel(3, "DG in data, no DG in MC");
+ fHist->Add(fGapResponseMCv0fmdspdtpcDg);
+ }
+
+ //= FAST-OR MULTIPLICITY STUDY =
+ if (fAnalysisStatus & AliCDMesonBase::kBitFastORmultStudy) {
+ fhspdV0dg =
+ new TH1I("b21_spd_V0dg",
+ "V0 double gap;fired SPD FastOR chips; counts", 50, 0., 50.);
+ fHist->Add(fhspdV0dg);
+ fhspdV0FMDdg =
+ new TH1I("b22_spd_V0FMDdg",
+ "V0-FMD double gap;fired SPD FastOR chips; counts", 50, 0., 50.);
+ fHist->Add(fhspdV0FMDdg);
+ fhspdV0FMDSPDdg =
+ new TH1I("b23_spd_V0FMDSPDdg",
+ "V0-FMD-SPD double gap;fired SPD FastOR chips; counts",
+ 50, 0., 50.);
+ fHist->Add(fhspdV0FMDSPDdg);
+ fhspdV0FMDSPDTPCdg =
+ new TH1I("b24_spd_V0FMDSPDTPCdg",
+ "V0-FMD-SPD-TPC double gap;fired SPD FastOR chips; counts",
+ 50, 0., 50.);
+ fHist->Add(fhspdV0FMDSPDTPCdg);
+ fhspdAfterCuts =
+ new TH1I("b25_spd_afterCuts",
+ "after cuts;fired SPD FastOR chips; counts",
+ 50, 0., 50.);
+ fHist->Add(fhspdAfterCuts);
+ }
+
+ //= STATISTICS FLOW =
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ fhspd = new TH1I("a00_spd", ";fired SPD FastOR chips; counts", 50, 0., 50.);
+ fHist->Add(fhspd);
+
+ if (fAnalysisStatus & AliCDMesonBase::kBitFastORStudy) {
+ fhfo =
+ new TH2I("a05_fo", ";fired chips inner layer;fired chips outer layer",
+ 50, 0., 50., 50, 0., 50.);
+ fHist->Add(fhfo);
+
+ fhFOchans =
+ new TH1F("a06_foChans", ";channel key; counts", 1200, 0., 1200.);
+ fHist->Add(fhFOchans);
+ }
+
+ fhpriVtxDist = new TH1I("a08_priVtxDist", ";vertex position (cm); counts",
+ 601, -30.05, 30.05);
+ fHist->Add(fhpriVtxDist);
+
+ fhpriVtxPos = new TH1I("a09_priVtxPos", ";vertex position (boolean);counts",
+ 2, 0., 2.);
+ fhpriVtxPos->GetXaxis()->SetBinLabel(1,"inside");
+ fhpriVtxPos->GetXaxis()->SetBinLabel(2,"outside");
+ fHist->Add(fhpriVtxPos);
+
+
+ fhpriv = new TH1I("a10_priv", ";vertex quality (boolean);counts",
+ 2, 0., 2.);
+ fhpriv->GetXaxis()->SetBinLabel(1,"bad vertex");
+ fhpriv->GetXaxis()->SetBinLabel(2,"good vertex");
+ fHist->Add(fhpriv);
+
+ fhntrk = new TH1I("a20_ntrk", ";track multiplicity after cuts;counts",
+ 10, 0., 10.);
+ fHist->Add(fhntrk);
+
+ fhStatsFlow = AliCDMesonBase::GetHistStatsFlow();
+ fHist->Add(fhStatsFlow);
+ }
+
+ //= ETA-PHI MAPS FOR TRACKS =
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitEtaPhiMaps)) {
+ fv0Map = new TH2F("d0_v0Map", ";#eta;#phi", 40, -2., 2.,
+ 60, 0., TMath::TwoPi());
+ fHist->Add(fv0Map);
+
+ fv0fmdMap = new TH2F("d0_v0fmdMap", ";#eta;#phi", 40, -2., 2.,
+ 60, 0., TMath::TwoPi());
+ fHist->Add(fv0fmdMap);
+
+ fv0fmdspdMap = new TH2F("d0_v0fmdspdMap", ";#eta;#phi", 40, -2., 2.,
+ 60, 0., TMath::TwoPi());
+ fHist->Add(fv0fmdspdMap);
+
+ fv0fmdspdtpcMap = new TH2F("d0_v0fmdspdtpcMap", ";#eta;#phi", 30, -1.5, 1.5,
+ 60, 0., TMath::TwoPi());
+ fHist->Add(fv0fmdspdtpcMap);
+ }
+
+ //= ETA-PHI MAPS FOR CUTTED TRACKS =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitEtaPhiMapsWithCuts)) {
+ fv0MapCutted = new TH2F("d0_v0MapCutted", ";#eta;#phi", 24, -1.2, 1.2,
+ 60, 0., TMath::TwoPi());
+ fHist->Add(fv0MapCutted);
+
+ fv0fmdMapCutted = new TH2F("d0_v0fmdMapCutted", ";#eta;#phi", 24, -1.2, 1.2,
+ 60, 0., TMath::TwoPi());
+ fHist->Add(fv0fmdMapCutted);
+
+ fv0fmdspdMapCutted = new TH2F("d0_v0fmdspdMapCutted", ";#eta;#phi",
+ 24, -1.2, 1.2, 60, 0., TMath::TwoPi());
+ fHist->Add(fv0fmdspdMapCutted);
+
+ fv0fmdspdtpcMapCutted = new TH2F("d0_v0fmdspdtpcMapCutted", ";#eta;#phi",
+ 24, -1.2, 1.2, 60, 0., TMath::TwoPi());
+ fHist->Add(fv0fmdspdtpcMapCutted);
+ }
+
+ //= SPD HIT MAPS =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitHitMapSPD) ||
+ (fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) {
+ fHitMapSPDinner = new TH2F("d1_HitMapSPDinner", ";#eta;#phi", 72, -3.6, 3.6,
+ 20, 0., TMath::TwoPi());
+ fHist->Add(fHitMapSPDinner);
+
+ fHitMapSPDouter = new TH2F("d2_HitMapSPDouter", ";#eta;#phi", 48, -2.4, 2.4,
+ 40, 0., TMath::TwoPi());
+ fHist->Add(fHitMapSPDouter);
+
+ fHitMapSPDtrklt = new TH2F("d3_HitMapSPDtrklt", ";#eta;#phi", 60, -3., 3.,
+ 72, 0., TMath::TwoPi());
+ fHist->Add(fHitMapSPDtrklt);
+ }
+
+ //= FMD HIT MAPS =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitHitMapFMD) ||
+ (fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) {
+ fHitMapFMDa = new TH2F("d4_HitMapFMDa", ";#eta;#phi",
+ 36, 1.5, 5.1, 40, 0., TMath::TwoPi());
+ fHist->Add(fHitMapFMDa);
+ fHitMapFMDc = new TH2F("d4_HitMapFMDc", ";#eta;#phi",
+ 20, -3.5, -1.5, 40, 0., TMath::TwoPi());
+ fHist->Add(fHitMapFMDc);
+ }
+
+ //= FMD SUMMATION =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitFMDsum) ||
+ (fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) {
+ fFMDsum1I = new TH1F("d5_FMDsum1I", ";sum;counts", 50, 0., 50.);
+ fHist->Add(fFMDsum1I);
+ fFMDsum2I = new TH1F("d5_FMDsum2I", ";sum;counts", 50, 0., 50.);
+ fHist->Add(fFMDsum2I);
+ fFMDsum2O = new TH1F("d5_FMDsum2O", ";sum;counts", 50, 0., 50.);
+ fHist->Add(fFMDsum2O);
+ fFMDsum3I = new TH1F("d5_FMDsum3I", ";sum;counts", 50, 0., 50.);
+ fHist->Add(fFMDsum3I);
+ fFMDsum3O = new TH1F("d5_FMDsum3O", ";sum;counts", 50, 0., 50.);
+ fHist->Add(fFMDsum3O);
+ }
+
+ //= VERTEX STUDIES =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitVtxStudies)) {
+ fPriVtxX = new TH1I("e0_PriVtxX", ";x (cm);counts", 100, -10., 10.);
+ fHist->Add(fPriVtxX);
+ fPriVtxY = new TH1I("e0_PriVtxY", ";y (cm);counts", 100, -10., 10.);
+ fHist->Add(fPriVtxY);
+ fPriVtxZ = new TH1I("e0_PriVtxZ", ";z (cm);counts", 100, -10., 10.);
+ fHist->Add(fPriVtxZ);
+
+ fPriVtxDst = new TH1I("e1_PriVtxDst", ";dst (cm);counts", 101, 0., 2.525);
+ fHist->Add(fPriVtxDst);
+
+ fPriVtxDstV0dg =
+ new TH1I("e2_PriVtxDstV0dg", ";dst (cm);counts", 101, 0., 2.525);
+ fHist->Add(fPriVtxDstV0dg);
+
+ fPriVtxDstV0FMDdg =
+ new TH1I("e3_PriVtxDstV0FMDdg", ";dst (cm);counts", 101, 0., 2.525);
+ fHist->Add(fPriVtxDstV0FMDdg);
+
+ fPriVtxDstV0FMDSPDdg =
+ new TH1I("e4_PriVtxDstV0FMDSPDdg", ";dst (cm);counts", 101, 0., 2.525);
+ fHist->Add(fPriVtxDstV0FMDSPDdg);
+
+ fPriVtxDstV0FMDSPDTPCdg =
+ new TH1I("e5_PriVtxDstV0FMDSPDTPCdg", ";dst (cm);counts", 101, 0., 2.525);
+ fHist->Add(fPriVtxDstV0FMDSPDTPCdg);
+ }
+
+ //= TPC GAP STUDY =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitTPCGapStudy)) {
+ if (!fDoAOD) {
+ fTPCGapDCAaSide = new TH2F("i00_TPCGapDCAaSide", ";dca_{xy};dca_{z}",
+ 2000, 0., 1., 1000, 0., 5.);
+ fHist->Add(fTPCGapDCAaSide);
+ fTPCGapDCAcSide = new TH2F("i01_TPCGapDCAcSide", ";dca_{xy};dca_{z}",
+ 2000, 0., 1., 1000, 0., 5.);
+ fHist->Add(fTPCGapDCAcSide);
+ }
+ }
+
+ //= VZERO TRIGGER STUDY =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitVZEROStudy)) {
+ if (!fDoAOD) { // not yet possible for AODs
+ fVZEROhists = AliCDMesonBase::GetHistVZEROStudies(fHist);
+ }
+ }
+
+ //= PID MATRIX FOR DATA AND MC =
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitPIDStudy)) {
+ fComb2trkPIDuls = AliCDMesonBase::GetHistPIDStudies("f0_Comb2trkPIDuls");
+ fHist->Add(fComb2trkPIDuls);
+
+ fComb2trkPIDls = AliCDMesonBase::GetHistPIDStudies("f1_Comb2trkPIDls");
+ fHist->Add(fComb2trkPIDls);
+
+ fComb2trkPIDulsDG =
+ AliCDMesonBase::GetHistPIDStudies("f2_Comb2trkPIDulsDG");
+ fHist->Add(fComb2trkPIDulsDG);
+
+ fComb2trkPIDlsDG = AliCDMesonBase::GetHistPIDStudies("f3_Comb2trkPIDlsDG");
+ fHist->Add(fComb2trkPIDlsDG);
+
+ if (fAnalysisStatus & AliCDMesonBase::kBitTHnMC) {
+ fComb2trkPIDulsMC =
+ AliCDMesonBase::GetHistPIDStudies("f4_Comb2trkPIDulsMC");
+ fHist->Add(fComb2trkPIDulsMC);
+
+ fComb2trkPIDlsMC = AliCDMesonBase::GetHistPIDStudies("f5_Comb2trkPIDlsMC");
+ fHist->Add(fComb2trkPIDlsMC);
+
+ fComb2trkPIDulsDGmc =
+ AliCDMesonBase::GetHistPIDStudies("f6_Comb2trkPIDulsDGmc");
+ fHist->Add(fComb2trkPIDulsDGmc);
+
+ fComb2trkPIDlsDGmc =
+ AliCDMesonBase::GetHistPIDStudies("f7_Comb2trkPIDlsDGmc");
+ fHist->Add(fComb2trkPIDlsDGmc);
+ }
+ }
+
+ //= MC PROCESS IDs =
+ if (!fAnalysisStatus || fAnalysisStatus & AliCDMesonBase::kBitMCProcess) {
+ fMCProcessUls = new TH2F("g0_MCProcessUls",
+ ";MC Process ID;double-gap topology",
+ 105, 1., 106., 12, 0., 12.);
+ TAxis* axis = fMCProcessUls->GetYaxis();
+ axis->SetBinLabel(1, "V0");
+ axis->SetBinLabel(2, "FMD");
+ axis->SetBinLabel(3, "SPD");
+ axis->SetBinLabel(4, "TPC");
+ axis->SetBinLabel(5, "V0-FMD");
+ axis->SetBinLabel(6, "V0-FMD-SPD");
+ axis->SetBinLabel(7, "V0-FMD-SPD-TPC");
+ axis->SetBinLabel(8, "TPC-SPD");
+ axis->SetBinLabel(9, "TPC-SPD-FMD");
+ axis->SetBinLabel(10, "TPC-SPD-FMD-V0");
+ axis->SetBinLabel(11, "SPD-FMD");
+ axis->SetBinLabel(12, "SPD-FMD-V0");
+ fHist->Add(fMCProcessUls);
+
+ fMCProcessLs = new TH2F("g1_MCProcessLs",
+ ";MC Process ID;double-gap topology",
+ 105, 1., 106., 12, 0., 12.);
+ axis = fMCProcessLs->GetYaxis();
+ axis->SetBinLabel(1, "V0");
+ axis->SetBinLabel(2, "FMD");
+ axis->SetBinLabel(3, "SPD");
+ axis->SetBinLabel(4, "TPC");
+ axis->SetBinLabel(5, "V0-FMD");
+ axis->SetBinLabel(6, "V0-FMD-SPD");
+ axis->SetBinLabel(7, "V0-FMD-SPD-TPC");
+ axis->SetBinLabel(8, "TPC-SPD");
+ axis->SetBinLabel(9, "TPC-SPD-FMD");
+ axis->SetBinLabel(10, "TPC-SPD-FMD-V0");
+ axis->SetBinLabel(11, "SPD-FMD");
+ axis->SetBinLabel(12, "SPD-FMD-V0");
+ fHist->Add(fMCProcessLs);
+ }
+
+ //= ALL TRACK INVARIANT MASS =================================================
+ if (fAnalysisStatus & AliCDMesonBase::kBitAllTrackMass) {
+ fMAllTrackMass = new TH1F("j00_AllTrackInvMass",
+ ";all track invariant mass (GeV);counts",
+ 95, 0.25, 5.);
+ fHist->Add(fMAllTrackMass);
+ }
+ //= PWA TREE =================================================================
+ if (!fAnalysisStatus || fAnalysisStatus & AliCDMesonBase::kBitPWAtree) {
+ fPWAtree = new TTree("freidt_pwa", "pwa");
+ fPWAtree->Branch("theta", &fTheta);
+ fPWAtree->Branch("phi", &fPhi);
+ fPWAtree->Branch("m", &fMass);
+ fPWAtree->Branch("theta", &fTheta);
+ fPWAtree->Branch("pX", &fMomentum[0]);
+ fPWAtree->Branch("pY", &fMomentum[1]);
+ fPWAtree->Branch("pZ", &fMomentum[2]);
+ fPWAtree->Branch("gapCondition", &fCurrentGapCondition);
+ }
+
+ PostOutputs();
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::UserExec(Option_t *)
+{
+ //
+ // Executed for every event which passed the physics selection
+ //
+
+ //printf("Entry: %ld\n", (Long_t)Entry()); // print current event number
+
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinTotalInput); // stats flow
+ }
+
+ //= INPUT DATA SANITY TESTS ==================================================
+ if (!CheckInput()) {
+ PostOutputs();
+ return;
+ }
+
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinGoodInput); // stats flow
+ }
+
+ //= ANALYZE ONLY PHOJET CD EVENTS ============================================
+ if (fAnalysisStatus & AliCDMesonBase::kBitCDonly) {
+ if (fMCEvent) {
+ AliGenEventHeader* header = fMCEvent->GenEventHeader();
+ if (header) {
+ // Phojet
+ if (TString(header->IsA()->GetName()) == "AliGenDPMjetEventHeader") {
+ fMCprocess = ((AliGenDPMjetEventHeader*)header)->ProcessType();
+ if (fMCprocess != 4) { // no central diffraction
+ PostOutputs();
+ return;
+ }
+ else{ // CD found
+ fhStatsFlow->Fill(AliCDMesonBase::kBinCDonlyEvents); // statsflow
+ }
+ }
+ }
+ }
+ }
+
+
+ //= V0-AND/OR ================================================================
+ if (!fDoAOD) {
+ if (AliCDMesonUtils::V0AND(fESDEvent))
+ fhStatsFlow->Fill(AliCDMesonBase::kBinV0OR);
+ if (AliCDMesonUtils::V0OR(fESDEvent))
+ fhStatsFlow->Fill(AliCDMesonBase::kBinV0AND);
+ }
+
+ //= EMPTY EVENT STUDY ========================================================
+ // only implemented on ESDs!
+ if (fAnalysisStatus & AliCDMesonBase::kBitEEStudy) {
+ DoEmptyEventStudy();
+ PostOutputs();
+ return;
+ }
+
+ //= MC TRUTH =================================================================
+ // for now only implemented on ESDs
+ Int_t nMCprimaries = 0;
+ DetermineMCprocessType();
+ if (fAnalysisStatus & AliCDMesonBase::kBitTHnMC) {
+ nMCprimaries = DoMCTruth();
+ }
+
+ //= EVENT SELECTION ==========================================================
+ Int_t kfo =
+ (fAnalysisStatus & AliCDMesonBase::kBitFastORStudy & !fDoAOD) ? 1 : 0;
+ Int_t ninnerp=-999, nouterp=-999;
+ Bool_t eventIsValid = (fDoAOD) ?
+ AliCDMesonUtils::CutEvent(fAODEvent, fhpriv, fPriVtxX, fPriVtxY, fPriVtxZ,
+ fhpriVtxPos, fhpriVtxDist) :
+ AliCDMesonUtils::CutEvent(fESDEvent, fhspd, fhpriv, fhpriVtxPos,
+ fhpriVtxDist, fhfo, fhFOchans, kfo, ninnerp,
+ nouterp, fPriVtxX, fPriVtxY, fPriVtxZ);
+ if (!eventIsValid) {
+ PostOutputs();
+ return;
+ }
+
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinEventsAfterCuts); // stats flow
+ }
+
+ //= PILE UP ==================================================================
+ const Bool_t isPileup = (fDoAOD) ?
+ fAODEvent->IsPileupFromSPD(2, 0.8, 3., 2., 5.) :
+ fESDEvent->IsPileupFromSPD(2, 0.8, 3., 2., 5.);
+ // using only 2 instead of three contributors
+
+ if (isPileup) {
+ PostOutputs();
+ return;
+ }
+
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinEventsWithOutPileUp); // stats flow
+ }
+
+ //= VZERO TRIGGER STUDY ======================================================
+ // dtermine the VZERO signal distributions
+ if (fVZEROhists) {
+ AliCDMesonUtils::DoVZEROStudy(fESDEvent, fVZEROhists, fRun);
+ }
+
+ //= GAP ======================================================================
+ // determine the complete gap configuration (for all gap-tagging detectors)
+ if (!DetermineGap()) {
+ PostOutputs();
+ return;
+ }
+
+ // fill gaprun (determines the gap fraction per run)
+ const Double_t x[2] = {(Double_t)fRun, (Double_t)fCurrentGapCondition};
+ if (fGapRun) fGapRun->Fill(x);
+
+ //= VERTEX COINCIDENCE AND POSITION ==========================================
+ AnalyzeVtx();
+
+ //= TRACK CUTS ===============================================================
+ Bool_t doSoft = !fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitSoftTracks);
+
+ fTracks->ProcessEvent(fAODEvent, fESDEvent, doSoft); // apply cuts
+ DoMultiplicityStudy(nMCprimaries); // fill corresponding histograms
+ if (fMAllTrackMass &&
+ (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG)) {
+ // calculate the invariant mass of all tracks in the event including soft
+ // ones if the event is a full double gap event
+ fMAllTrackMass->Fill(fTracks->GetInvariantMass());
+ }
+
+ // is multiplicity within the desired range of 2 to 3?
+ Int_t nch = fTracks->GetTracks();
+ Int_t ncombined = fTracks->GetCombinedTracks();
+ Bool_t wSoft = ((ncombined >= 2) && (ncombined <=3)); // including soft tracks
+ Bool_t woSoft = ((nch >= 2) && (nch <= 3)); // excluding soft tracks
+
+ if (!wSoft && !woSoft) {
+ // multiplicity out of range (both with soft and without soft track)
+ PostOutputs();
+ return;
+ }
+
+ //= TRACK PAIRS ==============================================================
+ // loop over all track combinations
+ for(Int_t ii=0; ii<ncombined; ii++){
+ for(Int_t jj=ii+1; jj<ncombined; jj++){
+ // assign current tracks
+ fTrkPair[0] = fTracks->GetTrack(ii);
+ fTrkPair[1] = fTracks->GetTrack(jj);
+
+ // analyze track pairs
+ if (wSoft) DoTrackPair(kTRUE);
+ if (woSoft && (ii < nch) && (jj < nch)) DoTrackPair(kFALSE);
+ }
+ }
+
+ //============================================================================
+ PostOutputs();
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::PostOutputs()
+{
+ //
+ // PostData
+ //
+
+ Int_t iOutputSlot = 1; // dynamic output slot number handling
+ PostData(iOutputSlot++, fHist);
+ if (!(fAnalysisStatus & AliCDMesonBase::kBitEEStudy)) { // normal operation
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitTHnMother)) {
+ // THnSparse for invariant mass distributions
+ PostData(iOutputSlot++, fThnMother);
+ }
+ if (!fAnalysisStatus
+ || ((fAnalysisStatus & AliCDMesonBase::kBitSoftTracks)
+ && (fAnalysisStatus & AliCDMesonBase::kBitTHnMother))) {
+ // same including soft tracks
+ PostData(iOutputSlot++, fThnMotherSoft);
+ }
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitMultStudy)) {
+ // multiplicity study
+ PostData(iOutputSlot++, fThnMultiplicity);
+ }
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitTHnMC)) {
+ // MC truth study is active
+ PostData(iOutputSlot++, fThnMC);
+ }
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitPWAtree)) {
+ // PWA tree active
+ PostData(iOutputSlot++, fPWAtree);
+ }
+ }
+ else if (fAnalysisStatus & AliCDMesonBase::kBitEEStudy) {
+ // empty event study
+ PostData(iOutputSlot++, fThnEmptyEvents);
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::FillEtaPhiMaps() const
+{
+ //
+ // Controls which event contributes to which map
+ //
+
+ AliVEvent* event = (fDoAOD) ? (AliVEvent*)fAODEvent : (AliVEvent*)fESDEvent;
+
+ if (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0fmdspdtpcMap,
+ fv0fmdspdtpcMapCutted);
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0fmdspdMap,
+ fv0fmdspdMapCutted);
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0fmdMap, fv0fmdMapCutted);
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0Map, fv0MapCutted);
+ }
+ else if (fGapInformation[kV0FMDSPD] == AliCDMesonBase::kBinDG) {
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0fmdspdMap,
+ fv0fmdspdMapCutted);
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0fmdMap, fv0fmdMapCutted);
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0Map, fv0MapCutted);
+ }
+ else if (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG) {
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0fmdMap, fv0fmdMapCutted);
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0Map, fv0MapCutted);
+ }
+ else if (fGapInformation[kV0] == AliCDMesonBase::kBinDG) {
+ AliCDMesonUtils::FillEtaPhiMap(event, fTracks, fv0Map, fv0MapCutted);
+ }
+}
+
+
+//------------------------------------------------------------------------------
+Bool_t AliAnalysisTaskCDMeson::CheckInput()
+{
+ //
+ //general protection
+ //
+ if (const AliESDInputHandler *esdH =
+ dynamic_cast<AliESDInputHandler*>(fInputHandler)){
+ fESDEvent = esdH->GetEvent();
+ }
+ else if (const AliAODInputHandler *aodH =
+ dynamic_cast<AliAODInputHandler*>(fInputHandler)){
+ fAODEvent = aodH->GetEvent();
+ fDoAOD = kTRUE;
+ }
+ fNoGapFraction = (fDoAOD) ? 1. : fNoGapFraction; // process all running on AOD
+ fPIDResponse = (AliPIDResponse*)fInputHandler->GetPIDResponse();
+
+ if(!fESDEvent && !fAODEvent){
+ printf("freidtlog No valid event\n");
+ return kFALSE;
+ }
+
+ if(!fPIDResponse){
+ printf("freidtlog No PIDd\n");
+ // PID is fixed to unknown
+ //return kFALSE;
+ }
+
+ if(fDoAOD && fabs(fAODEvent->GetMagneticField())<1){
+ printf("freidtlog strange Bfield! %f\n", fAODEvent->GetMagneticField());
+ return kFALSE;
+ }
+ else if((!fDoAOD) && fabs(fESDEvent->GetMagneticField())<1){
+ printf("freidtlog strange Bfield! %f\n", fESDEvent->GetMagneticField());
+ return kFALSE;
+ }
+
+ Int_t tmprun = 0;
+ if (fDoAOD) {
+ tmprun = fAODEvent->GetRunNumber();
+ }
+ else {
+ tmprun = fESDEvent->GetRunNumber();
+ }
+
+ if(fRun!=tmprun){
+ fRun = tmprun;
+ AliCDMesonUtils::SPDLoadGeom(fRun);
+
+ // change PID cuts if necessary
+ if ((fRun >= 162933) && (fRun <= 165462)) {
+ fPIDmode = 1;
+ }
+ else {
+ fPIDmode = 0;
+ }
+ }
+
+ // get MC event
+ fMCEvent = MCEvent();
+
+ return kTRUE;
+}
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::DoEmptyEventStudy() {
+ // Analyses the observables needed for the gap determination in empty events
+ // for cross checks same observables are stored for A/C/I(B)-triggers as well
+ //
+
+
+ if (!fAnalysisStatus & AliCDMesonBase::kBitEEStudy) {
+ // check whether empty event analysis is activated
+ return;
+ }
+
+ if(!fESDEvent || fDoAOD) { // ensure that we are running on ESDs
+ return;
+ }
+
+ Int_t eventType = AliCDMesonUtils::GetEventType(fESDEvent);
+ TRandom3 rnd(0);
+ if ((eventType != AliCDMesonBase::kBinEventUnknown) &&
+ ((rnd.Rndm() > 0.95) ||
+ (AliCDMesonBase::kBinEventE == eventType))) {
+ Int_t fmdA = 0, fmdC = 0, spdIA = 0, spdIC = 0, spdOA = 0, spdOC = 0,
+ spdTrklSum = 0, spdTrklCentral = 0, spdTrklForwardA = 0,
+ spdTrklForwardC = 0;
+
+ Float_t fmdSums[] = { 0., 0., 0., 0., 0. };
+ AliCDMesonUtils::GetMultFMD(fESDEvent, fmdA, fmdC, fmdSums);
+ // Obtain FMD multiplicity using FMDHitCombinations
+ AliCDMesonUtils::GetMultSPD(fESDEvent, spdIA, spdIC, spdOA, spdOC);
+ // Obtain SPD multiplicity using FastOR information
+ AliCDMesonUtils::GetSPDTrackletMult(fESDEvent, spdTrklSum,
+ spdTrklForwardA, spdTrklForwardC,
+ spdTrklCentral);
+ // obtain SPD multiplicity using tracklets via AliMultiplicity
+
+ AliCDMesonBase::FillThnEmptyEvents(fThnEmptyEvents, eventType, fmdA, fmdC,
+ spdIA, spdIC, spdOA, spdOC,
+ spdTrklForwardA, spdTrklForwardC,
+ (Int_t)fmdSums[0], (Int_t)fmdSums[1],
+ (Int_t)fmdSums[2], (Int_t)fmdSums[3],
+ (Int_t)fmdSums[4]);
+ }
+ // the following code is only executed for empty events, this ensures, that
+ // all the histograms contain only information on empty events, while the
+ // THnSparse contains all information on A/C/I events
+ if (AliCDMesonBase::kBinEventE == eventType) {
+ TH1F* fmdSums[] = { fFMDsum1I, fFMDsum2I, fFMDsum2O, fFMDsum3I,
+ fFMDsum3O };
+ AliCDMesonUtils::GetGapConfig(fESDEvent, fHitMapSPDinner, fHitMapSPDouter,
+ fHitMapSPDtrklt, fHitMapFMDa, fHitMapFMDc,
+ (TH1**)fmdSums, fTPCGapDCAaSide,
+ fTPCGapDCAcSide);
+ // fill hit maps
+ }
+ return;
+}
+
+
+//------------------------------------------------------------------------------
+Bool_t AliAnalysisTaskCDMeson::DetermineGap()
+{
+ // determines the gap configuration for all gap tagging detectors based on the
+ // data set which is available
+ //
+
+ fCurrentGapCondition = 0x0; // initialize gap condition
+
+ if (fDoAOD) {
+ AliAODHandler* aodHandler =
+ (AliAODHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+ TTree *aodTree = aodHandler->GetTree();
+ if (aodTree) {
+ aodTree->SetBranchAddress("gapCondition", &fCurrentGapCondition);
+ aodTree->GetEvent(Entry()); // seems to be needed! (loads current event)
+ }
+ if (!fCurrentGapCondition) {
+ fCurrentGapCondition = 0xfffe;
+ puts("freidtlog - error while gap condition determination using AODs\n");
+ return kFALSE;
+ }
+ // DEBUGGING STUFF
+ //aodTree->ls();
+ //aodTree->MakeCode();
+ //aodTree->MakeClass();
+ //aodTree->GetEvent(Entry());
+ //aodTree->Show();
+ //TBranch* branch = aodTree->GetBranch("gapCondition");
+ //printf("freidtlog - branch=x%x\n",branch);
+ //if (!fCurrentGapCondition) {
+ //branch->SetAddress(&fCurrentGapCondition);
+ //}
+ //Int_t entry = Entry();
+ //printf("Entry()=%d\n", entry);
+ //branch->GetEvent(entry);
+ //printf("freidtlog - gapcondition=0x%x\n", fCurrentGapCondition);
+ }
+ else {
+ if (fAnalysisStatus & AliCDMesonBase::kBitReadPreprocessedGap) {
+ // retrieve preprocessed gap information
+ AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
+ TTree* esdTree = am->GetInputEventHandler()->GetTree(); // get ESD tree
+ if (esdTree) {
+ esdTree->SetBranchAddress("gapCondition", &fCurrentGapCondition);
+ esdTree->GetEvent(Entry()); // seems to be needed! (loads current event)
+ }
+ }
+ else {
+ // determine gaps from ESDs
+ TH1F* fmdSums[] = {fFMDsum1I, fFMDsum2I, fFMDsum2O, fFMDsum3I, fFMDsum3O};
+ fCurrentGapCondition = AliCDMesonUtils::GetGapConfig(fESDEvent,
+ fHitMapSPDinner,
+ fHitMapSPDouter,
+ fHitMapSPDtrklt,
+ fHitMapFMDa,
+ fHitMapFMDc,
+ (TH1**)fmdSums,
+ fTPCGapDCAaSide,
+ fTPCGapDCAcSide);
+ }
+ if (!fCurrentGapCondition) {
+ fCurrentGapCondition = 0xfffe;
+ puts("freidtlog - error while gap condition determination using ESDs\n");
+ return kFALSE;
+ }
+ }
+ //printf("freidtlog - Event: %ld, gapCondition: 0x%x\n", (Long_t)Entry(),
+ // fCurrentGapCondition);
+
+
+ // disentagle the contributions to the gap conditions of different "tightness"
+ fGapInformation[kV0] =
+ AliCDMesonBase::GetGapBin("V0", fCurrentGapCondition, kFALSE);
+ fGapInformation[kV0FMD] =
+ AliCDMesonBase::GetGapBin("V0FMD", fCurrentGapCondition, kFALSE);
+ fGapInformation[kV0FMDSPD] =
+ AliCDMesonBase::GetGapBin("V0FMDSPD", fCurrentGapCondition, kFALSE);
+ fGapInformation[kV0FMDSPDTPC] =
+ AliCDMesonBase::GetGapBin("V0FMDSPDTPC", fCurrentGapCondition, kFALSE);
+ fGapInformation[kFMD] =
+ AliCDMesonBase::GetGapBin("FMD",fCurrentGapCondition, kFALSE);
+ fGapInformation[kSPD] =
+ AliCDMesonBase::GetGapBin("SPD",fCurrentGapCondition, kFALSE);
+ fGapInformation[kTPC] =
+ AliCDMesonBase::GetGapBin("TPC",fCurrentGapCondition, kFALSE);
+ fGapInformation[kTPCSPD] =
+ AliCDMesonBase::GetGapBin("TPCSPD",fCurrentGapCondition, kFALSE);
+ fGapInformation[kTPCSPDFMD] =
+ AliCDMesonBase::GetGapBin("TPCSPDFMD",fCurrentGapCondition, kFALSE);
+ fGapInformation[kTPCSPDFMDV0] =
+ AliCDMesonBase::GetGapBin("TPCSPDFMDV0",fCurrentGapCondition, kFALSE);
+ fGapInformation[kSPDFMD] =
+ AliCDMesonBase::GetGapBin("SPDFMD",fCurrentGapCondition, kFALSE);
+ fGapInformation[kSPDFMDV0] =
+ AliCDMesonBase::GetGapBin("SPDFMDV0",fCurrentGapCondition, kFALSE);
+
+ fGapInformationWCent[kV0] = AliCDMesonBase::GetGapBin("V0", fCurrentGapCondition);
+ fGapInformationWCent[kV0FMD] =
+ AliCDMesonBase::GetGapBin("V0FMD", fCurrentGapCondition);
+ fGapInformationWCent[kV0FMDSPD] =
+ AliCDMesonBase::GetGapBin("V0FMDSPD", fCurrentGapCondition);
+ fGapInformationWCent[kV0FMDSPDTPC] =
+ AliCDMesonBase::GetGapBin("V0FMDSPDTPC", fCurrentGapCondition);
+ fGapInformationWCent[kFMD] =
+ AliCDMesonBase::GetGapBin("FMD",fCurrentGapCondition);
+ fGapInformationWCent[kSPD] =
+ AliCDMesonBase::GetGapBin("SPD",fCurrentGapCondition);
+ fGapInformationWCent[kTPC] =
+ AliCDMesonBase::GetGapBin("TPC",fCurrentGapCondition);
+ fGapInformationWCent[kTPCSPD] =
+ AliCDMesonBase::GetGapBin("TPCSPD",fCurrentGapCondition);
+ fGapInformationWCent[kTPCSPDFMD] =
+ AliCDMesonBase::GetGapBin("TPCSPDFMD",fCurrentGapCondition);
+ fGapInformationWCent[kTPCSPDFMDV0] =
+ AliCDMesonBase::GetGapBin("TPCSPDFMDV0",fCurrentGapCondition);
+ fGapInformationWCent[kSPDFMD] =
+ AliCDMesonBase::GetGapBin("SPDFMD",fCurrentGapCondition);
+ fGapInformationWCent[kSPDFMDV0] =
+ AliCDMesonBase::GetGapBin("SPDFMDV0",fCurrentGapCondition);
+
+ return kTRUE;
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::DoMultiplicityStudy(Int_t nMCprimaries)
+{
+ // stores the multiplicity distributions for different gap conditions and
+ // compares the multiplicity distribution for "normal primary tracks (ITSTPC)"
+ // and the multiplicity distribution furthermore taking also soft tracks into
+ // account
+ //
+
+ // retrieve values from the track object
+ Int_t ntrk0 = fTracks->GetTracksBeforeCuts();
+ Int_t nch = fTracks->GetTracks();
+ Int_t ncombined = fTracks->GetCombinedTracks();
+ Int_t nITSpureSA = fTracks->GetITSpureSACount();
+
+ // determine the residual tracks / tracklets
+ fResidualTracks = ntrk0 - ncombined - nITSpureSA;
+ fResidualTracklets = fTracks->GetRemainingTrackletsCentralBarrel();
+
+ // soft track specific stuff
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitSoftTracks)) {
+ fMultStudy->Fill((Double_t)nch, (Double_t)ncombined);
+ if (fGapInformation[kV0] == AliCDMesonBase::kBinDG) {
+ fMultStudyV0dg->Fill((Double_t)nch, (Double_t)ncombined);
+ }
+ if (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG) {
+ fMultStudyV0FMDdg->Fill((Double_t)nch, (Double_t)ncombined);
+ }
+ if (fGapInformation[kV0FMDSPD] == AliCDMesonBase::kBinDG) {
+ fMultStudyV0FMDSPDdg->Fill((Double_t)nch, (Double_t)ncombined);
+ }
+ if (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ fMultStudyV0FMDSPDTPCdg->Fill((Double_t)nch, (Double_t)ncombined);
+ }
+ }
+
+
+ // multiplicity distributions for different gaps
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitMultPerGapHists)) {
+ fv0ntrk->Fill(ncombined, fGapInformation[kV0]);
+ fv0fmdntrk->Fill(ncombined, fGapInformation[kV0FMD]);
+ fv0fmdspdntrk->Fill(ncombined, fGapInformation[kV0FMDSPD]);
+ fv0fmdspdtpcntrk->Fill(ncombined, fGapInformation[kV0FMDSPDTPC]);
+ }
+
+ // multiplicity removed by cuts for different gaps
+ if (!fAnalysisStatus ||
+ (fAnalysisStatus & AliCDMesonBase::kBitRmMultPerGapHists)) {
+ fv0Rmntrk->Fill(ntrk0-ncombined-nITSpureSA, fGapInformation[kV0]);
+ fv0fmdRmntrk->Fill(ntrk0-ncombined-nITSpureSA, fGapInformation[kV0FMD]);
+ fv0fmdspdRmntrk->Fill(ntrk0-ncombined-nITSpureSA,
+ fGapInformation[kV0FMDSPD]);
+ fv0fmdspdtpcRmntrk->Fill(ntrk0-ncombined-nITSpureSA,
+ fGapInformation[kV0FMDSPDTPC]);
+ }
+
+
+ // fill maps with eta phi information for QA purposes
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitEtaPhiMaps) ||
+ (fAnalysisStatus & AliCDMesonBase::kBitEtaPhiMapsWithCuts)) {
+ FillEtaPhiMaps();
+ }
+
+ // fill stats and general multiplicity histograms
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ fhntrk->Fill(ncombined); // multiplicity distribution
+
+ if (fGapInformationWCent[kV0] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinv0Gap);
+ }
+ if (fGapInformationWCent[kV0FMD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinv0fmdGap);
+ }
+ if (fGapInformationWCent[kV0FMDSPD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinv0fmdspdGap);
+ }
+ if (fGapInformationWCent[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinv0fmdspdtpcGap);
+ }
+ if (fGapInformationWCent[kTPC] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBintpcGap);
+ }
+ if (fGapInformationWCent[kTPCSPD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBintpcspdGap);
+ }
+ if (fGapInformationWCent[kTPCSPDFMD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBintpcspdfmdGap);
+ }
+ if (fGapInformationWCent[kTPCSPDFMDV0] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBintpcspdfmdv0Gap);
+ }
+ if (fGapInformationWCent[kFMD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinfmdGap);
+ }
+ if (fGapInformationWCent[kSPD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinspdGap);
+ }
+ if (fGapInformationWCent[kSPDFMD] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinspdfmdGap);
+ }
+ if (fGapInformationWCent[kSPDFMDV0] == AliCDMesonBase::kBinDG) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinspdfmdv0Gap);
+ }
+ if (nch == 2) fhStatsFlow->Fill(AliCDMesonBase::kBinTwoTrackEvents);
+ else if (nch == 3) fhStatsFlow->Fill(AliCDMesonBase::kBinThreeTrackEvents);
+ }
+
+ // fill multiplicity response histograms (data vs. MC)
+ if (fAnalysisStatus & AliCDMesonBase::kBitMultResponseMC) {
+ fMultResponseMC->Fill(nMCprimaries, ncombined);
+ if (fGapInformation[kV0] == AliCDMesonBase::kBinDG) {
+ fMultResponseV0dgMC->Fill(nMCprimaries, ncombined);
+ }
+ if (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG) {
+ fMultResponseV0FMDdgMC->Fill(nMCprimaries, ncombined);
+ }
+ if (fGapInformation[kV0FMDSPD] == AliCDMesonBase::kBinDG) {
+ fMultResponseV0FMDSPDdgMC->Fill(nMCprimaries, ncombined);
+ }
+ if (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ fMultResponseV0FMDSPDTPCdgMC->Fill(nMCprimaries, ncombined);
+ }
+ }
+
+ // event cleanliness
+ if (fResidualTracks == 0 && fhStatsFlow) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinResidualTracks);
+ }
+ if (fResidualTracklets == 0 && fhStatsFlow) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinResidualTracklets);
+ }
+
+
+ // multiplicity THnSparse
+ if (!fAnalysisStatus
+ || (fAnalysisStatus & AliCDMesonBase::kBitMultStudy)) {
+ AliCDMesonBase::FillThnMultiplicity(fThnMultiplicity, nch,
+ ncombined-nch, ncombined,
+ fGapInformation[kV0],
+ fGapInformation[kFMD],
+ fGapInformation[kSPD],
+ fGapInformation[kTPC],
+ fResidualTracks, fResidualTracklets,
+ fVtxZ, fVtxDst, fMCprocessType);
+ }
+
+ if (fAnalysisStatus & AliCDMesonBase::kBitFastORmultStudy) {
+ if (fhspdV0dg && fhspdV0FMDdg && fhspdV0FMDSPDdg && fhspdV0FMDSPDTPCdg &&
+ !fDoAOD) {
+ Int_t mult = AliCDMesonUtils::GetFastORmultiplicity(fESDEvent);
+
+ if (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ fhspdV0FMDSPDTPCdg->Fill(mult);
+ fhspdV0FMDSPDdg->Fill(mult);
+ fhspdV0FMDdg->Fill(mult);
+ fhspdV0dg->Fill(mult);
+ fhspdAfterCuts->Fill(mult);
+ }
+ else if (fGapInformation[kV0FMDSPD] == AliCDMesonBase::kBinDG) {
+ fhspdV0FMDSPDdg->Fill(mult);
+ fhspdV0FMDdg->Fill(mult);
+ fhspdV0dg->Fill(mult);
+ fhspdAfterCuts->Fill(mult);
+ }
+ else if (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG) {
+ fhspdV0FMDdg->Fill(mult);
+ fhspdV0dg->Fill(mult);
+ fhspdAfterCuts->Fill(mult);
+ }
+ else if (fGapInformation[kV0] == AliCDMesonBase::kBinDG) {
+ fhspdV0dg->Fill(mult);
+ fhspdAfterCuts->Fill(mult);
+ }
+ else {
+ fhspdAfterCuts->Fill(mult);
+ }
+ }
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::DoTrackPair(Bool_t soft /* = kFALSE */)
+{
+ //
+ // process a two track pair
+ //
+
+ AliCDMesonUtils::SwapTrack((const AliVTrack**)fTrkPair); // random order
+
+ const Int_t combCh =
+ AliCDMesonUtils::GetCombCh((const AliVTrack**)fTrkPair); // uls or ls?
+
+ const Int_t combPID = DoPID(combCh); // retrieve PID information
+
+ // analyze vertex quality
+ Int_t vtxInRng = 0;
+ Int_t vtxCoincidence = 0;
+
+ vtxInRng = (fabs(fVtxZ) < 4.) ? 1 : 0;
+ // 4cm range, determined from detector geometry
+
+ vtxCoincidence = (fVtxDst < fMaxVtxDst) ? 1 : 0;
+ // vertex coincidence of the track and SPD vertex
+
+ // initialize kinematic values
+ Double_t mass = -999., pt = -999., cts = -999., oa = -999.;
+ const TVector3 p1(fTrkPair[0]->Px(), fTrkPair[0]->Py(), fTrkPair[0]->Pz());
+ const TVector3 p2(fTrkPair[1]->Px(), fTrkPair[1]->Py(), fTrkPair[1]->Pz());
+ const TVector3* momenta[2] = { &p1, &p2 };
+ AliCDMesonUtils::GetMassPtCtsOA(combPID, momenta, mass, pt, cts, oa);
+
+ //= THnMother ================================================================
+ if (!fAnalysisStatus ||
+ ((fAnalysisStatus & AliCDMesonBase::kBitTHnMother)
+ && (!soft || fAnalysisStatus & AliCDMesonBase::kBitSoftTracks))) {
+ TRandom3 rnd(0);
+ if (((fGapInformation[kV0] == AliCDMesonBase::kBinDG)
+ && (!(fAnalysisStatus & AliCDMesonBase::kBitReduceGapEvents)
+ || (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG)
+ || (rnd.Rndm() > (1. - fReducedGapFraction))))
+ || (rnd.Rndm() > (1. - fNoGapFraction)) || fDoAOD
+ || (vtxInRng && vtxCoincidence && (fResidualTracks == 0)
+ && (fResidualTracklets == 0))) {
+ // reduce amount of no-gap events by a factor of 10 in order to spare
+ // memory and shrink the final data volume
+
+ // do switching needed to process either events including soft or without
+ THnSparseD* thn = (soft) ? fThnMotherSoft : fThnMother; // output variable
+ Int_t mult = (soft) ? fTracks->GetCombinedTracks() : fTracks->GetTracks();
+ // correct multiplicity value
+ Int_t residualTracks = (soft) ?
+ fResidualTracks : fResidualTracks + fTracks->GetSoftTracks();
+ // if the analysis is done without soft tracks, they will be added to the
+ // residualTracks
+
+ // fill MC histograms
+ if (fMCEvent && (!fAnalysisStatus ||
+ fAnalysisStatus & AliCDMesonBase::kBitMCProcess)) {
+ if (mult == 2) {
+ FillMChists(combCh);
+ }
+ }
+
+ AliCDMesonBase::FillThnMother(thn, (Double_t)mult, (Double_t)combCh,
+ (Double_t)combPID,
+ (Double_t)fGapInformation[kV0],
+ (Double_t)fGapInformation[kFMD],
+ (Double_t)fGapInformation[kSPD],
+ (Double_t)fGapInformation[kTPC],
+ mass, pt, cts, oa, fTrkPair[0]->Pt(),
+ (residualTracks == 0) ? 0. : 1.,
+ (Double_t)vtxInRng,
+ (Double_t)fMCprocessType,
+ (Double_t)vtxCoincidence,
+ (fResidualTracklets == 0) ? 0. : 1.);
+ }
+ }
+
+ //= PWA ======================================================================
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitPWAtree)) {
+ if ((soft && (fTracks->GetCombinedTracks() == 2))
+ || ((fTracks->GetTracks() == 2)
+ && !(!fAnalysisStatus
+ || (fAnalysisStatus & AliCDMesonBase::kBitSoftTracks)))) {
+ // TODO refine this condition
+ // TODO add flag for 3 track events and soft tracks and the other quality
+ // flags
+ // PWA is only done for two track events
+ AliCDMesonUtils::GetPWAinfo(combPID, (const AliVTrack**)fTrkPair, fTheta,
+ fPhi, fMass, fMomentum);
+ TRandom3 rnd(0);
+ if (((fGapInformation[kV0] == AliCDMesonBase::kBinDG)
+ && (!(fAnalysisStatus & AliCDMesonBase::kBitReduceGapEvents)
+ || (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG)
+ || (rnd.Rndm() > (1. - fReducedGapFraction)))) ||
+ (rnd.Rndm() > (1. - fNoGapFraction))) {
+ fPWAtree->Fill();
+ }
+ }
+ }
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliAnalysisTaskCDMeson::DoPID(Int_t combCh)
+{
+ // determine the PID of the two tracks for full double gap events, store
+ // the results in another histogram than for the remaining events
+ //
+
+ Int_t combPID = 0;
+ if (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ // double-gap events
+ combPID = (combCh == AliCDMesonBase::kBinPM) ? // assignment to uls/ls histo
+ AliCDMesonUtils::GetCombPID(fPIDResponse, (const AliVTrack**)fTrkPair,
+ fPIDmode, fComb2trkPIDulsDG) :
+ AliCDMesonUtils::GetCombPID(fPIDResponse, (const AliVTrack**)fTrkPair,
+ fPIDmode, fComb2trkPIDlsDG);
+ }
+ else { // non full double-gap events
+ combPID = (combCh == AliCDMesonBase::kBinPM) ? // assignment to uls/ls histo
+ AliCDMesonUtils::GetCombPID(fPIDResponse, (const AliVTrack**)fTrkPair,
+ fPIDmode, fComb2trkPIDuls) :
+ AliCDMesonUtils::GetCombPID(fPIDResponse, (const AliVTrack**)fTrkPair,
+ fPIDmode, fComb2trkPIDls);
+ }
+
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitStatsFlow)) {
+ if (combPID == AliCDMesonBase::kBinPion ||
+ combPID == AliCDMesonBase::kBinPionE) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinPionEvents);
+ }
+ else if (combPID == AliCDMesonBase::kBinKaon ||
+ combPID == AliCDMesonBase::kBinKaonE) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinKaonEvents);
+ }
+ else if (combPID == AliCDMesonBase::kBinProton ||
+ combPID == AliCDMesonBase::kBinProtonE) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinProtonEvents);
+ }
+ else if (combPID == AliCDMesonBase::kBinElectron ||
+ combPID == AliCDMesonBase::kBinElectronE) {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinElectronEvents);
+ }
+ else {
+ fhStatsFlow->Fill(AliCDMesonBase::kBinUnknownPIDEvents);
+ }
+ }
+
+ return combPID;
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::DetermineMCprocessType()
+{
+ //
+ // retrieves the MC process type from the AliGenEventHeader and classifies
+ // them
+ //
+
+ // get MC information
+ fMCprocess = -1; //detailed MC sub process information
+ fMCprocessType = AliCDMesonBase::kBinND; // ND is default, also for data
+
+ if (fMCEvent) {
+ AliGenEventHeader* header = fMCEvent->GenEventHeader();
+ if (header) {
+ // Pythia6
+ if (TString(header->IsA()->GetName()) == "AliGenPythiaEventHeader") {
+ fMCprocess = ((AliGenPythiaEventHeader*)header)->ProcessType();
+ switch(fMCprocess) {
+ case 92:
+ case 93:
+ case 103:
+ case 104: fMCprocessType = AliCDMesonBase::kBinSD; break;
+ case 94:
+ case 105: fMCprocessType = AliCDMesonBase::kBinDD; break;
+ default: fMCprocessType = AliCDMesonBase::kBinND; break;
+ }
+ }
+ // Phojet
+ else if (TString(header->IsA()->GetName()) == "AliGenDPMjetEventHeader") {
+ fMCprocess = ((AliGenDPMjetEventHeader*)header)->ProcessType();
+ switch(fMCprocess) {
+ case 5:
+ case 6: fMCprocessType = AliCDMesonBase::kBinSD; break;
+ case 7: fMCprocessType = AliCDMesonBase::kBinDD; break;
+ case 4: fMCprocessType = AliCDMesonBase::kBinCD; break;
+ default: fMCprocessType = AliCDMesonBase::kBinND; break;
+ }
+ }
+ }
+ }
+}
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::FillMChists(Int_t combCh)
+{
+ // fill the MC histograms for different gap conditions. Only real two track
+ // events are taken into account
+ //
+
+ if (fMCprocess > 0) { // MC process
+ for (Int_t iGap = kV0; iGap < kMax; ++iGap) {
+ if (fGapInformation[iGap] == AliCDMesonBase::kBinDG) {
+ if (combCh == AliCDMesonBase::kBinPM) {
+ fMCProcessUls->Fill(fMCprocess, iGap);
+ }
+ else if (combCh == AliCDMesonBase::kBinPPMM) {
+ fMCProcessLs->Fill(fMCprocess, iGap);
+ }
+ }
+ }
+ }
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliAnalysisTaskCDMeson::DoMCTruth()
+{
+ //
+ // analyses the MC truth and does a gap selection based on the eta values of
+ // the primaries
+ //
+
+ if (!fMCEvent) return -1;
+ AliStack* stack = fMCEvent->Stack();
+ if (!stack) return -1;
+
+ //= Multiplicity =============================================================
+ // determine number of charged physical primaries on the stack
+ Int_t nPhysicalPrimaries = 0;
+ Int_t nPrimaries = stack->GetNprimary();
+ for (Int_t iTracks = 0; iTracks < nPrimaries; ++iTracks) {
+ TParticle* part = stack->Particle(iTracks);
+ if (stack->IsPhysicalPrimary(iTracks) && (part->GetPDG()->Charge() != 0.) &&
+ (part->Eta() < 0.9) && (part->Eta() > -0.9)) { // TODO add parameter
+ ++nPhysicalPrimaries;
+ }
+ }
+
+
+ //= Track Gap && Multiplicity in Region Bins =================================
+ Int_t gapA = 0;
+ Int_t gapAv0 = 0;
+ Int_t gapAv0fmd = 0;
+ Int_t gapC = 0;
+ Int_t gapCv0 = 0;
+ Int_t gapCv0fmd = 0;
+ Int_t central = 0;
+ for (Int_t iTracks = 0; iTracks < nPrimaries; ++iTracks) {
+ TParticle* part = (TParticle*)stack->Particle(iTracks);
+ if (part && stack->IsPhysicalPrimary(iTracks) &&
+ (part->GetPDG()->Charge() != 0.)) {
+ if (part->Eta() > -0.9 && part->Eta() < 0.9) central++;
+ if (part->Eta() > 0.9 && part->Eta() < 5.1) gapA++;
+ if (part->Eta() > 2.8 && part->Eta() < 5.1) gapAv0++;
+ if (part->Eta() > 1.7 && part->Eta() < 5.1) gapAv0fmd++;
+ if (part->Eta() < -0.9 && part->Eta() > -3.7) gapC++;
+ if (part->Eta() < -1.9 && part->Eta() > -3.7) gapCv0++;
+ if (part->Eta() < -1.9 && part->Eta() > -3.7) gapCv0fmd++;
+ }
+ }
+
+ if (fMultRegionsMC) {
+ // multiplicity distribution separated in A-side, central barrel and C-side
+ fMultRegionsMC->Fill(gapA, 0);
+ fMultRegionsMC->Fill(central, 1);
+ fMultRegionsMC->Fill(gapC, 2);
+ fMultRegionsMC->Fill(gapA+gapC, 3);
+ }
+
+ // WARNING: gap response determination based on primaries only, method too
+ // simple to obtain reliable results; at least decays should be taken into
+ // account properly
+ if (fGapResponseMCv0Dg) {
+ Bool_t gapFromMC = (gapAv0 == 0) && (gapCv0 == 0);
+ Bool_t gapFromData = (fGapInformation[kV0] == AliCDMesonBase::kBinDG);
+
+ if (gapFromMC && !gapFromData) {
+ fGapResponseMCv0Dg->Fill(0);
+ }
+ else if (gapFromMC && gapFromData) {
+ fGapResponseMCv0Dg->Fill(1);
+ }
+ else if (!gapFromMC && gapFromData) {
+ fGapResponseMCv0Dg->Fill(2);
+ }
+ }
+ if (fGapResponseMCv0fmdDg) {
+ Bool_t gapFromMC = (gapAv0fmd == 0) && (gapCv0fmd == 0);
+ Bool_t gapFromData = (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG);
+
+ if (gapFromMC && !gapFromData) {
+ fGapResponseMCv0fmdDg->Fill(0);
+ }
+ else if (gapFromMC && gapFromData) {
+ fGapResponseMCv0fmdDg->Fill(1);
+ }
+ else if (!gapFromMC && gapFromData) {
+ fGapResponseMCv0fmdDg->Fill(2);
+ }
+ }
+ if (fGapResponseMCv0fmdspdDg) {
+ Bool_t gapFromMC = (gapA == 0) && (gapC == 0);
+ Bool_t gapFromData = (fGapInformation[kV0FMDSPD] == AliCDMesonBase::kBinDG);
+
+ if (gapFromMC && !gapFromData) {
+ fGapResponseMCv0fmdspdDg->Fill(0);
+ }
+ else if (gapFromMC && gapFromData) {
+ fGapResponseMCv0fmdspdDg->Fill(1);
+ }
+ else if (!gapFromMC && gapFromData) {
+ fGapResponseMCv0fmdspdDg->Fill(2);
+ }
+ }
+ if (fGapResponseMCv0fmdspdtpcDg) {
+ Bool_t gapFromMC = (gapA == 0) && (gapC == 0);
+ Bool_t gapFromData =
+ (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG);
+
+ if (gapFromMC && !gapFromData) {
+ fGapResponseMCv0fmdspdtpcDg->Fill(0);
+ }
+ else if (gapFromMC && gapFromData) {
+ fGapResponseMCv0fmdspdtpcDg->Fill(1);
+ }
+ else if (!gapFromMC && gapFromData) {
+ fGapResponseMCv0fmdspdtpcDg->Fill(2);
+ }
+ }
+
+ if ((nPhysicalPrimaries < 2 ) || (nPhysicalPrimaries >= 4 ))
+ return nPhysicalPrimaries;
+
+ Int_t gapCond = AliCDMesonBase::kBinNG;
+ if ((gapA == 0) && (gapC == 0) && (central > 0))
+ gapCond = AliCDMesonBase::kBinDG;
+ else if ((gapA == 0) && (gapC > 0)) gapCond = AliCDMesonBase::kBinGA;
+ else if ((gapA > 0) && (gapC == 0)) gapCond = AliCDMesonBase::kBinGC;
+ else if ((gapA > 0) && (gapC > 0)) gapCond = AliCDMesonBase::kBinNG;
+
+ //= Two Track ================================================================
+ for (Int_t iTracks = 0; iTracks < stack->GetNprimary(); ++iTracks) {
+ TParticle* part1 =0x0;
+ part1 = (TParticle*)stack->Particle(iTracks);
+ if (part1) {
+ if (stack->IsPhysicalPrimary(iTracks) &&
+ (part1->GetPDG()->Charge() != 0.) && (part1->Eta() < 0.9) &&
+ (part1->Eta() > -0.9)) {
+ for (Int_t jTracks = iTracks; jTracks < stack->GetNprimary(); ++jTracks) {
+ TParticle* part2 = 0x0;
+ part2 = (TParticle*)stack->Particle(jTracks);
+ if (part2) {
+ if (stack->IsPhysicalPrimary(jTracks) &&
+ (part2->GetPDG()->Charge() != 0.) && (part2->Eta() < 0.9) &&
+ (part2->Eta() > -0.9)) {
+ TParticle* particles[2] = { part1, part2 };
+ DoMCTrackPair(particles, gapCond, nPhysicalPrimaries);
+ }
+ }
+ }
+ }
+ }
+ }
+ return nPhysicalPrimaries;
+}
+
+
+//------------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::DoMCTrackPair(TParticle* particles[],
+ Int_t gapCond, Int_t multiplicity)
+{
+ //
+ // analyses a two track pair from MC truth
+ //
+
+ // swap tracks randomly
+ TRandom3 rnd(0);
+ if (rnd.Rndm() >= 0.5) {
+ TParticle* t = particles[1];
+ particles[1] = particles[0];
+ particles[0] = t;
+ }
+
+ // determine whether the pair is like sign or unlike sign
+ const Int_t combCh =
+ (particles[0]->GetPDG()->Charge()*particles[1]->GetPDG()->Charge() > 0.) ?
+ AliCDMesonBase::kBinPPMM : AliCDMesonBase::kBinPM;
+
+ // determine PID
+ Int_t combPID = 0;
+ if (gapCond == AliCDMesonBase::kBinDG) {
+ combPID = (combCh == AliCDMesonBase::kBinPM) ? // assignment to uls/ls histo
+ AliCDMesonUtils::GetCombPID((const TParticle**)particles,
+ fComb2trkPIDulsDGmc) :
+ AliCDMesonUtils::GetCombPID((const TParticle**)particles,
+ fComb2trkPIDlsDGmc);
+ }
+ else {
+ combPID = (combCh == AliCDMesonBase::kBinPM) ? // assignment to uls/ls histo
+ AliCDMesonUtils::GetCombPID((const TParticle**)particles,
+ fComb2trkPIDulsMC) :
+ AliCDMesonUtils::GetCombPID((const TParticle**)particles,
+ fComb2trkPIDlsMC);
+ }
+
+ // vertex quality is always nice for MC truth ;)
+ Int_t vtxInRng = 1;
+ Int_t vtxCoincidence = 1;
+
+ // initialize kinematic values
+ Double_t mass = -999., pt = -999., cts = -999., oa = -999.;
+ const TVector3 p1(particles[0]->Px(), particles[0]->Py(), particles[0]->Pz());
+ const TVector3 p2(particles[1]->Px(), particles[1]->Py(), particles[1]->Pz());
+ const TVector3* momenta[2] = { &p1, &p2 };
+ AliCDMesonUtils::GetMassPtCtsOA(combPID, momenta, mass, pt, cts, oa);
+
+ //= THnMother ================================================================
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitTHnMother)) {
+ if ((gapCond != AliCDMesonBase::kBinNG) ||
+ (rnd.Rndm() > (1. - fNoGapFraction)) || fDoAOD) {
+
+ AliCDMesonBase::FillThnMother(fThnMC, multiplicity, (Double_t)combCh,
+ (Double_t)combPID, (Double_t)gapCond,
+ (Double_t)gapCond, (Double_t)gapCond,
+ (Double_t)gapCond, mass, pt, cts, oa,
+ particles[0]->Pt(), 0., (Double_t)vtxInRng,
+ (Double_t)fMCprocessType,
+ (Double_t)vtxCoincidence, 1.);
+ }
+ }
+
+ //= PWA ======================================================================
+ // TODO ? should this be implemented???
+}
+
+
+//--------------------------------------------------------------------------
+void AliAnalysisTaskCDMeson::AnalyzeVtx()
+{
+ // calculates the distance between the vertex obtain from tracks and the
+ // vertex obtain from spd tracklets
+ // stores the z position of the primary vertex from tracks
+
+ fVtxDst = 0.; // reset distance
+
+ // retrieve the pointers of the current primary vertices
+ AliVVertex* trackVtx = (fDoAOD) ?
+ (AliVVertex*)fAODEvent->GetPrimaryVertex() :
+ (AliVVertex*)fESDEvent->GetPrimaryVertexTracks();
+ AliVVertex* spdVtx = (fDoAOD) ?
+ (AliVVertex*)fAODEvent->GetPrimaryVertexSPD() :
+ (AliVVertex*)fESDEvent->GetPrimaryVertexSPD();
+
+ fVtxZ = trackVtx->GetZ(); // store the vertex z position
+
+ if (fDoAOD && (trackVtx == spdVtx)) { // no primary track vertex in the AOD
+ fVtxDst = -5; // set arbitrary distance (counted in underflow bin!)
+ }
+ else { // do proper calculation of the geometrical distance
+ fVtxDst += (trackVtx->GetX()-spdVtx->GetX())
+ * (trackVtx->GetX()-spdVtx->GetX());
+ fVtxDst += (trackVtx->GetY()-spdVtx->GetY())
+ * (trackVtx->GetY()-spdVtx->GetY());
+ fVtxDst += (trackVtx->GetZ()-spdVtx->GetZ())
+ * (trackVtx->GetZ()-spdVtx->GetZ());
+ fVtxDst = TMath::Sqrt(fVtxDst);
+ }
+
+ if (!fAnalysisStatus || (fAnalysisStatus & AliCDMesonBase::kBitVtxStudies)) {
+ fPriVtxDst->Fill(fVtxDst);
+ if (fGapInformation[kV0FMDSPDTPC] == AliCDMesonBase::kBinDG) {
+ fPriVtxDstV0FMDSPDTPCdg->Fill(fVtxDst);
+ fPriVtxDstV0FMDSPDdg->Fill(fVtxDst);
+ fPriVtxDstV0FMDdg->Fill(fVtxDst);
+ fPriVtxDstV0dg->Fill(fVtxDst);
+ }
+ else if (fGapInformation[kV0FMDSPD] == AliCDMesonBase::kBinDG) {
+ fPriVtxDstV0FMDSPDdg->Fill(fVtxDst);
+ fPriVtxDstV0FMDdg->Fill(fVtxDst);
+ fPriVtxDstV0dg->Fill(fVtxDst);
+ }
+ else if (fGapInformation[kV0FMD] == AliCDMesonBase::kBinDG) {
+ fPriVtxDstV0FMDdg->Fill(fVtxDst);
+ fPriVtxDstV0dg->Fill(fVtxDst);
+ }
+ else if (fGapInformation[kV0] == AliCDMesonBase::kBinDG) {
+ fPriVtxDstV0dg->Fill(fVtxDst);
+ }
+ }
+}
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+//
+// Select events according to gap conditions, analyze two track events in pp
+// collisions
+//
+// Author:
+// Xianguo Lu <lu@physi.uni-heidelberg.de>
+// continued by
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#ifndef ALIANALYSISTASKCDMESON_H
+#define ALIANALYSISTASKCDMESON_H
+
+#ifndef ALIANALYSISTASK_H
+#include "AliAnalysisTaskSE.h"
+#endif
+
+#include "THnSparse.h" // forward declaration is not possible
+
+class AliESDEvent;
+class AliVTrack;
+class AliPIDResponse;
+class AliPhysicsSelection;
+
+class TH1I;
+class TH1F;
+class TH2I;
+class TH2F;
+class TH2D;
+class TList;
+class TParticle;
+class THnSparse;
+class TObjString;
+class TObjArray;
+
+class AliCDMesonTracks;
+
+class AliAnalysisTaskCDMeson : public AliAnalysisTaskSE
+{
+public:
+ AliAnalysisTaskCDMeson(const char* name, Long_t state = 0x0);
+ AliAnalysisTaskCDMeson();
+ virtual ~AliAnalysisTaskCDMeson();
+
+ virtual void UserCreateOutputObjects();
+ virtual void UserExec(Option_t *);
+
+private:
+ enum { //
+ kV0 = 0,
+ kFMD,
+ kSPD,
+ kTPC,
+ kV0FMD,
+ kV0FMDSPD,
+ kV0FMDSPDTPC,
+ kTPCSPD,
+ kTPCSPDFMD,
+ kTPCSPDFMDV0,
+ kSPDFMD,
+ kSPDFMDV0,
+ kMax
+ };
+
+ AliAnalysisTaskCDMeson(const AliAnalysisTaskCDMeson &p);
+ AliAnalysisTaskCDMeson& operator=(const AliAnalysisTaskCDMeson &p);
+
+ void FillEtaPhiMaps() const; // controls which event contributes to which map
+
+
+ // functions called by the UserExec(...), not to be called elsewhere!
+ //-------------------------------------------------------------------
+ Bool_t CheckInput();
+ void PostOutputs(); // cares about posting the output before exiting UserExec,
+ // WARNING: PostOutputs should only be used directly in the UserExec!!
+ void DoEmptyEventStudy(); // does the empty event study
+ Bool_t DetermineGap(); // determines the gap of all available detectors
+ void AnalyzeVtx(); // calcs the distance of the pri. vertex from tracks and SPD
+ void DoMultiplicityStudy(Int_t nMCprimaries);
+ // compares ITSTPC and ITSTPC+soft distribution and the multiplicity response
+ void DoTrackPair(Bool_t soft = kFALSE);
+ // analyses track pairs and prepares THnMother(Soft) input
+ Int_t DoMCTruth(); // analyses the MCtruth for corrections, returns #Primaries
+ void DoMCTrackPair(TParticle* particles[], Int_t gapCond, Int_t multiplicity);
+ // analyzes track pairs in MC truth
+ void DetermineMCprocessType(); // determines the MC process ID
+
+ // functions called by the DoTrackPair(...), not to be called elsewhere!
+ // ---------------------------------------------------------------------
+ Int_t DoPID(Int_t combCh); //retrieves the PID information for the current two
+ void FillMChists(Int_t combCh); //fills the MC histograms for two track events
+
+ // analysis task status
+ //---------------------
+ Bool_t fDoAOD; // true for running on AODs
+ Long_t fAnalysisStatus; // stores the status bits used to distinguish
+ // which histograms should be generated
+ TObjString *fAnalysisStatusString; // to be stored in the output list
+ Double_t fNoGapFraction; // fraction of no-gap events to be processed only
+ // used in case of ESD input
+ Double_t fReducedGapFraction; // fraction of events stored altough they
+ // show only the requested and not the refined
+ // gap
+ Double_t fMaxVtxDst; // maximum distance of the track and SPD vertex
+
+ // event information
+ //------------------
+ AliESDEvent *fESDEvent; // esd event object
+ AliAODEvent *fAODEvent; // esd event object
+ AliPIDResponse *fPIDResponse; // esd pid object
+ AliPhysicsSelection *fPhysicsSelection; // physics selection object
+ AliCDMesonTracks *fTracks; // object taking care about the track cuts
+ Double_t fVtxDst; // distance of the primary vertex from tracks and from SPD
+ Double_t fVtxZ; // z-position of the primary vertex from tracks
+ Int_t fResidualTracks; // tracks rejected by cuts within the event
+ Int_t fResidualTracklets; // SPD tracklets not assigned to tracks
+ Int_t fMCprocessType; // MC process type, 0 for data
+ Int_t fMCprocess; // detailed MC sub process information
+
+ // information about the trackpair which is currently processed
+ AliVTrack* fTrkPair[2]; // track objects
+
+ Int_t fRun; // number of the run which is about to be processed
+ Int_t fPIDmode; // selects set of PID cuts, 0 for 3sigma standard cuts,
+ //1 for LHC11f
+
+ // information stored for the PWA (addresses needed for the tree)
+ Float_t fTheta; // theta angle in the helicity frame
+ Float_t fPhi; // phi angle in the helicity frame
+ Float_t fMass; // mass of the two track system
+ Float_t fMomentum[3]; // momentum of the two track system
+ Int_t fCurrentGapCondition; // gap condition of the current event
+ Int_t fGapInformationWCent[kMax]; // gap condition for different detectors
+ // individually and their combinations
+ // requiring central act. from SPD FastOR
+ Int_t fGapInformation[kMax]; // same as above, without the requirement
+ // on central activity
+
+ // output objects
+ //---------------
+ THnSparseI *fGapRun; //contains the gap information of the different detectors
+ // per run (x-axis: run number, y-axis gap-state (bit encoded)
+
+ TList *fHist; // output list (contains all histograms)
+ THnSparseD *fThnMother; // THnSparse for mother pt and mass
+ THnSparseD *fThnMotherSoft; // same as above, using soft tracks as well
+ THnSparseD *fThnMultiplicity; // ThnSparse for multiplicity studies
+ THnSparseD *fThnMC; // THnSparse for MC mother pt and mass
+
+ THnSparseD *fThnEmptyEvents; // THnSparse used to analyse empty events
+ TTree *fPWAtree; // Tree containing information needed for PWA
+
+ // Multiplicity distributions for the different gap conditions
+ TH2D *fv0ntrk; //v0bit vs. nch
+ TH2D *fv0fmdntrk; //v0fmdbit vs. nch
+ TH2D *fv0fmdspdntrk; //v0fmdspdbit vs. nch
+ TH2D *fv0fmdspdtpcntrk; //v0fmdspdtpcbit vs. nch
+
+ // How many tracks are removed by the current cuts
+ TH2F *fv0Rmntrk; // V0 gap events vs. of removed tracks
+ TH2F *fv0fmdRmntrk; // V0-FMD gap events vs. of removed tracks
+ TH2F *fv0fmdspdRmntrk; // V0-FMD-SPD gap events vs. of removed tracks
+ TH2F *fv0fmdspdtpcRmntrk; // V0-FMD-SPD-TPC gap events vs. of removed tracks
+
+ TH2F *fMultStudy; // multiplicity after standard cuts and soft cuts vs. mult
+ TH2F *fMultStudyV0dg; // same for with events V0 Double-Gap
+ TH2F *fMultStudyV0FMDdg; // same for events with V0-FMD DG
+ TH2F *fMultStudyV0FMDSPDdg; // same for events with V0-FMD-SPD DG
+ TH2F *fMultStudyV0FMDSPDTPCdg; // same for events with V0-FMD-SPD-TPC DG
+
+ TH2F *fMultResponseMC; // multiplicity response MC to reconstruction
+ TH2F *fMultResponseV0dgMC; // same for events with V0 double gap
+ TH2F *fMultResponseV0FMDdgMC; // with V0-FMD DG
+ TH2F *fMultResponseV0FMDSPDdgMC; // with V0-FMD-SPD DG
+ TH2F *fMultResponseV0FMDSPDTPCdgMC; // with V0-FMD-SPD-TPC DG
+
+ TH2F *fMultRegionsMC; // multiplicity in the 3 regions of interest
+
+ TH1I *fhspdV0dg; // FastOR multiplicity in SPD FastOR with V0 double gap
+ TH1I *fhspdV0FMDdg; // with V0-FMD double gap
+ TH1I *fhspdV0FMDSPDdg; // with V0-FMD-SPD double gap
+ TH1I *fhspdV0FMDSPDTPCdg; // with V0-FMD-SPD double gap
+ TH1I *fhspdAfterCuts; // fast OR multiplicity after the event (vertex) cuts
+
+ TH1I *fGapResponseMCv0Dg; // gap status from MC and reconstruction for V0
+ TH1I *fGapResponseMCv0fmdDg; // for V0-FMD
+ TH1I *fGapResponseMCv0fmdspdDg; // for V0-FMD-SPD
+ TH1I *fGapResponseMCv0fmdspdtpcDg; // for V0-FMD-SPD-TPC
+
+ // Statistics flow diagrams
+ TH1I *fhspd; // distribution of fastOR-multiplicity
+ //SPD fastOR offline esdEvent->GetMultiplicity->TestFastOrFiredChips()
+ TH2I *fhfo; // histogram containing the SPD FastOR information ??
+ TH1I *fhpriVtxDist; // primary vertex distribution
+ TH1I *fhpriVtxPos; // position of the primary vertex (0 within range, 1 else)
+ TH1I *fhpriv; //primary vertex cut effect
+ TH1I *fhntrk; //n-trk-after-cut effect w/o requiring a gap-condition
+ TH1F *fhStatsFlow; // stepwise statistics flow
+ TH1F *fhFOchans; // hits in the different fastOR channels
+
+ TObjArray* fVZEROhists; // histograms for the VZERO trigger study
+
+ // eta-phi maps containing the ESDtracks
+ TH2F *fv0Map; // V0 gap events
+ TH2F *fv0fmdMap; // V0-FMD gap events
+ TH2F *fv0fmdspdMap; // V0-FMD-SPD gap events
+ TH2F *fv0fmdspdtpcMap; // V0-FMD-SPD-TPC gap events
+
+ // eta-phi maps containing the ESDtracks after the track cut
+ TH2F *fv0MapCutted; // V0 gap events
+ TH2F *fv0fmdMapCutted; // V0-FMD gap events
+ TH2F *fv0fmdspdMapCutted; // V0-FMD-SPD gap events
+ TH2F *fv0fmdspdtpcMapCutted; // V0-FMD-SPD-TPC gap events
+
+ TH2F *fHitMapSPDinner; // eta-phi map containing SPD FastOR hits (inner layer)
+ TH2F *fHitMapSPDouter; // eta-phi map containing SPD FastOR hits (outer layer)
+ TH2F *fHitMapSPDtrklt; // eta-phi map containing SPD hits based on tracklets
+ TH2F *fHitMapFMDa; // eta-phi-multiplicity map containing FMD hits at A side
+ TH2F *fHitMapFMDc; // eta-phi-multiplicity map containing FMD hits at C side
+
+ // Summed FMD
+ TH1F *fFMDsum1I; // FMD 1
+ TH1F *fFMDsum2I; // FMD 2 inner ring
+ TH1F *fFMDsum2O; // FMD 2 outer ring
+ TH1F *fFMDsum3I; // FMD 3 inner ring
+ TH1F *fFMDsum3O; // FMD 3 outer ring
+
+ // Vertex
+ TH1I *fPriVtxX; // x distribution of the primary vertex
+ TH1I *fPriVtxY; // y distribution of the primary vertex
+ TH1I *fPriVtxZ; // z distribution of the primary vertex
+ TH1I *fPriVtxDst; // distance distribution of the pri. vtx from SPD and tracks
+ TH1I *fPriVtxDstV0dg; // with V0 DG
+ TH1I *fPriVtxDstV0FMDdg; // with V0 FMD DG
+ TH1I *fPriVtxDstV0FMDSPDdg; // with V0 FMD SPD DG
+ TH1I *fPriVtxDstV0FMDSPDTPCdg; // with V0 FMD SPD TPC DG
+
+ // TPC Gap
+ TH2F *fTPCGapDCAaSide; // dca distribution for tracks destroying the TPC gap
+ TH2F *fTPCGapDCAcSide; // for V0-FMD-SPD gap events on A and C side
+
+ // PID
+ // All events but double gap events
+ TH2F *fComb2trkPIDuls; // PID matrix on unlike-sign two particles events
+ TH2F *fComb2trkPIDls; // PID matrix on like-sign about two particles
+ // Full-gap events (V0-FMD-SPD-TPC)
+ TH2F *fComb2trkPIDulsDG; // PID matrix on unlike-sign two particles events
+ TH2F *fComb2trkPIDlsDG; // PID matrix on like-sign about two particles
+ // 'PID' in MC
+ // All events but double gap events
+ TH2F *fComb2trkPIDulsMC; // PID matrix on unlike-sign two particles events
+ TH2F *fComb2trkPIDlsMC; // PID matrix on like-sign about two particles
+ // Full-gap events (V0-FMD-SPD-TPC)
+ TH2F *fComb2trkPIDulsDGmc; // PID matrix on unlike-sign two particles events
+ TH2F *fComb2trkPIDlsDGmc; // PID matrix on like-sign about two particles
+
+ // MC Process
+ TH2F *fMCProcessUls; // Process fulfilling a gap condition with two
+ // unlike-sign particles
+ TH2F *fMCProcessLs; // same for like-sign particles
+
+ TH1F *fMAllTrackMass; // histogram for the invariant mass dist. of all tracks
+
+ ClassDef(AliAnalysisTaskCDMeson, 1);
+};
+
+
+#endif
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+//
+// $Id$
+//
+// Task to skim ESD files, basic idea found at PWGGA/EMCALTasks as
+// AliEsdSkimTask. This version is modified to produce a exact copy of the ESD
+// stream containing only events which are of interest for the central
+// diffractive analysis (PWGUD) and a control sample. In order to be able to
+// process MC data as well, the usual concept using an output slot is for
+// output data is violated and an output file is directly created without
+// a corresponding output slot.
+//
+// Author:
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+// ROOT headers
+#include <TClonesArray.h>
+#include <TFile.h>
+#include <TTree.h>
+#include <TH1.h>
+#include <TAxis.h>
+#include <TRandom3.h>
+//#include <TKey.h>
+//#include <TROOT.h>
+
+// AliRoot headers
+#include "AliAnalysisManager.h"
+#include "AliCentrality.h"
+#include "AliESDEvent.h"
+#include "AliESDFMD.h"
+#include "AliEventplane.h"
+#include "AliInputEventHandler.h"
+#include "AliMultiplicity.h"
+#include "AliVEventHandler.h"
+#include "AliMCEventHandler.h"
+//#include "AliMCEvent.h"
+//#include "AliStack.h"
+//#include "AliGenEventHeader.h"
+//#include "AliRunLoader.h"
+
+// own classes
+#include "AliCDMesonBase.h"
+#include "AliCDMesonUtils.h"
+
+// own header
+#include "AliAnalysisTaskCDskimESD.h"
+
+ClassImp(AliAnalysisTaskCDskimESD)
+
+//______________________________________________________________________________
+AliAnalysisTaskCDskimESD::AliAnalysisTaskCDskimESD(const char *opt,
+ Bool_t reduceGapEvents)
+ : AliAnalysisTaskSE(opt)
+ , fInputTree(0x0)
+ , fInputTreeCopy(0x0)
+ , fOutputTree(0x0)
+ , fDoMC(kFALSE)
+ , fWrkDir(0x0)
+ , fKinematicsFile(0x0)
+ , fTrackRefsFile(0x0)
+ , fTEinput(0x0)
+ , fTE(0x0)
+ , fgaliceFile(0x0)
+ , fDoGapCond(kTRUE)
+ , fCurrentGapCondition(0)
+ , fRequestedGapCondition(AliCDMesonBase::kBitV0A + AliCDMesonBase::kBitV0C)
+ // only events containing a V0 double gap or
+ , fNoGapFraction(0.006) // 0.6% of the minimum bias events are stored
+ , fReduceGapEvents(reduceGapEvents)
+ , fRefinedGapCondition(AliCDMesonBase::kBitV0A + AliCDMesonBase::kBitV0C +
+ AliCDMesonBase::kBitFMDA + AliCDMesonBase::kBitFMDC)
+ , fReducedGapFraction(0.02) // 2% of the events with a loose gap are stored
+ , fStatsFlow(0x0)
+ , fSkimmingList(0x0)
+ , fFileName()
+ , fEventNumberInFile(-1)
+ , fRunNumber(-1)
+ , fEventTime(0)
+{
+ //
+ // Constructor.
+ //
+
+ if (!opt)
+ return;
+
+ //fBranchNames = "ESD:AliESDHeader.,AliESDRun."; // TODO don't we need all?
+
+ // check whether we are running on MC data
+ AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
+ AliMCEventHandler* mcH =
+ dynamic_cast<AliMCEventHandler*>(am->GetMCtruthEventHandler());
+ if(mcH) {
+ fDoMC = kTRUE;
+
+ fWrkDir = gDirectory; // save the old directory
+ fKinematicsFile = new TFile("Kinematics.root", "recreate");
+ fKinematicsFile->Close();
+ delete fKinematicsFile;
+ fKinematicsFile = 0x0;
+ gDirectory = fWrkDir; // restore the old directory
+
+ fTrackRefsFile = new TFile("TrackRefs.root", "recreate");
+ fTrackRefsFile->Close();
+ delete fTrackRefsFile;
+ fTrackRefsFile = 0x0;
+ gDirectory = fWrkDir; // restore the old directory
+ }
+
+ DefineOutput(1, TTree::Class());
+ DefineOutput(2, TH1::Class());
+ DefineOutput(3, TTree::Class());
+ if (fDoMC) DefineOutput(4, TTree::Class());
+}
+
+
+//______________________________________________________________________________
+AliAnalysisTaskCDskimESD::AliAnalysisTaskCDskimESD()
+ : AliAnalysisTaskSE()
+ , fInputTree(0x0)
+ , fInputTreeCopy(0x0)
+ , fOutputTree(0x0)
+ , fDoMC(kFALSE)
+ , fWrkDir(0x0)
+ , fKinematicsFile(0x0)
+ , fTrackRefsFile(0x0)
+ , fTEinput(0x0)
+ , fTE(0x0)
+ , fgaliceFile(0x0)
+ , fDoGapCond(kTRUE)
+ , fCurrentGapCondition(0)
+ , fRequestedGapCondition(AliCDMesonBase::kBitV0A + AliCDMesonBase::kBitV0C)
+ // only events containing a V0 double gap or
+ , fNoGapFraction(0.006) // 0.6% of the minimum bias events are stored
+ , fReduceGapEvents(kFALSE)
+ , fRefinedGapCondition(AliCDMesonBase::kBitV0A + AliCDMesonBase::kBitV0C +
+ AliCDMesonBase::kBitFMDA + AliCDMesonBase::kBitFMDC)
+ , fReducedGapFraction(0.02) // 2% of the events with a loose gap are stored
+ , fStatsFlow(0x0)
+ , fSkimmingList(0x0)
+ , fFileName()
+ , fEventNumberInFile(-1)
+ , fRunNumber(-1)
+ , fEventTime(0)
+{
+ //
+ // Default Constructor
+ //
+}
+
+
+//______________________________________________________________________________
+AliAnalysisTaskCDskimESD::~AliAnalysisTaskCDskimESD()
+{
+ //
+ // deconstructor
+ //
+
+ /* // this delete lead to chrases caused by the garbage collection ...
+ if (fOutputTree) {
+ delete fOutputTree;
+ fOutputTree = 0x0;
+ }
+ if (fInputTreeCopy) {
+ delete fInputTreeCopy;
+ fInputTreeCopy = 0x0;
+ }
+ if (fStatsFlow) {
+ delete fStatsFlow;
+ fStatsFlow = 0x0;
+ } */
+ if (fDoMC && fKinematicsFile) {
+ fKinematicsFile->Close();
+ delete fKinematicsFile;
+ fKinematicsFile = 0x0;
+ }
+ if (fDoMC && fTrackRefsFile) {
+ fTrackRefsFile->Close();
+ delete fTrackRefsFile;
+ fTrackRefsFile = 0x0;
+ }
+}
+
+
+//______________________________________________________________________________
+void AliAnalysisTaskCDskimESD::UserExec(Option_t */*opt*/)
+{
+ //
+ // Process event.
+ //
+
+ fStatsFlow->Fill(0.); // UserExec(...) <= bin name
+
+ if (fDoMC) {
+ // replace dummy tree by the correct tree
+ AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
+ if(!fTEinput) {
+ fTEinput =
+ ((AliMCEventHandler*)am->GetMCtruthEventHandler())->GetTree();
+ if (fTEinput) {
+ delete fTE;
+ fTE = fTEinput->CloneTree(0);
+ fTE->SetDirectory(fgaliceFile);
+ fTE->SetAutoFlush(-10*1024*1024);
+ //CopygALICE(am);
+ }
+ }
+ else if (fTEinput != am->GetMCtruthEventHandler()->GetTree()) {
+ fTEinput =
+ ((AliMCEventHandler*)am->GetMCtruthEventHandler())->GetTree();
+ fTEinput->CopyAddresses(fTE);
+ }
+ }
+
+ // check whether the ESD input tree is still the same - NOT NECESSARY
+ //TTree* currentInput =
+ // AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->GetTree();
+ //if (fInputTree != currentInput) { // ESD input tree has changed
+ // puts("CHANGED INPUT TREE");
+ // fInputTree = currentInput;
+ // fInputTreeCopy = (TTree*)fInputTree->Clone("esdTree");
+ // fInputTreeCopy->CopyAddresses(fOutputTree);
+ //}
+
+ AliESDEvent *esdin = dynamic_cast<AliESDEvent*>(InputEvent());
+ if (!esdin || !fInputTree || !fInputTreeCopy || !fOutputTree || !fTEinput
+ || !fTE) {
+ // check whether everything is ready
+ PostData(1, fOutputTree);
+ PostData(2, fStatsFlow);
+ PostData(3, fSkimmingList);
+ if (fDoMC) PostData(4, fTE);
+ return;
+ }
+
+ fStatsFlow->Fill(1.); // Ready <= bin name
+
+ //============================================================================
+ // event selection, every event passing the next lines is written to the new
+ // ESD stream
+ //
+ // modify this lines in order to adjust the ESD selection to your needs
+
+ // determine current gap condition
+ fCurrentGapCondition = AliCDMesonUtils::GetGapConfig(esdin, 0x0, 0x0, 0x0,
+ 0x0, 0x0, 0x0, 0x0, 0x0);
+
+ TRandom3 rnd(0);
+ if (!(AliCDMesonBase::kBitBaseLine & fCurrentGapCondition) ||
+ (fRequestedGapCondition & fCurrentGapCondition) || fReduceGapEvents) {
+ // check whether the event satifies the required gap condition and whether
+ // all double-gap events will be stored
+
+ if (fReduceGapEvents &&
+ (AliCDMesonBase::kBitBaseLine & fCurrentGapCondition) &&
+ !(fRequestedGapCondition & fCurrentGapCondition)) {
+ // not all double gap events can be stored, but the gap determination
+ // was successful and at least a according to the less stringent condition
+ // a double gap was found
+ if (!(fRefinedGapCondition & fCurrentGapCondition)) {
+ // event fulfilled the more stringent gap and hence is will be stored
+ fStatsFlow->Fill(5.);
+ }
+ else if (rnd.Rndm() > (1.-fReducedGapFraction)) {
+ // randomly select a fraction of
+ fStatsFlow->Fill(4.); // Double Gap Event <= bin name
+ }
+ else {
+ // reject event (no refined gap and not randomly selected)
+ fStatsFlow->Fill(3.); // Event Rejected <= bin name
+ PostData(1, fOutputTree);
+ PostData(2, fStatsFlow);
+ PostData(3, fSkimmingList);
+ if (fDoMC) PostData(4, fTE);
+ return;
+ }
+ }
+ else if (rnd.Rndm() > (1.-fNoGapFraction)) { // randomly selected
+ fStatsFlow->Fill(2.); // Control Sample Event <= bin name
+ }
+ else {
+ fStatsFlow->Fill(3.); // Event Rejected <= bin name
+ PostData(1, fOutputTree);
+ PostData(2, fStatsFlow);
+ PostData(3, fSkimmingList);
+ if (fDoMC) PostData(4, fTE);
+ return;
+ }
+ }
+ else {
+ fStatsFlow->Fill(4.); // Double Gap Event <= bin name
+ }
+ // end of event selection
+ //============================================================================
+
+ // load the current event in the cloned input tree
+ //
+ // unfortunately the Entry() gives the entry number with in the current input
+ // tree, not within the chain for MC data
+ // => use GetReadEntry() of the ESD tree (value return by Entry() for real
+ //data)
+ fInputTreeCopy->GetEvent(fInputTree->GetReadEntry());
+ fOutputTree->Fill(); // fill the current event into the new ESD stream
+ //printf("Entry()=%lld\nEntries=%lld\nChainOffset=%lld\nReadEntry=%lld\n",
+ // Entry(), fInputTreeCopy->GetEntries(),
+ // fInputTreeCopy->GetChainOffset(), fInputTree->GetReadEntry());
+
+ // MC specific stuff ---------------------------------------------------------
+ AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
+
+ if (fDoMC) { // do the MC related part of the skimming
+ Long64_t iEntry = fOutputTree->GetEntries() - 1;
+ AliMCEventHandler* mcHandler =
+ (AliMCEventHandler*)am->GetMCtruthEventHandler();
+ TTree* treeK = mcHandler->TreeK();
+ TTree* treeTR = mcHandler->TreeTR();
+ if (!treeK) AliFatal("TreeK not found!");
+ if (!treeTR) AliFatal("TreeTR not found!");
+
+ fWrkDir = gDirectory;
+ gDirectory = fKinematicsFile;
+ TString dirName = TString(Form("Event%ld", (Long_t)iEntry));
+ gDirectory->mkdir(dirName);
+ gDirectory->cd(dirName);
+ TTree* outputTreeK = treeK->CloneTree(); // copy the whole tree
+ outputTreeK->Write();
+ treeK->CopyAddresses(outputTreeK, kTRUE); // separate clone again
+ treeK->GetListOfClones()->Remove((TObject*)outputTreeK);
+ outputTreeK->Delete();
+ outputTreeK = 0x0;
+
+ gDirectory = fTrackRefsFile;
+ gDirectory->mkdir(dirName);
+ gDirectory->cd(dirName);
+ TTree* outputTreeTR = treeTR->CloneTree(); // copy the whole tree
+ outputTreeTR->Write();
+ treeTR->CopyAddresses(outputTreeTR, kTRUE); // separate clone again
+ treeTR->GetListOfClones()->Remove((TObject*)outputTreeTR);
+ outputTreeTR->Delete();
+ outputTreeTR = 0x0;
+ gDirectory = fWrkDir;
+
+ TTree* currentTEinput =
+ ((AliMCEventHandler*)am->GetMCtruthEventHandler())->GetTree();
+ if (fTEinput != currentTEinput) {
+ // when the input file is changed, the kinematics tree address changes
+ // hence our output tree fTE has to be linked to the new one!
+ fTEinput = currentTEinput;
+ fTEinput->CopyAddresses(fTE);
+ }
+ fTE->Fill(); // fill MC event headers
+ }
+
+ fStatsFlow->Fill(6.); // Tree Filled <= bin name
+
+ // remember which event we stored
+ fEventNumberInFile = esdin->GetEventNumberInFile();
+ fRunNumber = esdin->GetRunNumber();
+ fFileName = fInputTree->GetCurrentFile()->GetName();
+ fEventTime = esdin->GetTimeStamp();
+ fSkimmingList->Fill();
+
+ PostData(1, fOutputTree);
+ PostData(2, fStatsFlow);
+ PostData(3, fSkimmingList);
+ if (fDoMC) PostData(4, fTE);
+}
+
+//______________________________________________________________________________
+void AliAnalysisTaskCDskimESD::UserCreateOutputObjects()
+{
+ // Create output objects.
+ fStatsFlow = (TH1*)(new TH1F("skimmingStatsFlow", "", 7, 0., 7.));
+ TAxis* x = fStatsFlow->GetXaxis();
+ x->SetBinLabel(1, "UserExec(...)");
+ x->SetBinLabel(2, "Ready");
+ x->SetBinLabel(3, "Control Sample Event");
+ x->SetBinLabel(4, "Rejected Event");
+ x->SetBinLabel(5, "Double Gap Event");
+ x->SetBinLabel(6, "Refined Double Gap Event");
+ x->SetBinLabel(7, "Tree Filled");
+
+ PostData(2, fStatsFlow);
+
+ // TODO implement fSkimmingList
+ fSkimmingList = new TTree("SkimmingList", "SkimmingList");
+ fSkimmingList->Branch("FileName", &fFileName);
+ fSkimmingList->Branch("EventNumberInFile", &fEventNumberInFile);
+ fSkimmingList->Branch("RunNumber", &fRunNumber);
+ fSkimmingList->Branch("TimeStamp", &fEventTime);
+ PostData(3, fSkimmingList);
+
+ // Get input information
+ AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
+ fInputTree = am->GetInputEventHandler()->GetTree();
+ if (!fInputTree) {
+ puts("AliAnalysisTaskCDskimESD: input tree not found\n");
+ return;
+ }
+ fInputTreeCopy = (TTree*)fInputTree->Clone("esdTree");
+
+
+ // prevent the task from being run on proof
+ if (AliAnalysisManager::GetAnalysisManager()->GetAnalysisType() ==
+ AliAnalysisManager::kProofAnalysis) {
+ AliFatal("AliAnalysisTaskCDskimESD: cannot be run on PROOF!");
+ }
+
+ TFile *file = 0x0;
+ file = OpenFile(1); // open file for the first output slot
+ if (!file) {
+ puts("AliAnalysisTaskCDskimESD: file not opened\n");
+ return;
+ }
+
+ fOutputTree = fInputTreeCopy->CloneTree(0);
+ if (!fOutputTree) {
+ puts("AliAnalysisTaskCDskimESD: cloning tree not successful\n");
+ return;
+ }
+
+ file->SetCompressionLevel(5); // caused a crash
+ fOutputTree->SetDirectory(file);
+ fOutputTree->SetAutoFlush(-10*1024*1024);
+
+ if (fDoGapCond) {
+ if(!fOutputTree->GetBranch("gapCondition")) {
+ fOutputTree->Branch("gapCondition", &fCurrentGapCondition);
+ }
+ }
+
+ PostData(1, fOutputTree);
+
+ // process the mc information
+ if(fDoMC) {
+ fWrkDir = gDirectory; // save the old directory
+ fKinematicsFile = new TFile("Kinematics.root", "update");
+ gDirectory = fWrkDir; // restore the old directory
+
+ fTrackRefsFile = new TFile("TrackRefs.root", "update");
+ gDirectory = fWrkDir; // restore the old directory
+
+ fTEinput =
+ ((AliMCEventHandler*)am->GetMCtruthEventHandler())->GetTree();
+ if (!fTEinput) {
+ // this trick is done in order to post an output tree, although the tree
+ // which is cloned is not available during the UserCreateOutputObjects()
+ fTE = new TTree("TE", "TE");
+ }
+ else {
+ fTE = fTEinput->CloneTree(0);
+ //CopygALICE(am); // not needed, not properly working/tested
+ }
+ fgaliceFile = OpenFile(4);
+ if (!file) {
+ puts("AliAnalysisTaskCDskimESD: galice.root file not opened\n");
+ return;
+ }
+ fgaliceFile->SetCompressionLevel(5);
+ fTE->SetDirectory(fgaliceFile);
+ fTE->SetAutoFlush(-10*1024*1024);
+
+
+ PostData(4, fTE);
+ }
+}
+
+
+//______________________________________________________________________________
+/* // not needed => not yet working / properly tested
+void AliAnalysisTaskCDskimESD::CopygALICE(AliAnalysisManager* am)
+{
+ //
+ // copy all objects contained in the galice root excluding the TE tree
+ //
+
+ TString galiceInput =
+ *(((AliMCEventHandler*)am->GetMCtruthEventHandler())->GetInputPath());
+ galiceInput += "galice.root";
+ printf("galiceInput=%s\n", galiceInput.Data());
+ fWrkDir = gDirectory;
+ TFile* input = new TFile(galiceInput.Data(), "READ");
+
+ input->cd();
+ gDirectory->ls();
+
+ //loop on all entries of of the input directory
+ TKey *key;
+ TIter nextkey(input->GetListOfKeys());
+ input->GetListOfKeys()->ls();
+ while ((key = (TKey*)nextkey())) {
+ TString name = key->GetName();
+ printf("name=%s\n", name.Data());
+ const TString classname = key->GetClassName();
+ printf("classname=%s\n", classname.Data());
+ TClass *cl = gROOT->GetClass(classname);
+ if (!cl) continue;
+ if (!name.Contains("TE") && !classname.Contains("AliRun")) {
+ TObject *obj = key->ReadObj();
+ fgaliceFile->cd();
+ obj->Write();
+ delete obj;
+ input->cd();
+ }
+ else if (cl->InheritsFrom("AliRunLoader")) {
+ // TODO write the runloader stuff
+ //AliRunLoader* rl = AliRunLoader::Instance();
+ //printf("FILENAME=%s\n", rl->GetFileName().Data());
+ }
+ fgaliceFile->ls();
+ }
+ fgaliceFile->SaveSelf(kTRUE);
+ fWrkDir->cd();
+}
+*/
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+//
+// $Id$
+//
+// Task to skim ESD files, basic idea found at PWGGA/EMCALTasks as
+// AliEsdSkimTask. This version is modified to produce a exact copy of the ESD
+// stream containing only events which are of interest for the central
+// diffractive analysis (PWGUD) and a control sample. In order to be able to
+// process MC data as well, the usual concept using an output slot is for
+// output data is violated and an output file is directly created without
+// a corresponding output slot.
+//
+// Author:
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#ifndef ALIANALYSISTASKCDSKIMESD_H
+#define ALIANALYSISTASKCDSKIMESD_H
+
+#include "AliAnalysisTaskSE.h"
+
+// $Id$
+
+class TTree;
+class TH1;
+class TFile;
+class TDirectory;
+
+class AliAnalysisTaskCDskimESD : public AliAnalysisTaskSE {
+public:
+ AliAnalysisTaskCDskimESD(const char *opt, Bool_t reduceGapEvents);
+ AliAnalysisTaskCDskimESD();
+ virtual ~AliAnalysisTaskCDskimESD();
+
+ void UserExec(Option_t *opt);
+ void UserCreateOutputObjects();
+
+ void SetGapCondition(Int_t cond) { fRequestedGapCondition=cond; }
+ void SetRefinedGapCondition(Int_t cond) {fRefinedGapCondition=cond; }
+
+protected:
+ // variables concerning the reconstructed information in the ESDs
+ TTree* fInputTree; // input tree
+ TTree* fInputTreeCopy; // input tree copy
+ TTree* fOutputTree; // output tree
+
+ // variables needed for the reproduction of the MC information
+ Bool_t fDoMC; // do we process MC information?
+ TDirectory* fWrkDir; // directory the analysis task used before
+ // it was manipulated to store the
+ // kinemactics truth
+ TFile* fKinematicsFile; // output file containing the skimmed
+ // kinematics information
+ TFile* fTrackRefsFile; // output file containing the skimmed
+ // TrackRefs information
+ TTree* fTEinput; // input tree containing MC eventt headers
+ TTree* fTE; // tree containing the MC event headers
+ TFile* fgaliceFile; // pointer to the galice.root output file
+
+ // variables concerning the gaps used for the selection
+ Bool_t fDoGapCond; // if true the gap condition for central
+ // diffractive analysis are stored
+ Int_t fCurrentGapCondition; // gap condition
+ Int_t fRequestedGapCondition; // gap condition required for selection
+ Double_t fNoGapFraction; // fraction of no-gap events to be stored
+ Bool_t fReduceGapEvents; // do not store all V0 double-gap events
+ // necessary due to high double gap rates
+ // in MC productions (LHC10d{1,2,4,4a})
+ Int_t fRefinedGapCondition; // gap condition with more stringent
+ // selection applied in case
+ // fReduceDgEvents is active
+ Double_t fReducedGapFraction; // fraction of events stored altough they
+ // show only the requested and not the
+ // refined gap
+
+ // selection monitoring
+ TH1* fStatsFlow; // histogram indicating the number of
+ // input events
+ TTree* fSkimmingList; // which events of the intial ones are
+ // skimmed?
+ TString fFileName; // name of the input file
+ Int_t fEventNumberInFile; // event index in the ESD file
+ Int_t fRunNumber; // run number of the event
+ UInt_t fEventTime; // event timestamp
+
+private:
+ //void CopygALICE(AliAnalysisManager* am); // copy objects in galice.root
+ // not implemented:
+ AliAnalysisTaskCDskimESD(const AliAnalysisTaskCDskimESD&);
+ AliAnalysisTaskCDskimESD &operator=(const AliAnalysisTaskCDskimESD&);
+ ClassDef(AliAnalysisTaskCDskimESD, 1); // Esd skimming task
+};
+#endif // ALIANALYSISTASKCDSKIMESD_H
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+//
+// AliDDMesonBase
+// for
+// AliAnalysisTaskDDMeson
+//
+// Author:
+// Xianguo Lu <lu@physi.uni-heidelberg.de>
+// continued by
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#include "TH2.h"
+#include "TH2D.h"
+#include "TH1I.h"
+#include "TAxis.h"
+#include "THnSparse.h"
+#include "TString.h"
+#include "TList.h"
+
+#include "AliCDMesonBase.h"
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonBase::GetGapBin(TString tag, const Int_t gapcg,
+ Bool_t checkCentralActivity /* = kTRUE */)
+{
+ //
+ // retrieve gap topology for a given string
+ //
+
+ tag.ToUpper();
+
+ Bool_t ka = kFALSE, kc = kFALSE;
+
+ if(tag.Contains("V0")){
+ ka = ka || (gapcg & kBitV0A);
+ kc = kc || (gapcg & kBitV0C);
+ }
+ if(tag.Contains("FMD")){
+ ka = ka || (gapcg & kBitFMDA);
+ kc = kc || (gapcg & kBitFMDC);
+ }
+ if(tag.Contains("SPD")){
+ ka = ka || (gapcg & kBitSPDA);
+ kc = kc || (gapcg & kBitSPDC);
+ }
+ if(tag.Contains("TPC")){
+ ka = ka || (gapcg & kBitTPCA);
+ kc = kc || (gapcg & kBitTPCC);
+ }
+ if(tag.Contains("ZDC")){
+ ka = ka || (gapcg & kBitZDCA);
+ kc = kc || (gapcg & kBitZDCC);
+ }
+
+ if(ka && kc)
+ return kBinNG;
+ else{
+ if(!ka && !kc)
+ if (((gapcg & kBitCentAct) && checkCentralActivity) ||
+ !checkCentralActivity) {
+ return kBinDG; // central activity seen (or not required)
+ }
+ else {
+ return kBinNG; // no central activity
+ }
+ else if(!kc)
+ return kBinGC;
+ else
+ return kBinGA;
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonBase::CheckRange(Double_t &var, const Double_t min,
+ const Double_t max)
+{
+ //
+ // check whether the value is with in the range specified with min and max
+ //
+
+ const Double_t eps = 1e-3;
+ if( var >= max ) var = max - eps;
+ if( var <= min ) var = min + eps;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonBase::GetAxis(TString thntit, TString name)
+{
+ //
+ // get the number of an axis, derived from the ThnSparse title (list of axes)
+ //
+
+ thntit.ToUpper();
+ thntit.ReplaceAll(" ","");
+ const Int_t nmax = 20;
+ TString tits[nmax];
+ Int_t counter = 0;
+ while(thntit.Contains(",")){
+ const Int_t pos = thntit.First(",");
+ tits[counter] = thntit(0, pos);
+ thntit = thntit(pos+1, thntit.Length()-pos);
+ counter++;
+ if(counter>=nmax-1){
+ printf("freidtlog AliCDMesonBase::GetAxis too small nmax! %d %d\n",
+ counter, nmax);
+ return -1; //exit(1); // TODO
+ }
+ }
+ tits[counter++] = thntit;
+
+ //----------------
+
+ name.ToUpper();
+
+ for(Int_t ii=0; ii<counter; ii++){
+ if( tits[ii] == name )
+ return ii;
+ }
+ printf("freidtlog AliCDMesonBase::GetAxis !%s! %s not found!\n", name.Data(),
+ thntit.Data());
+ for(Int_t ii=0; ii<counter; ii++){
+ printf("*************** AliCDMesonBase::GetAxis *****************\n");
+ printf("freidtlog AliCDMesonBase::GetAxis %d !%s!\n", ii, tits[ii].Data());
+ printf("\n");
+ }
+ return -1; //exit(1); // TODO
+}
+
+
+//------------------------------------------------------------------------------
+TString AliCDMesonBase::GetTitleMother()
+{
+ //
+ // get the title of ThnMother (= list of axes)
+ //
+
+ TString title = "Ncombined, CombCh, CombPID, V0, FMD, SPD, TPC,";
+ title += " Mass, Pt, CTS, OA, DaughterPt, TrackResiduals, VertexZinRng,";
+ title += " ProcessType, VertexCoincidence, TrackletResiduals";
+ return title;
+}
+
+
+//------------------------------------------------------------------------------
+THnSparseD * AliCDMesonBase::GetThnMother(TString name /* = "CDMeson_Mother" */)
+{
+ // creates the THnSparse called Mother (used for the generation of invariant
+ // mass distributions)
+ //
+
+ const Int_t nbin[] = {
+ fgkNNcombined, fgkNCombCh, fgkNCombPID, fgkNGapConfig, fgkNGapConfig,
+ fgkNGapConfig, fgkNGapConfig, fgkNMass, fgkNMotherPt, fgkNCTS, fgkNOA,
+ fgkNDaughterPt, fgkNTrackResiduals, fgkNVertexZinRng, fgkNProcessType,
+ fgkNVertexCoincidence, fgkNTrackletResiduals
+ };
+ const Double_t binmin[] = {
+ fgkMinNcombined, fgkMinCombCh, fgkMinCombPID, fgkMinGapConfig,
+ fgkMinGapConfig, fgkMinGapConfig, fgkMinGapConfig, fgkMinMass,
+ fgkMinMotherPt, fgkMinCTS, fgkMinOA, fgkMinDaughterPt,
+ fgkMinTrackResiduals, fgkMinVertexZinRng, fgkMinProcessType,
+ fgkMinVertexCoincidence, fgkMinTrackletResiduals
+ };
+ const Double_t binmax[] = {
+ fgkMaxNcombined, fgkMaxCombCh, fgkMaxCombPID, fgkMaxGapConfig,
+ fgkMaxGapConfig, fgkMaxGapConfig, fgkMaxGapConfig, fgkMaxMass,
+ fgkMaxMotherPt, fgkMaxCTS, fgkMaxOA, fgkMaxDaughterPt,
+ fgkMaxTrackResiduals, fgkMaxVertexZinRng, fgkMaxProcessType,
+ fgkMaxVertexCoincidence, fgkMaxTrackletResiduals
+ };
+
+ const Int_t npar = sizeof(nbin)/sizeof(Int_t);
+
+ return new THnSparseD(name.Data(), GetTitleMother(), npar, nbin, binmin,
+ binmax);
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonBase::FillThnMother(THnSparseD *thn, const Double_t vNch,
+ const Double_t vCombCh,
+ const Double_t vCombPID, const Double_t vV0,
+ const Double_t vFMD, const Double_t vSPD,
+ const Double_t vTPC, const Double_t vMass,
+ const Double_t vPt, const Double_t vOA,
+ const Double_t vCTS,
+ const Double_t vDaughterPt,
+ const Double_t vTrackResiduals,
+ const Double_t vVertexZ,
+ const Double_t vProcessType,
+ const Double_t vVertexCoincidence,
+ const Double_t vTrkltResiduals)
+{
+ //
+ // fill ThnMother
+ //
+
+ Double_t var[]={
+ vNch, vCombCh, vCombPID, vV0, vFMD, vSPD, vTPC, vMass, vPt, vOA, vCTS,
+ vDaughterPt, vTrackResiduals, vVertexZ, vProcessType, vVertexCoincidence,
+ vTrkltResiduals
+ };
+ const Int_t nv = sizeof(var)/sizeof(Double_t);
+ if(nv!=thn->GetNdimensions()){
+ printf("AliCDMesonBase::FillThnMother nv error!! %d %d\n", nv,
+ thn->GetNdimensions());
+ return; //exit(1);
+ }
+
+ CheckRange(var[7], fgkMinMass, fgkMaxMass);
+ CheckRange(var[8], fgkMinMotherPt, fgkMaxMotherPt);
+ CheckRange(var[11], fgkMinDaughterPt, fgkMaxDaughterPt);
+
+ thn->Fill(var);
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonBase::GetAxisMother(const TString name)
+{
+ //
+ // return axis number corresponding to the name
+ //
+
+ return GetAxis(GetTitleMother(), name);
+}
+
+
+//------------------------------------------------------------------------------
+TString AliCDMesonBase::GetTitleEmptyEvents()
+{
+ //
+ // title / list of axes of the empty event study THnSparse
+ //
+
+ TString title = "EventType, FMD-A, FMD-C, SPD-I-A, SPD-I-C, SPD-O-A, SPD-O-C";
+ title += ", SPDtrkltA, SPDtrkltC";
+ title += ", fmdSum1I, fmdSum2I, fmdSum2O, fmdSum3I, fmdSum3O";
+ //title += ", TPC_all, TPC_diffVertex";
+ return title;
+}
+
+
+//------------------------------------------------------------------------------
+THnSparseD* AliCDMesonBase::GetThnEmptyEvents()
+{
+ // creates the THnSparse called Empty Events
+ // EventType, FMD-A, FMD-C, SPD-I-A, SPD-I-C, SPD-O-A, SPD-O-C, FMD1I..FMD3O
+ // TPC_all, TPC_diffVertex (not used so far)
+
+ const Int_t nbin[] = {
+ fgkNEventType, fgkNMultW, fgkNMultW, fgkNMult, fgkNMult, fgkNMult, fgkNMult,
+ fgkNMult, fgkNMult, fgkNMult, fgkNMult, fgkNMult, fgkNMult, fgkNMult
+ //, fgkNMult, fgkNMult
+ };
+ const Double_t binmin[] = {
+ fgkMinEventType, fgkMinMultW, fgkMinMultW, fgkMinMult, fgkMinMult,
+ fgkMinMult, fgkMinMult, fgkMinMult, fgkMinMult, fgkMinMult, fgkMinMult,
+ fgkMinMult, fgkMinMult, fgkMinMult //, fgkMinMult, fgkMinMult
+ };
+ const Double_t binmax[] = {
+ fgkMaxEventType, fgkMaxMultW, fgkMaxMultW, fgkMaxMult, fgkMaxMult,
+ fgkMaxMult, fgkMaxMult, fgkMaxMult, fgkMaxMult, fgkMaxMult, fgkMaxMult,
+ fgkMaxMult, fgkMaxMult, fgkMaxMult //, fgkMaxMult, fgkMaxMult
+ };
+
+ const Int_t npar = sizeof(nbin)/sizeof(Int_t);
+
+ return new THnSparseD("CDMeson_EmptyEvents", GetTitleEmptyEvents(), npar,
+ nbin, binmin, binmax);
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonBase::FillThnEmptyEvents(THnSparseD * thn, const Int_t eventType,
+ const Int_t multFMDA,
+ const Int_t multFMDC,
+ const Int_t multSPDIA,
+ const Int_t multSPDIC,
+ const Int_t multSPDOA,
+ const Int_t multSPDOC,
+ const Int_t multSPDtrkltA,
+ const Int_t multSPDtrkltC,
+ const Int_t fmdSum1I,
+ const Int_t fmdSum2I,
+ const Int_t fmdSum2O,
+ const Int_t fmdSum3I,
+ const Int_t fmdSum3O/*,
+ const Int_t multTPC,
+ const Int_t multTPCdiffVertex */)
+{
+ //
+ // Fill ThnEmptyEvents
+ //
+
+ Double_t var[]={
+ eventType, multFMDA, multFMDC, multSPDIA, multSPDIC, multSPDOA,
+ multSPDOC, multSPDtrkltA, multSPDtrkltC, fmdSum1I, fmdSum2I, fmdSum2O,
+ fmdSum3I, fmdSum3O//, multTPC, multTPCdiffVertex
+ };
+ const Int_t nv = sizeof(var)/sizeof(Double_t);
+ if(nv!=thn->GetNdimensions()){
+ printf("AliCDMesonBase::FillThnEmptyEvents nv error!! %d %d\n", nv,
+ thn->GetNdimensions());
+ return; //exit(1); // TODO
+ }
+
+ thn->Fill(var);
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonBase::GetAxisEmptyEvents(const TString name)
+{
+ //
+ // return axis number corresponding to the name
+ //
+
+ return GetAxis(GetTitleEmptyEvents(), name);
+}
+
+
+//------------------------------------------------------------------------------
+TString AliCDMesonBase::GetTitleMultiplicity()
+{
+ //
+ // get title of the multiplicity study ThnSparse
+ //
+
+ TString title = "Nch, Nsoft, Ncombined, V0, FMD, SPD, TPC, NresidualTracks";
+ title += ", NresidualTracklets, VertexZ, VerticesDistance, ProcessType";
+ return title;
+}
+
+
+//------------------------------------------------------------------------------
+THnSparseD* AliCDMesonBase::GetThnMultiplicity()
+{
+ //
+ // creates the THnSparse called Multiplicity
+ //
+
+ const Int_t nbin[] = {
+ fgkNNch, fgkNNsoft, fgkNNcomb,fgkNGapConfig, fgkNGapConfig,
+ fgkNGapConfig, fgkNGapConfig, fgkNNresidualTracks, fgkNNresidualTracklets,
+ fgkNVertexZ, fgkNVerticesDistance, fgkNProcessType
+ };
+ const Double_t binmin[] = {
+ fgkMinNch, fgkMinNsoft, fgkMinNcomb, fgkMinGapConfig, fgkMinGapConfig,
+ fgkMinGapConfig, fgkMinGapConfig, fgkMinNresidualTracks,
+ fgkMinNresidualTracklets, fgkMinVertexZ, fgkMinVerticesDistance,
+ fgkMinProcessType
+ };
+ const Double_t binmax[] = {
+ fgkMaxNch, fgkMaxNsoft, fgkMaxNcomb, fgkMaxGapConfig, fgkMaxGapConfig,
+ fgkMaxGapConfig, fgkMaxGapConfig, fgkMaxNresidualTracks,
+ fgkMaxNresidualTracklets, fgkMaxVertexZ, fgkMaxVerticesDistance,
+ fgkMaxProcessType
+ };
+
+ const Int_t npar = sizeof(nbin)/sizeof(Int_t);
+
+ return new THnSparseD("CDMeson_Multiplicity", GetTitleMultiplicity(), npar,
+ nbin, binmin, binmax);
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonBase::FillThnMultiplicity(THnSparseD *thn, const Double_t vNch,
+ const Double_t vNsoft,
+ const Double_t vNcombined,
+ const Double_t vV0,
+ const Double_t vFMD,
+ const Double_t vSPD,
+ const Double_t vTPC,
+ const Double_t vNresidualTracks,
+ const Double_t vNresidualTracklets,
+ const Double_t vVertexZ,
+ const Double_t vVerticesDistance,
+ const Double_t vProcessType)
+{
+ // fill ThnMultiplicity
+ // input list copied from GetTitle
+ // var[] copied from GetTitle
+
+ Double_t var[]={
+ vNch, vNsoft, vNcombined, vV0, vFMD, vSPD, vTPC, vNresidualTracks,
+ vNresidualTracklets, vVertexZ, vVerticesDistance, vProcessType
+ };
+ const Int_t nv = sizeof(var)/sizeof(Double_t);
+ if(nv!=thn->GetNdimensions()){
+ printf("AliCDMesonBase::FillThnMultiplicity nv error!! %d %d\n", nv,
+ thn->GetNdimensions());
+ return; //exit(1);
+ }
+
+ thn->Fill(var);
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonBase::GetAxisMultiplicity(const TString name)
+{
+ //
+ // return axis number corresponding to the name
+ //
+
+ return GetAxis(GetTitleMultiplicity(), name);
+}
+
+
+//------------------------------------------------------------------------------
+TH1F* AliCDMesonBase::GetHistStatsFlow()
+{
+ //
+ // setup the stats flow histogram
+ //
+
+ TH1F *hist = new TH1F("c00_statsFlow", "", AliCDMesonBase::kBinLastValue,
+ 0, AliCDMesonBase::kBinLastValue);
+ TAxis* axis = hist->GetXaxis();
+ axis->SetBinLabel(AliCDMesonBase::kBinTotalInput+1, "total Input");
+ axis->SetBinLabel(AliCDMesonBase::kBinGoodInput+1, "good ESDs");
+ axis->SetBinLabel(AliCDMesonBase::kBinV0OR+1, "V0-OR");
+ axis->SetBinLabel(AliCDMesonBase::kBinV0AND+1, "V0-AND");
+ axis->SetBinLabel(AliCDMesonBase::kBinEventsAfterCuts+1, "after cuts");
+ axis->SetBinLabel(AliCDMesonBase::kBinEventsWithOutPileUp+1, "w/o pile up");
+ axis->SetBinLabel(AliCDMesonBase::kBinv0Gap+1, "with V0 DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinv0fmdGap+1, "with V0-FMD DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinv0fmdspdGap+1,
+ "with V0-FMD-SPD DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinv0fmdspdtpcGap+1,
+ "with V0-FMD-SPD-TPC DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinv0fmdspdtpczdcGap+1,
+ "with V0-FMD-SPD-TPC-ZDC DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinfmdGap+1, "with FMD DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinspdGap+1, "with SPD DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBintpcGap+1, "with TPC DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBintpcspdGap+1, "width TPC-SPD DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBintpcspdfmdGap+1,
+ "width TPC-SPD-FMD DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBintpcspdfmdv0Gap+1,
+ "width TPC-SPD-FMD-V0 DG gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinspdfmdGap+1, "with SPD FMD gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinspdfmdv0Gap+1, "with SPD FMD V0 gap");
+ axis->SetBinLabel(AliCDMesonBase::kBinTwoTrackEvents+1, "with two tracks");
+ axis->SetBinLabel(AliCDMesonBase::kBinThreeTrackEvents+1,
+ "with three tracks");
+ axis->SetBinLabel(AliCDMesonBase::kBinPionEvents+1, "with two pions");
+ axis->SetBinLabel(AliCDMesonBase::kBinKaonEvents+1, "with two kaons");
+ axis->SetBinLabel(AliCDMesonBase::kBinProtonEvents+1, "with two proton");
+ axis->SetBinLabel(AliCDMesonBase::kBinElectronEvents+1, "with two electron");
+ axis->SetBinLabel(AliCDMesonBase::kBinUnknownPIDEvents+1, "with unknown PID");
+ axis->SetBinLabel(AliCDMesonBase::kBinResidualTracks+1,
+ "without residual tracks");
+ axis->SetBinLabel(AliCDMesonBase::kBinResidualTracklets+1,
+ "without residual tracklets");
+ axis->SetBinLabel(AliCDMesonBase::kBinCDonlyEvents+1, "CD only events");
+ return hist;
+}
+
+
+//------------------------------------------------------------------------------
+TH2F* AliCDMesonBase::GetHistPIDStudies(TString name)
+{
+ //
+ // setup the PID studies histogram
+ //
+
+ TH2F *hist = new TH2F(name.Data(), ";particle 1;particle 2",
+ AliCDMesonBase::kBinPIDUnknown,
+ AliCDMesonBase::kBinPionE,
+ AliCDMesonBase::kBinPIDUnknown+1,
+ AliCDMesonBase::kBinPIDUnknown,
+ AliCDMesonBase::kBinPionE,
+ AliCDMesonBase::kBinPIDUnknown+1);
+ TAxis* x = hist->GetXaxis();
+ TAxis* y = hist->GetYaxis();
+ x->SetBinLabel(AliCDMesonBase::kBinPionE, "#pi (ex)");
+ x->SetBinLabel(AliCDMesonBase::kBinPion, "#pi");
+ x->SetBinLabel(AliCDMesonBase::kBinSinglePion, "-");
+ x->SetBinLabel(AliCDMesonBase::kBinKaonE, "K (ex)");
+ x->SetBinLabel(AliCDMesonBase::kBinKaon, "K");
+ x->SetBinLabel(AliCDMesonBase::kBinSingleKaon, ",");
+ x->SetBinLabel(AliCDMesonBase::kBinProtonE, "p (ex)");
+ x->SetBinLabel(AliCDMesonBase::kBinProton, "p");
+ x->SetBinLabel(AliCDMesonBase::kBinSingleProton, "_");
+ x->SetBinLabel(AliCDMesonBase::kBinElectronE, "e (ex)");
+ x->SetBinLabel(AliCDMesonBase::kBinElectron, "e");
+ x->SetBinLabel(AliCDMesonBase::kBinSingleElectron, ".");
+ x->SetBinLabel(AliCDMesonBase::kBinPIDUnknown, "X");
+ y->SetBinLabel(AliCDMesonBase::kBinPionE, "#pi (ex)");
+ y->SetBinLabel(AliCDMesonBase::kBinPion, "#pi");
+ y->SetBinLabel(AliCDMesonBase::kBinSinglePion, "-");
+ y->SetBinLabel(AliCDMesonBase::kBinKaonE, "K (ex)");
+ y->SetBinLabel(AliCDMesonBase::kBinKaon, "K");
+ y->SetBinLabel(AliCDMesonBase::kBinSingleKaon, ",");
+ y->SetBinLabel(AliCDMesonBase::kBinProtonE, "p (ex)");
+ y->SetBinLabel(AliCDMesonBase::kBinProton, "p");
+ y->SetBinLabel(AliCDMesonBase::kBinSingleProton, "_");
+ y->SetBinLabel(AliCDMesonBase::kBinElectronE, "e (ex)");
+ y->SetBinLabel(AliCDMesonBase::kBinElectron, "e");
+ y->SetBinLabel(AliCDMesonBase::kBinSingleElectron, ".");
+ y->SetBinLabel(AliCDMesonBase::kBinPIDUnknown, "X");
+ return hist;
+}
+
+
+//------------------------------------------------------------------------------
+TObjArray* AliCDMesonBase::GetHistVZEROStudies(TList* l)
+{
+ //
+ // Create Histograms for the VZERO trigger studies
+ //
+
+ TObjArray* arr = new TObjArray(130);
+
+ TObject* o = 0x0;
+ for (Int_t iPMT = 0; iPMT < 64; ++iPMT) {
+ o = (TObject*)new TH2F(Form("h00_%02d_ADC_TriggerThr", iPMT),
+ ";ADC Counts;Trigger Threshold (ADC Counts)",
+ 400, 0., 4000., 48, 2., 50.);
+ arr->Add(o);
+ l->Add(o);
+ }
+ for (Int_t iPMT = 0; iPMT < 64; ++iPMT) {
+ o = (TObject*)new TH2F(Form("h01_%02d_ADC_Multiplicity", iPMT),
+ ";ADC Counts;Multiplicity",
+ 400, 0., 4000., 250, 0., 250.);
+ arr->Add(o);
+ l->Add(o);
+ }
+ /*
+ // not of use for pp data in 2010
+ o = (TObject*)new TH2F("h02_TriggerChargeA_Trigger",
+ ";Trigger Charge A;Trigger Decision A",
+ 250, 0., 5000., 2, 0., 2.);
+ arr->Add(o);
+ l->Add(o);
+ o = (TObject*)new TH2F("h03__TriggerChargeC_Trigger",
+ ";Trigger Charge C;Trigger Decision C",
+ 250, 0., 5000., 2, 0., 2.);
+ arr->Add(o);
+ l->Add(o);
+ */
+
+ return arr;
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonBase::GetGapTriggers(THnSparseI* gaprun, Int_t gapCondition,
+ const Int_t run, Double_t& triggers,
+ Double_t& total)
+{
+ // determine the number of certain triggers with a gap in the detectors
+ // specified by gapCondition in run and the total number of events
+ // surviving all cuts (including pile-up rejection)
+
+ triggers = 0;
+ total = 0;
+ Int_t nTuple[] = { gaprun->GetAxis(0)->FindBin(run), 0 };
+ for (Int_t i = 0; i< kBitGapMax; ++i) {
+ nTuple[1] = i;
+ Double_t temp = gaprun->GetBinContent(nTuple);
+ if (!(i & gapCondition)) {
+ triggers += temp;
+ }
+ total += temp;
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonBase::GetNoGapTriggers(THnSparseI* gaprun, Int_t gapCondition,
+ const Int_t run, Double_t& triggers,
+ Double_t& total)
+{
+ // determine the number of certain triggers with a NO gap in the detectors
+ // specified by the gapCondition in run and the total number of events
+ // surviving all cuts (including pile-up rejection)
+
+ triggers = 0;
+ total = 0;
+ Int_t nTuple[] = { gaprun->GetAxis(0)->FindBin(run), 0 };
+ for (Int_t i = 0; i< kBitGapMax; ++i) {
+ nTuple[1] = i;
+ Double_t temp = gaprun->GetBinContent(nTuple);
+ if (i & gapCondition) {
+ triggers += temp;
+ }
+ total += temp;
+ }
+}
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+//
+// AliCDMesonBase
+// for
+// AliAnalysisTaskCDMeson
+//
+//
+// Author:
+// Xianguo Lu <lu@physi.uni-heidelberg.de>
+// continued by
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#ifndef ALICDMESONBASE_H
+#define ALICDMESONBASE_H
+
+#include "THnSparse.h"
+
+class TH2;
+class TH2F;
+class TH1F;
+class TList;
+class TObjArray;
+
+class AliCDMesonBase
+{
+public:
+ enum{
+ // combined charged
+ kBinPM = 1, // unlike-sign
+ kBinPPMM, // like-sign
+
+ // pid results (used for single particles and pairs)
+ kBinPionE = 1, // pions with exclusion cut
+ kBinPion, // pions without exclusion cut but TPC and TOF
+ kBinSinglePion, // one pion identified, the other particle unidentified
+ kBinKaonE, // same for kaons ...
+ kBinKaon,
+ kBinSingleKaon,
+ kBinProtonE, // protons
+ kBinProton,
+ kBinSingleProton,
+ kBinElectronE, // and electrons
+ kBinElectron,
+ kBinSingleElectron,
+ kBinPIDUnknown, // could not be identified
+
+ // gap Conditions
+ kBinDG = 1, // double gap
+ kBinGC, // single gap c side
+ kBinGA, // single gap a side
+ kBinNG, // no gap
+
+
+ // event types used for the empty event studies
+ kBinEventI = 1, // Beam-Beam Interaction
+ kBinEventA, // Beam from A-side
+ kBinEventC, // Beam from C-side
+ kBinEventAC, // Beam from one side (there is no separation in 2011)
+ kBinEventE, // no beam (empty event)
+ kBinEventUnknown,
+
+ // StatsFlow histogram entries
+ // names for the bins are specified in the .cxx-file
+ kBinTotalInput = 0,
+ kBinGoodInput,
+ kBinV0OR,
+ kBinV0AND,
+ kBinEventsAfterCuts,
+ kBinEventsWithOutPileUp,
+ kBinv0Gap,
+ kBinv0fmdGap,
+ kBinv0fmdspdGap,
+ kBinv0fmdspdtpcGap,
+ kBinv0fmdspdtpczdcGap,
+ kBinfmdGap,
+ kBinspdGap,
+ kBintpcGap,
+ kBintpcspdGap,
+ kBintpcspdfmdGap,
+ kBintpcspdfmdv0Gap,
+ kBinspdfmdGap,
+ kBinspdfmdv0Gap,
+ kBinTwoTrackEvents,
+ kBinThreeTrackEvents,
+ kBinPionEvents,
+ kBinKaonEvents,
+ kBinProtonEvents,
+ kBinElectronEvents,
+ kBinUnknownPIDEvents,
+ kBinResidualTracks,
+ kBinResidualTracklets,
+ kBinCDonlyEvents,
+ kBinLastValue, // used to specify the correct histogram width
+
+
+ // MC event/process types
+ kBinND = 0, // used for data as well, "Non-Diffractive"
+ kBinCD, // central diffractive
+ kBinSD, // single diffractive
+ kBinDD, // double diffractive
+
+
+ // gap-condition bits used in AliAnalysisTaskCDMeson::fGapRun
+ kBitBaseLine = (1<<0),
+
+ kBitV0A = (1<<1),
+ kBitV0C = (1<<2),
+ kBitFMDA = (1<<3),
+ kBitFMDC = (1<<4),
+
+ kBitSPDA = (1<<5),
+ kBitSPDC = (1<<6),
+ kBitTPCA = (1<<7),
+ kBitTPCC = (1<<8),
+
+ kBitZDCA = (1<<9),
+ kBitZDCC = (1<<10),
+ kBitCentAct = (1<<11),
+ kBitGapMax = (1<<12),
+
+
+ // analysis task status bits
+ // do not change the order in order to be backward compatible!
+ kBitConfigurationSet = (1<<0), // if not set everything is active
+ kBitImposeZDCGap = (1<<1), // excluding the ZDC gap
+ kBitEtaPhiMaps = (1<<2),
+ kBitEtaPhiMapsWithCuts = (1<<3),
+ kBitStatsFlow = (1<<4),
+ kBitMultPerGapHists = (1<<5),
+ kBitRmMultPerGapHists = (1<<6),
+ kBitTHnMother = (1<<7),
+ kBitFastORStudy = (1<<8), // also not active by default
+ kBitHitMapSPD = (1<<9),
+ kBitHitMapFMD = (1<<10),
+ kBitVtxStudies = (1<<11),
+ kBitEEStudy = (1<<12),
+ kBitPIDStudy = (1<<13),
+ kBitMCProcess = (1<<14),
+ kBitFMDsum = (1<<15),
+ kBitSoftTracks = (1<<16),
+ kBitPWAtree = (1<<17),
+ kBitTHnMC = (1<<18),
+ kBitMultResponseMC = (1<<19),
+ kBitMultStudy = (1<<20),
+ kBitReadPreprocessedGap = (1<<21), // not active by default as well
+ kBitReduceGapEvents = (1<<22),
+ kBitVZEROStudy = (1<<23),
+ kBitTPCGapStudy = (1<<24),
+ kBitAllTrackMass = (1<<25),
+ kBitCDonly = (1<<26),
+ kBitFastORmultStudy = (1<<27),
+ kBitConfigurationVersion = (1<<28) // always set, last bit
+ };
+
+ static Int_t GetGapBin(const TString tag, const Int_t gapcg,
+ Bool_t checkCentralActivity = kTRUE);
+
+ static THnSparseD* GetThnMother(TString name = "CDMeson_Mother");
+ static void FillThnMother(THnSparseD *thn, const Double_t vNch,
+ const Double_t vCombCh, const Double_t vCombPID,
+ const Double_t vV0, const Double_t vFMD,
+ const Double_t vSPD, const Double_t vTPC,
+ const Double_t vMass, const Double_t vPt,
+ const Double_t vOA, const Double_t vCTS,
+ const Double_t vDaughterPt,
+ const Double_t vTrackResiduals,
+ const Double_t vVertexZ,
+ const Double_t vProcessType,
+ const Double_t vVertexCoincidence,
+ const Double_t vTrkltResiduals);
+ static Int_t GetAxisMother(const TString name);
+
+ static THnSparseD* GetThnEmptyEvents();
+ static void FillThnEmptyEvents(THnSparseD *thn, const Int_t eventType,
+ const Int_t multFMDA, const Int_t multFMDC,
+ const Int_t multSPDIA, const Int_t multSPDIC,
+ const Int_t multSPDOA, const Int_t multSPDOC,
+ const Int_t multSPDtrklA,
+ const Int_t multSPDtrklC, const Int_t fmdSum1I,
+ const Int_t fmdSum2I, const Int_t fmdSum20,
+ const Int_t fmdSum3I, const Int_t fmdSum3O/*,
+ const Int_t multTPC,
+ const Int_t multTPCdiffVertex */);
+ static Int_t GetAxisEmptyEvents(const TString name);
+
+ static THnSparseD* GetThnMultiplicity();
+ static void FillThnMultiplicity(THnSparseD *thn, const Double_t vNch,
+ const Double_t vNsoft,
+ const Double_t vNcombined,
+ const Double_t vV0, const Double_t vFMD,
+ const Double_t vSPD, const Double_t vTPC,
+ const Double_t vNresidualTracks,
+ const Double_t vNresidualTracklets,
+ const Double_t vVertexZ,
+ const Double_t vVerticesDistance,
+ const Double_t vProcessType);
+ static Int_t GetAxisMultiplicity(const TString name);
+
+ static TH1F* GetHistStatsFlow();
+ static TH2F* GetHistPIDStudies(TString name);
+ static TObjArray* GetHistVZEROStudies(TList* l);
+
+ static void GetGapTriggers(THnSparseI* gaprun, Int_t gapCondition,
+ const Int_t run, Double_t& triggers,
+ Double_t& total);
+
+ static void GetNoGapTriggers(THnSparseI* gaprun, Int_t gapCondition,
+ const Int_t run, Double_t& triggers,
+ Double_t& total);
+private:
+ static void CheckRange(Double_t &var, const Double_t min, const Double_t max);
+ static Int_t GetAxis(TString thntit, TString name);
+ static TString GetTitleMother();
+ static TString GetTitleEmptyEvents();
+ static TString GetTitleMultiplicity();
+
+ //== INVARIANT MASS DISTRIBUTIONS (ThnMother) ================================
+ //-- event characteristics
+ // number of charged primary particles - combined => can be w/ soft tracks
+ static const Int_t fgkNNcombined = 2; // number of bins
+ static const Double_t fgkMinNcombined = 2; // lower border
+ static const Double_t fgkMaxNcombined = 4; // upper border (#bins + lower border)
+
+ // unlike sign or like sign charged?
+ static const Int_t fgkNCombCh = 2; // kBinPPMM = number of bins
+ static const Double_t fgkMinCombCh = 1.; // kBinPM = lower border
+ static const Double_t fgkMaxCombCh = 3.; // kBinPPMM + kBinPM = upper border
+
+ // two track PID, take care the enum changed from time to time
+ static const Int_t fgkNCombPID = 13; // kBinPIDUnknown = number of bins
+ static const Double_t fgkMinCombPID = 1.; // kBinPionE = lower border
+ static const Double_t fgkMaxCombPID = 14.;// ... = upper border
+
+ // gap configuration, used for different detectors
+ static const Int_t fgkNGapConfig = 4; // kBinNG = number of bins
+ static const Double_t fgkMinGapConfig = 1.; // kBinDG = lower border
+ static const Double_t fgkMaxGapConfig = 5.; // kBinNG + kBinDG = upper border
+
+ //-- mother kinematics
+ // invariant mass of the two-track system
+ static const Int_t fgkNMass = 1024; // number of bins
+ static const Double_t fgkMinMass = 0.; // lower border
+ static const Double_t fgkMaxMass = 5.12; // upper border
+
+ // transverse momentum of the two-track system
+ static const Int_t fgkNMotherPt = 128;//10; // number of bins
+ static const Double_t fgkMinMotherPt = 0.; // lower border
+ static const Double_t fgkMaxMotherPt = .64; //6.4 //2; // upper border
+
+ // cosine theta* (opening angle of the two daugther tracks in the
+ // centre-of-mass system of the two-track/mother system)
+ // **no meaning full variable**
+ static const Int_t fgkNCTS = 2; // number of bins
+ static const Double_t fgkMinCTS = -1.; // lower border
+ static const Double_t fgkMaxCTS = -0.9; // upper border
+
+ // opening angle in the lab frame
+ static const Int_t fgkNOA = 20; // number of bins
+ static const Double_t fgkMinOA = -1.; // lower border
+ static const Double_t fgkMaxOA = 1.; // upper border
+
+ //-- daughter kinematics
+ // transverse momentum of one of the two daughter particles
+ // (randomly selected)
+ static const Int_t fgkNDaughterPt = 128; // number of bins
+ static const Double_t fgkMinDaughterPt = 0.; // lower border
+ static const Double_t fgkMaxDaughterPt = 6.4; // upper border
+
+ // pseudo rapidity of one of the two daughter particles
+ // (randomly selected)
+ //static const Int_t fgkNDaughterEta = 64; // number of bins
+ //static const Double_t fgkMinDaughterEta = -1.28; // lower border
+ //static const Double_t fgkMaxDaughterEta = 1.28; // upper border
+
+ //-- Event quality information
+ // boolean values to reduce output size
+
+ // are there tracks in addition to the ones selected using AliCDMesonTracks
+ // (ITSTPC, ITSsa, ITSpureSA) (0 = no, 1 = yes)
+ static const Int_t fgkNTrackResiduals = 2; // number of bins
+ static const Double_t fgkMinTrackResiduals = 0.; // lower border
+ static const Double_t fgkMaxTrackResiduals = 2.; // upper border
+
+ // vertex with in +/-4cm (0 = no, 1 = yes)
+ static const Int_t fgkNVertexZinRng = 2; // number of bins
+ static const Double_t fgkMinVertexZinRng = 0.; // lower border
+ static const Double_t fgkMaxVertexZinRng = 2.; // upper border
+
+ // are the vertices from SPD and tracks within 0.5cm? (0 = no, 1 = yes)
+ static const Int_t fgkNVertexCoincidence = 2; // number of bins
+ static const Double_t fgkMinVertexCoincidence = 0.; // lower border
+ static const Double_t fgkMaxVertexCoincidence = 2.; // upper border
+
+ // are there SPD tracklets which are not assigned to tracks? (0 = no, 1 = yes)
+ static const Int_t fgkNTrackletResiduals = 2; // number of bins
+ static const Double_t fgkMinTrackletResiduals = 0.; // lower border
+ static const Double_t fgkMaxTrackletResiduals = 2.; // upper border
+
+ //-- MC event information
+ static const Int_t fgkNProcessType = 4; // kBinDD = number of bins
+ static const Double_t fgkMinProcessType = 0.; // kBinND = lower border
+ static const Double_t fgkMaxProcessType = 4.; // kBinDD = upper border
+
+
+ //== EMPTY EVENT STUDY (ThnEmptyEvents) ======================================
+ // event type
+ static const Int_t fgkNEventType = 5; // kBinEventE = number of bins (5)
+ static const Double_t fgkMinEventType = 1.; // kBinEventI = lower border (1)
+ static const Double_t fgkMaxEventType = 6.; // kBinEventE+kBinEventI = u. b.
+
+ // multiplicities (reused for different detectors and ways of counting)
+ static const Int_t fgkNMult = 32; // number of bins
+ static const Double_t fgkMinMult = 0.; // lower border
+ static const Double_t fgkMaxMult = 31.; // upper border
+
+ // multplicities - extended range
+ // (reused for different detectors and ways of counting)
+ static const Int_t fgkNMultW = 64; // number of bins
+ static const Double_t fgkMinMultW = 0; // lower border
+ static const Double_t fgkMaxMultW = 63; // upper border
+
+ //== MULTIPLICITY STUDY (TnnMultiplicity) ====================================
+ // number of ITSTPC tracks in event
+ static const Int_t fgkNNch = 51; // number of bins
+ static const Double_t fgkMinNch = 0.; // lower border
+ static const Double_t fgkMaxNch = 51.; // upper border
+
+ // number of ITS standalone tracks in event
+ static const Int_t fgkNNsoft = 11; // number of bins
+ static const Double_t fgkMinNsoft = 0.; // lower border
+ static const Double_t fgkMaxNsoft = 11.; // upper border
+
+ // combined multiplicity
+ static const Int_t fgkNNcomb = 61; // number of bins
+ static const Double_t fgkMinNcomb = 0.; // lower border
+ static const Double_t fgkMaxNcomb = 61.; // upper border
+
+ // gap configuration is reused from THnMother
+
+ // number of residual tracks
+ static const Int_t fgkNNresidualTracks = 11; // number of bins
+ static const Double_t fgkMinNresidualTracks = 0.; // lower border
+ static const Double_t fgkMaxNresidualTracks = 11.; // upper border
+
+ // number of residual tracklets
+ static const Int_t fgkNNresidualTracklets = 21; // number of bins
+ static const Double_t fgkMinNresidualTracklets = 0.; // lower border
+ static const Double_t fgkMaxNresidualTracklets = 21.; // upper border
+
+ // vertex z-position
+ static const Int_t fgkNVertexZ = 20; // number of bins
+ static const Double_t fgkMinVertexZ = -10.; // lower border
+ static const Double_t fgkMaxVertexZ = 10.; // upper border
+
+ // SPD and track vertex distance
+ static const Int_t fgkNVerticesDistance = 10; // number of bins
+ static const Double_t fgkMinVerticesDistance = 0.; // lower border
+ static const Double_t fgkMaxVerticesDistance = 5.; // upper border
+
+ // MC process type is resued from THnMother
+};
+
+#endif
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+//
+// AliCDMesonUtils
+// for
+// AliAnalysisTaskCDMeson
+//
+// Author:
+// Xianguo Lu <lu@physi.uni-heidelberg.de>
+// continued by
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#include <TH1.h>
+#include <TH2.h>
+#include <TGeoMatrix.h>
+#include <THnSparse.h>
+#include <TLorentzVector.h>
+#include <TRandom3.h>
+#include <TParticle.h>
+#include <TObjArray.h>
+
+#include "AliCDBManager.h"
+#include "AliCDBEntry.h"
+#include "AliCDBStorage.h"
+#include "AliESDEvent.h"
+#include "AliPIDResponse.h"
+#include "AliVTrack.h"
+#include "AliVParticle.h"
+#include "AliESDtrack.h"
+#include "AliESDtrackCuts.h"
+#include "AliESDFMD.h"
+#include "AliESDVZERO.h"
+#include "AliGeomManager.h"
+#include "AliITSAlignMille2Module.h"
+#include "AliITSsegmentationSPD.h"
+#include "AliMultiplicity.h"
+#include "AliPIDResponse.h"
+#include "AliSPDUtils.h"
+#include "AliTriggerAnalysis.h"
+#include "AliVZEROTriggerData.h"
+#include "AliVZEROCalibData.h"
+
+#include "AliAODTracklets.h"
+#include "AliAODEvent.h"
+
+#include "AliCDMesonBase.h"
+#include "AliCDMesonTracks.h"
+
+#include "AliCDMesonUtils.h"
+
+
+//==============================================================================
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::GetMassPtCtsOA(const Int_t pid,
+ const TVector3* momenta[], Double_t & mass,
+ Double_t &pt, Double_t &cts, Double_t &oa)
+{
+ //
+ // Get Mass Pt Cts OA
+ //
+
+ Double_t tmpp1[3], tmpp2[3];
+ momenta[0]->GetXYZ(tmpp1);
+ momenta[1]->GetXYZ(tmpp2);
+
+ Double_t masstrk = -999;
+ if(pid==AliCDMesonBase::kBinPionE || pid==AliCDMesonBase::kBinPion ||
+ pid==AliCDMesonBase::kBinSinglePion)
+ masstrk = AliPID::ParticleMass(AliPID::kPion);
+ else if(pid==AliCDMesonBase::kBinKaonE || pid==AliCDMesonBase::kBinKaon ||
+ pid==AliCDMesonBase::kBinSingleKaon)
+ masstrk = AliPID::ParticleMass(AliPID::kKaon);
+ else if(pid==AliCDMesonBase::kBinProtonE || pid==AliCDMesonBase::kBinProton ||
+ pid==AliCDMesonBase::kBinSingleProton)
+ masstrk = AliPID::ParticleMass(AliPID::kProton);
+ else if(pid==AliCDMesonBase::kBinElectronE ||
+ pid==AliCDMesonBase::kBinElectron ||
+ pid==AliCDMesonBase::kBinSingleElectron)
+ masstrk = AliPID::ParticleMass(AliPID::kElectron);
+ else
+ masstrk = AliPID::ParticleMass(AliPID::kPion);
+
+ const TLorentzVector sumv = GetKinematics(tmpp1, tmpp2, masstrk, masstrk,
+ cts);
+ mass = sumv.M();
+ pt = sumv.Pt();
+
+ oa = GetOA(tmpp1, tmpp2);
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::GetPWAinfo(const Int_t pid, const AliVTrack *trks[],
+ Float_t& theta, Float_t& phi, Float_t& mass,
+ Float_t momentum[])
+{
+ // transforms the coordiante system to the helicity frame and determines
+ // invariant mass and 3-vector of the two track system
+ //
+
+ Double_t tmpp1[3], tmpp2[3]; // 3-vectors of the daughter tracks
+ trks[0]->GetPxPyPz(tmpp1);
+ trks[1]->GetPxPyPz(tmpp2);
+
+ // determine mass of the daugther tracks
+ Double_t masstrk = -999;
+ if(pid==AliCDMesonBase::kBinPion || pid==AliCDMesonBase::kBinSinglePion)
+ masstrk = AliPID::ParticleMass(AliPID::kPion);
+ else if(pid==AliCDMesonBase::kBinKaon || pid==AliCDMesonBase::kBinSingleKaon)
+ masstrk = AliPID::ParticleMass(AliPID::kKaon);
+ else if(pid==AliCDMesonBase::kBinProton ||
+ pid==AliCDMesonBase::kBinSingleProton)
+ masstrk = AliPID::ParticleMass(AliPID::kProton);
+ else if(pid==AliCDMesonBase::kBinElectron ||
+ pid==AliCDMesonBase::kBinSingleElectron)
+ masstrk = AliPID::ParticleMass(AliPID::kElectron);
+ else
+ masstrk = AliPID::ParticleMass(AliPID::kPion);
+
+ TLorentzVector va, vb; // 4-vectors of the two tracks
+ va.SetXYZM(tmpp1[0], tmpp1[1], tmpp1[2], masstrk);
+ vb.SetXYZM(tmpp2[0], tmpp2[1], tmpp2[2], masstrk);
+
+ TLorentzVector sumv = va+vb; // 4-vector of the pair
+ TVector3 sumXZ(sumv.X(), 0, sumv.Z()); // its projection to the xz-plane
+
+ // mass and momentum
+ mass = sumv.M();
+ momentum[0] = sumv.X();
+ momentum[1] = sumv.Y();
+ momentum[2] = sumv.Z();
+
+ // coordinate axis vectors
+ const TVector3 x(1.,0,0);
+ const TVector3 y(0,1.,0);
+ const TVector3 z(0,0,1.);
+
+ const Double_t alpha = -z.Angle(sumXZ);
+
+ sumXZ.Rotate(alpha, y);
+ sumv.Rotate(alpha, y);
+ va.Rotate(alpha, y);
+ vb.Rotate(alpha, y);
+
+ const Double_t beta = z.Angle(sumv.Vect());
+ sumv.Rotate(beta, x);
+ va.Rotate(beta, x);
+ vb.Rotate(beta, x);
+
+
+ const TVector3 bv = -sumv.BoostVector();
+ va.Boost(bv);
+ vb.Boost(bv);
+
+ // angles
+ theta = va.Theta();
+ phi = va.Phi();
+}
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetEventType(const AliESDEvent *ESDEvent)
+{
+ // checks of which type a event is:
+ // beam-beam interaction (I), beam-gas (A/C), empty (E)
+ //
+
+ TString firedTriggerClasses = ESDEvent->GetFiredTriggerClasses();
+
+ if (firedTriggerClasses.Contains("CINT1A-ABCE-NOPF-ALL")) { // A
+ return AliCDMesonBase::kBinEventA;
+ }
+ if (firedTriggerClasses.Contains("CINT1C-ABCE-NOPF-ALL")) { // C
+ return AliCDMesonBase::kBinEventC;
+ }
+ if (firedTriggerClasses.Contains("CINT1B-ABCE-NOPF-ALL")) { // I
+ return AliCDMesonBase::kBinEventI;
+ }
+ if (firedTriggerClasses.Contains("CINT1-E-NOPF-ALL")) { // E
+ return AliCDMesonBase::kBinEventE;
+ }
+ if (firedTriggerClasses.Contains("CDG5-E")) { // E
+ return AliCDMesonBase::kBinEventE;
+ }
+ if (firedTriggerClasses.Contains("CDG5-I")) { // I
+ return AliCDMesonBase::kBinEventI;
+ }
+ if (firedTriggerClasses.Contains("CDG5-AC")) { // AC
+ return AliCDMesonBase::kBinEventAC;
+ }
+ return AliCDMesonBase::kBinEventUnknown;
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::SwapTrack(const AliVTrack* trks[])
+{
+ //
+ // swap two esdtracks, needed since only the properties of one a stored and
+ // AliRoot sorts them
+ //
+
+ TRandom3 tmprd(0);
+ if(tmprd.Rndm()>0.5)
+ return;
+
+ const AliVTrack* tmpt = trks[0];
+ trks[0]=trks[1];
+ trks[1]=tmpt;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetCombCh(const AliVTrack * trks[])
+{
+ //
+ // get combination of charges
+ //
+
+ const Double_t ch1 = trks[0]->Charge();
+ const Double_t ch2 = trks[1]->Charge();
+
+ if(ch1*ch2<0){
+ return AliCDMesonBase::kBinPM;
+ }
+ else
+ return AliCDMesonBase::kBinPPMM;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetCombPID(AliPIDResponse* pid, const AliVTrack* trks[],
+ const Int_t mode, TH2* comb2trkPID /*= 0x0 */)
+{
+ //
+ // Get PID for a two track system (data)
+ //
+
+ if (!pid){
+ return AliCDMesonBase::kBinPIDUnknown;
+ }
+
+ // get PID for the single tracks
+ Int_t kpid[2];
+ for(Int_t ii=0; ii<2; ii++){
+ kpid[ii] = GetPID(pid, trks[ii], mode);
+ }
+
+ // PID QA histogram
+ if (comb2trkPID) {
+ comb2trkPID->Fill(kpid[0], kpid[1]);
+ }
+
+ return CombinePID(kpid);
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetCombPID(const TParticle* particles[],
+ TH2 *comb2trkPID /* = 0x0 */)
+{
+ //
+ // Get PID for a two track system (MC)
+ //
+
+ // get PID for the single tracks
+ Int_t kpid[2];
+ for(Int_t ii=0; ii<2; ii++){
+ kpid[ii] = GetPID(particles[ii]->GetPdgCode());
+ }
+
+ // PID QA histogram
+ if (comb2trkPID) {
+ comb2trkPID->Fill(kpid[0], kpid[1]);
+ }
+
+ return CombinePID(kpid);
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::GetSPDTrackletMult(const AliESDEvent *ESDEvent,
+ Int_t& sum, Int_t& forwardA,
+ Int_t& forwardC, Int_t& central)
+{
+ // obtain the multiplicity for eta < -.9 (forwardC), eta > 0.9 (fowardA), in
+ // the central barrel (central) and its sum from SPD tracklets with the
+ // AliMultiplicity class
+ // code simular to PWGPP/ITS/AliAnalysisTaskSPD.cxx
+
+ sum = forwardA = forwardC = central = 0; // initialize values
+
+ const AliMultiplicity *mult = ESDEvent->GetMultiplicity();
+ for (Int_t iTracklet = 0; iTracklet < mult->GetNumberOfTracklets();
+ iTracklet++) {
+ float_t eta = mult->GetEta(iTracklet);
+ if (eta < -0.9) {
+ forwardC++;
+ }
+ else if (eta < 0.9) {
+ central++;
+ }
+ else {
+ forwardA++;
+ }
+ sum++;
+ }
+}
+
+
+//------------------------------------------------------------------------------
+Bool_t AliCDMesonUtils::CutEvent(const AliESDEvent *ESDEvent, TH1 *hspd,
+ TH1 *hpriv, TH1 *hpriVtxPos, TH1 *hpriVtxDist,
+ TH2 *hfo, TH1* hfochans, Int_t &kfo,
+ Int_t &nip, Int_t &nop, TH1 *hpriVtxX,
+ TH1 *hpriVtxY, TH1 *hpriVtxZ)
+{
+ //
+ // CutEvent
+ //
+
+ AliTriggerAnalysis triggerAnalysis;
+ /*
+ //printf("freidtlog active triggers: %s\n",
+ // ESDEvent->GetHeader()->GetActiveTriggerInputs().Data());
+ //printf("freidtlog fired triggers: %s\n",
+ // ESDEvent->GetHeader()->GetFiredTriggerInputs().Data());
+ //*/
+ //http://alisoft.cern.ch/viewvc/trunk/ANALYSIS/AliTriggerAnalysis.cxx?view=markup&root=AliRoot
+ //Int_t AliTriggerAnalysis::SPDFiredChips(const AliESDEvent* aEsd, Int_t origin, Bool_t fillHists, Int_t layer)
+ // returns the number of fired chips in the SPD
+ // origin = 0 --> aEsd->GetMultiplicity()->GetNumberOfFiredChips() (filled from clusters)
+ // origin = 1 --> aEsd->GetMultiplicity()->TestFastOrFiredChips() (from hardware bits)
+
+
+ // SPD FastOR cut
+ const Int_t fastORhw = triggerAnalysis.SPDFiredChips(ESDEvent, 1);
+ if (hspd) hspd->Fill(fastORhw);
+ /*
+ if(fastORhw<1)
+ return kFALSE;
+ */
+
+ if (hfochans) {
+ const AliMultiplicity *mult = ESDEvent->GetMultiplicity();
+
+ for(Int_t iChipKey=0; iChipKey < 1200; iChipKey++){
+ if (mult->TestFastOrFiredChips(iChipKey)) {
+ hfochans->Fill((Double_t)iChipKey);
+ }
+ }
+ }
+
+ if(kfo){ // spd trigger study
+ Int_t nfoctr[10];
+ GetNFO(ESDEvent, "[1.0]", nfoctr, 0x0, 0x0);
+ nip = nfoctr[kInnerPixel];
+ nop = nfoctr[kOuterPixel];
+
+ if(AliCDMesonBase::GetGapBin("V0", GetV0(ESDEvent)) ==
+ AliCDMesonBase::kBinDG) {
+ hfo->Fill(nip, nop);
+ }
+
+ if(nop<2)
+ kfo = 1;
+ else
+ kfo = nop;
+
+ if(kfo>=10)
+ kfo=9;
+ }
+
+ // collision vertex cut
+ // A cut in XY is implicitly done during the reconstruction by constraining
+ // the vertex to the beam diamond.
+
+ // Primary vertex
+ Bool_t kpr0 = kTRUE;
+ const AliESDVertex *vertex = ESDEvent->GetPrimaryVertexTracks();
+ if(vertex->GetNContributors()<1) {
+ // SPD vertex
+ vertex = ESDEvent->GetPrimaryVertexSPD();
+ if(vertex->GetNContributors()<1) {
+ // NO VERTEX, SKIP EVENT
+ kpr0 = kFALSE;
+ }
+ }
+ const Bool_t kpriv = kpr0 && (fabs(vertex->GetZ()) < 10.);
+ // 10 is the common value, unit: cm
+ if (hpriv) hpriv->Fill(kpriv);
+ if (hpriVtxDist) hpriVtxDist->Fill(vertex->GetZ());
+ if (hpriVtxPos && kpr0) hpriVtxPos->Fill(!kpriv);
+ if(!kpriv)
+ return kFALSE;
+
+ if(hpriVtxX) hpriVtxX->Fill(vertex->GetX());
+ if(hpriVtxY) hpriVtxY->Fill(vertex->GetY());
+ if(hpriVtxZ) hpriVtxZ->Fill(vertex->GetZ());
+ return kTRUE;
+}
+
+//------------------------------------------------------------------------------
+Bool_t AliCDMesonUtils::CutEvent(const AliAODEvent *AODEvent, TH1 *hpriv,
+ TH1 *hpriVtxX, TH1 *hpriVtxY, TH1 *hpriVtxZ,
+ TH1 *hpriVtxPos, TH1 *hpriVtxDist)
+{
+ //
+ // Cut Event for AOD Events, to be combined with the ESD Track Cut
+ //
+
+ // TODO: no idea about fast or yet, to be thought of
+
+ // Primary vertex
+ Bool_t kpr0 = kTRUE;
+ const AliAODVertex *vertex = AODEvent->GetPrimaryVertex();
+ if(vertex->GetNContributors()<1) {
+ // SPD vertex
+ vertex = AODEvent->GetPrimaryVertexSPD();
+ if(vertex->GetNContributors()<1) {
+ // NO GOOD VERTEX, SKIP EVENT
+ kpr0 = kFALSE;
+ }
+ }
+ const Bool_t kpriv = kpr0 && (fabs(vertex->GetZ())<10.);
+ // 10 is the common value, unit: cm
+ hpriv->Fill(kpriv);
+ if (hpriVtxDist) hpriVtxDist->Fill(vertex->GetZ());
+ if (hpriVtxPos && kpr0) hpriVtxPos->Fill(!kpriv);
+ if(!kpriv)
+ return kFALSE;
+
+ if(hpriVtxX) hpriVtxX->Fill(vertex->GetX());
+ if(hpriVtxY) hpriVtxY->Fill(vertex->GetY());
+ if(hpriVtxZ) hpriVtxZ->Fill(vertex->GetZ());
+ return kTRUE;
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::DoVZEROStudy(const AliESDEvent *ESDEvent,
+ TObjArray* hists, const Int_t run)
+{
+ //
+ //
+ // IMPORTANT: the order of the histograms here and in
+ // AliCDMesonBase::GetHistVZEROStudies(..) has to match
+
+ const AliESDVZERO* esdV0 = ESDEvent->GetVZEROData();
+ if (!esdV0) {
+ Printf("ERROR: esd V0 not available");
+ return;
+ }
+
+ // determine trigger decision
+ AliTriggerAnalysis triggerAnalysis;
+ //const Bool_t khw = kFALSE;
+
+ //Float_t v0a = (triggerAnalysis.V0Trigger(ESDEvent, AliTriggerAnalysis::kASide, khw) ==
+ // AliTriggerAnalysis::kV0BB) ? 1. : 0.;
+ //Float_t v0c = (triggerAnalysis.V0Trigger(ESDEvent, AliTriggerAnalysis::kCSide, khw) ==
+ // AliTriggerAnalysis::kV0BB) ? 1. : 0.;
+
+ // obtain OCDB objects
+ AliCDBManager *man = AliCDBManager::Instance();
+ TString cdbpath;
+ if (man->IsDefaultStorageSet()) {
+ const AliCDBStorage *dsto = man->GetDefaultStorage();
+ cdbpath = TString(dsto->GetBaseFolder());
+ }
+ else { //should not be used!
+ man->SetDefaultStorage(gSystem->Getenv("TRAIN_CDB_PATH"));
+ cdbpath = TString(gSystem->Getenv("TRAIN_CDB_PATH"));
+ }
+ man->SetSpecificStorage("VZERO/Trigger/Data",cdbpath);
+ man->SetSpecificStorage("VZERO/Calib/Data",cdbpath);
+ man->SetRun(run);
+
+ AliCDBEntry *ent1 = man->Get("VZERO/Trigger/Data");
+ AliVZEROTriggerData *trigData = (AliVZEROTriggerData*)ent1->GetObject();
+ if (!trigData) {
+ printf("freidtlog failed loading VZERO trigger data from OCDB\n");
+ return;
+ }
+
+ AliCDBEntry *ent2 = man->Get("VZERO/Calib/Data");
+ AliVZEROCalibData *calData = (AliVZEROCalibData*)ent2->GetObject();
+ if (!calData) {
+ printf("freidtlog failed loading VZERO calibration data from OCDB\n");
+ return;
+ }
+
+ // fill histograms
+ Int_t pmtHist = 0;
+ Int_t iPMT = 0;
+ for (iPMT = 0; iPMT < 64; ++iPMT) {
+ Int_t board = AliVZEROCalibData::GetBoardNumber(iPMT);
+ Int_t channel = AliVZEROCalibData::GetFEEChannelNumber(iPMT);
+ ((TH2*)hists->At(iPMT+pmtHist))->Fill(esdV0->GetAdc(iPMT),
+ trigData->GetPedestalCut(0, board, channel));
+ // calData->GetDiscriThr(iPMT));
+ }
+ pmtHist = iPMT;
+ for (iPMT = 0; iPMT < 64; ++iPMT) {
+ ((TH2*)hists->At(iPMT+pmtHist))->Fill(esdV0->GetAdc(iPMT),
+ esdV0->GetMultiplicity(iPMT));
+ }
+ /*
+ // not used in 2010 for pp
+ ((TH2*)hists->At(pmtHist++))->Fill((Float_t)esdV0->GetTriggerChargeA(), v0a);
+ ((TH2*)hists->At(pmtHist++))->Fill((Float_t)esdV0->GetTriggerChargeC(), v0c);
+ */
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetGapConfig(const AliESDEvent *ESDEvent,
+ TH2 *hitMapSPDinner, TH2 *hitMapSPDouter,
+ TH2 *hitMapSPDtrklt, TH2 *hitMapFMDa,
+ TH2 *hitMapFMDc, TH1 **fmdSums,
+ TH2 *TPCGapDCAaSide, TH2 *TPCGapDCAcSide)
+{
+ //
+ // GetGapConfigAndTracks
+ //
+ // retrieves the gap configuration of a track and returns it as
+ // an bit vector
+ // kBaseLine ensures, that this event is valid
+ // + is equivalent to | in this case
+ Int_t gapConfig = AliCDMesonBase::kBitBaseLine + GetV0(ESDEvent)
+ + GetFMD(ESDEvent, hitMapFMDa, hitMapFMDc, fmdSums)
+ + GetSPD(ESDEvent, hitMapSPDinner, hitMapSPDouter, hitMapSPDtrklt);
+ if (gapConfig == AliCDMesonBase::kBitBaseLine) {
+ gapConfig += GetTPC(ESDEvent, TPCGapDCAaSide, TPCGapDCAcSide);
+ }
+ else {
+ gapConfig += GetTPC(ESDEvent, 0x0, 0x0);
+ }
+ if (GetFastORmultiplicity(ESDEvent) > 0) {
+ gapConfig += AliCDMesonBase::kBitCentAct;
+ }
+ return gapConfig; // + GetZDC(ESDEvent);
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::FillEtaPhiMap(const AliVEvent *event,
+ const AliCDMesonTracks* tracks, TH2 *map,
+ TH2 *map_c)
+{
+ //
+ // Fills the eta phi information about all tracks and about the tracks surving
+ // the track cuts (CutTrack) into two separate histograms
+ //
+
+ // all tracks
+ for (Int_t itrk = 0; itrk < event->GetNumberOfTracks(); ++itrk) {
+ if (AliVParticle *trk = event->GetTrack(itrk)) {
+ if (map) {
+ map->Fill(trk->Eta(), trk->Phi());
+ }
+ }
+ }
+
+ // tracks that survived the cuts
+ for (Int_t itrk = 0; itrk < tracks->GetCombinedTracks(); ++itrk) {
+ if (map_c) {
+ map_c->Fill(tracks->GetTrack(itrk)->Eta(), tracks->GetTrack(itrk)->Phi());
+ }
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::GetMultFMD(const AliESDEvent *ESDEvent, Int_t& fmdA,
+ Int_t& fmdC, Float_t *fmdSums)
+{
+ //
+ // Multiplicity seen by FMD
+ //
+ // WARNING: this function is only working with a modified AliRoot so far
+
+#ifdef STD_ALIROOT
+ fmdA = FMDHitCombinations(ESDEvent, 0);
+ fmdC = FMDHitCombinations(ESDEvent, 1);
+#else
+ AliTriggerAnalysis triggerAnalysis;
+ triggerAnalysis.SetFMDThreshold(0.3, 0.5); // parameters got from FMD
+ triggerAnalysis.FMDTrigger(ESDEvent, AliTriggerAnalysis::kASide, &fmdA);
+ triggerAnalysis.FMDTrigger(ESDEvent, AliTriggerAnalysis::kCSide, &fmdC);
+#endif
+
+ if (fmdSums) {
+ const AliESDFMD* fmdData =
+ (const_cast<AliESDEvent*>(ESDEvent))->GetFMDData();
+ for (UShort_t det=1; det<=3; det++) {
+ Int_t nRings = (det == 1 ? 1 : 2);
+ for (UShort_t ir = 0; ir < nRings; ir++) {
+ Char_t ring = (ir == 0 ? 'I' : 'O');
+ UShort_t nsec = (ir == 0 ? 20 : 40);
+ UShort_t nstr = (ir == 0 ? 512 : 256);
+ for (UShort_t sec =0; sec < nsec; sec++) {
+ for (UShort_t strip = 0; strip < nstr; strip++) {
+ Float_t mult = fmdData->Multiplicity(det,ring,sec,strip);
+ if (mult == AliESDFMD::kInvalidMult) continue;
+
+ if (det == 1 && ring == 'I') fmdSums[0] += mult;
+ else if (det == 2 && ring == 'I') fmdSums[1] += mult;
+ else if (det == 2 && ring == 'O') fmdSums[2] += mult;
+ else if (det == 3 && ring == 'I') fmdSums[3] += mult;
+ else if (det == 3 && ring == 'O') fmdSums[4] += mult;
+ }
+ }
+ }
+ }
+ }
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::GetMultSPD(const AliESDEvent * ESDEvent, Int_t& spdIA,
+ Int_t& spdIC, Int_t& spdOA, Int_t& spdOC)
+{
+ //
+ // Retrieves the multiplicity seen by the SPD FastOR
+ //
+
+ Int_t nfoctr[10];
+ GetNFO(ESDEvent, "]0.9[", nfoctr, 0x0, 0x0);
+
+ spdIA = nfoctr[kIPA];
+ spdIC = nfoctr[kIPC];
+ spdOA = nfoctr[kOPA];
+ spdOC = nfoctr[kOPC];
+}
+
+
+//==============================================================================
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetV0(const AliESDEvent * ESDEvent)
+{
+ //
+ //GetV0
+ //
+
+ AliTriggerAnalysis triggerAnalysis;
+ const Bool_t khw = kFALSE;
+ const Bool_t v0A =
+ (triggerAnalysis.V0Trigger(ESDEvent, AliTriggerAnalysis::kASide, khw) ==
+ AliTriggerAnalysis::kV0BB);
+ const Bool_t v0C =
+ (triggerAnalysis.V0Trigger(ESDEvent, AliTriggerAnalysis::kCSide, khw) ==
+ AliTriggerAnalysis::kV0BB);
+
+ return v0A * AliCDMesonBase::kBitV0A + v0C * AliCDMesonBase::kBitV0C;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetFMD(const AliESDEvent *ESDEvent, TH2 *hitMapFMDa,
+ TH2 *hitMapFMDc, TH1 **fmdSums)
+{
+ //
+ // GetFMD
+ //
+
+ AliTriggerAnalysis triggerAnalysis;
+ triggerAnalysis.SetFMDThreshold(0.3, 0.5); // parameters got from FMD
+ const Bool_t fmdA =
+ triggerAnalysis.FMDTrigger(ESDEvent, AliTriggerAnalysis::kASide);
+ const Bool_t fmdC =
+ triggerAnalysis.FMDTrigger(ESDEvent, AliTriggerAnalysis::kCSide);
+
+ //printf("FR - GetFMD\n");
+
+ // prepartions for a charge summation algorithm
+ Bool_t hitMaps = (Bool_t)(hitMapFMDa && hitMapFMDc);
+ Bool_t calcSum = fmdSums ?
+ (fmdSums[0] && fmdSums[1] && fmdSums[2] && fmdSums[3] && fmdSums[4]) :
+ kFALSE; // are the histograms defined?
+ Float_t sum[] = { 0., 0., 0., 0., 0. };
+ // summed multiplicity in the FMD detectors
+
+ // hit map generation
+ if (hitMaps || calcSum) {
+ const AliESDFMD* fmdData =
+ (const_cast<AliESDEvent*>(ESDEvent))->GetFMDData();
+
+ for (UShort_t det=1; det<=3;det++) {
+ Int_t nRings = (det == 1 ? 1 : 2);
+ for (UShort_t ir = 0; ir < nRings; ir++) {
+ Char_t ring = (ir == 0 ? 'I' : 'O');
+ UShort_t nsec = (ir == 0 ? 20 : 40);
+ UShort_t nstr = (ir == 0 ? 512 : 256);
+ for (UShort_t sec =0; sec < nsec; sec++) {
+ for (UShort_t strip = 0; strip < nstr; strip++) {
+ Float_t mult = fmdData->Multiplicity(det,ring,sec,strip);
+ if (mult == AliESDFMD::kInvalidMult) continue;
+
+ if (calcSum) {
+ if (det == 1 && ring == 'I') sum[0] += mult;
+ else if (det == 2 && ring == 'I') sum[1] += mult;
+ else if (det == 2 && ring == 'O') sum[2] += mult;
+ else if (det == 3 && ring == 'I') sum[3] += mult;
+ else if (det == 3 && ring == 'O') sum[4] += mult;
+ }
+
+ if (hitMaps) { // care about hit map specific information
+ const Float_t eta = fmdData->Eta(det,ring,sec,strip);
+ const Float_t phi =
+ fmdData->Phi(det,ring,sec,strip) / 180. * TMath::Pi();
+ //printf("FR - GetFMD: %f %f %f\n", eta, phi, mult);
+ if (eta != AliESDFMD::kInvalidEta) {
+ if ((-3.5 < eta) && (eta < -1.5)) {
+ hitMapFMDc->Fill(eta, phi, mult); // use mult as weight
+ }
+ else if ((1.5 < eta) && (eta < 5.5)) {
+ hitMapFMDa->Fill(eta, phi, mult); // use mult as weight
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+
+ if (calcSum) {
+ //printf("DEBUG -- SUM(%f,%f,%f,%f,%f)\n", sum[0], sum[1], sum[2], sum[3],
+ // sum[4]);
+ for (UInt_t i = 0; i < 5; i++) { //
+ fmdSums[i]->Fill(sum[i]);
+ }
+ }
+
+ return fmdA * AliCDMesonBase::kBitFMDA + fmdC * AliCDMesonBase::kBitFMDC;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetSPD(const AliESDEvent *ESDEvent, TH2 *hitMapSPDinner,
+ TH2 *hitMapSPDouter, TH2 *hitMapSPDtrklt)
+{
+ //
+ // GetSPD
+ //
+
+ Int_t nfoctr[10];
+ GetNFO(ESDEvent, "]0.9[", nfoctr, hitMapSPDinner, hitMapSPDouter);
+ // get multiplicity from fastOR and fill corresponding hit maps
+
+ if (hitMapSPDtrklt) FillSPDtrkltMap(ESDEvent, hitMapSPDtrklt);
+ // fill tracklet hit map
+
+ const Int_t ipA = nfoctr[kIPA]; // inner layer A side
+ const Int_t ipC = nfoctr[kIPC]; // inner layer C side
+ const Int_t opA = nfoctr[kOPA]; // outer layer A side
+ const Int_t opC = nfoctr[kOPC]; // outer layer C side
+
+ const Bool_t spdA = ipA + opA; // A side hit?
+ const Bool_t spdC = ipC + opC; // C side hit?
+
+ return spdA * AliCDMesonBase::kBitSPDA + spdC * AliCDMesonBase::kBitSPDC;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetTPC(const AliESDEvent * ESDEvent, TH2 *TPCGapDCAaSide,
+ TH2 *TPCGapDCAcSide)
+{
+ //
+ //GetTPC
+ //
+
+ const Double_t etacut = 0.9;
+ Int_t nA = 0;
+ Int_t nC = 0;
+
+ AliESDtrackCuts cuts;
+ cuts.SetMaxDCAToVertexXY(0.1);
+ cuts.SetMaxDCAToVertexZ(2.);
+
+ for(Int_t itrack = 0; itrack < ESDEvent->GetNumberOfTracks(); itrack++){
+ const AliESDtrack* esdtrack = ESDEvent->GetTrack(itrack);
+ Float_t b[2];
+ Float_t bCov[3];
+ esdtrack->GetImpactParameters(b,bCov);
+
+ if (bCov[0]<=0 || bCov[2]<=0) {
+ printf("freidtlog - Estimated b resolution lower or equal zero!\n");
+ bCov[0]=0;
+ bCov[2]=0;
+ }
+
+ Float_t dcaToVertexXY = b[0];
+ Float_t dcaToVertexZ = b[1];
+ if (esdtrack->Eta() > etacut) {
+ if (cuts.AcceptTrack(esdtrack)) nA++;
+ if (TPCGapDCAaSide) {
+ TPCGapDCAaSide->Fill(dcaToVertexXY, dcaToVertexZ);
+ }
+ }
+ else if (esdtrack->Eta() < -etacut) {
+ if (cuts.AcceptTrack(esdtrack)) nC++;
+ if (TPCGapDCAcSide) {
+ TPCGapDCAcSide->Fill(dcaToVertexXY, dcaToVertexZ);
+ }
+ }
+ }
+
+ const Bool_t tpcA = nA;
+ const Bool_t tpcC = nC;
+
+ return tpcA * AliCDMesonBase::kBitTPCA + tpcC * AliCDMesonBase::kBitTPCC;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetZDC(const AliESDEvent * ESDEvent)
+{
+ //
+ //GetZDC
+ //
+
+ const Int_t qa = ESDEvent->GetESDZDC()->GetESDQuality();
+ Bool_t zdcA = kFALSE, zdcC = kFALSE;
+ for(Int_t ii=0; ii<6; ii++){
+ if(qa & (1<<ii)){
+ if(ii<4) zdcA = kTRUE;
+ else zdcC = kTRUE;
+ }
+ }
+
+ return zdcA * AliCDMesonBase::kBitZDCA + zdcC * AliCDMesonBase::kBitZDCC;
+}
+
+
+//==============================================================================
+//------------------------------------------------------------------------------
+#ifdef STD_ALIROOT
+Int_t AliCDMesonUtils::FMDHitCombinations(const AliESDEvent* ESDEvent,
+ Int_t side)
+{
+ //
+ // copy of the FMDHitCombinations function originating from AliTriggerAnalysis
+ //
+ // side == 0 -> A side, side > 0 -> C side
+
+ // workaround for AliESDEvent::GetFMDData is not const!
+ const AliESDFMD* fmdData =
+ (const_cast<AliESDEvent*>(ESDEvent))->GetFMDData();
+ if (!fmdData)
+ {
+ puts("AliESDFMD not available");
+ return -1;
+ }
+
+ Int_t detFrom = (side == 0) ? 1 : 3;
+ Int_t detTo = (side == 0) ? 2 : 3;
+
+ Int_t triggers = 0;
+ Float_t totalMult = 0;
+ for (UShort_t det=detFrom;det<=detTo;det++) {
+ Int_t nRings = (det == 1 ? 1 : 2);
+ for (UShort_t ir = 0; ir < nRings; ir++) {
+ Char_t ring = (ir == 0 ? 'I' : 'O');
+ UShort_t nsec = (ir == 0 ? 20 : 40);
+ UShort_t nstr = (ir == 0 ? 512 : 256);
+ for (UShort_t sec =0; sec < nsec; sec++) {
+ for (UShort_t strip = 0; strip < nstr; strip++) {
+ Float_t mult = fmdData->Multiplicity(det,ring,sec,strip);
+ if (mult == AliESDFMD::kInvalidMult) continue;
+
+ if (mult > 0.3) // fFMDLowCut
+ totalMult = totalMult + mult;
+ else {
+ if (totalMult > 0.5) // fFMDHitCut
+ triggers++;
+ }
+ }
+ }
+ }
+ }
+ return triggers;
+}
+#endif
+
+
+//==============================================================================
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::SPDLoadGeom(const Int_t run)
+{
+ // method to get the gGeomanager
+ // it is called at the CreatedOutputObject stage
+ // to comply with the CAF environment
+
+ AliCDBManager *man = AliCDBManager::Instance();
+
+ TString cdbpath;
+ if (man->IsDefaultStorageSet()) {
+ const AliCDBStorage *dsto = man->GetDefaultStorage();
+ cdbpath = TString(dsto->GetBaseFolder());
+ //printf("default was set to: %s\n", cdbpath.Data());
+ }
+ else { //should not be used!
+ // man->SetDefaultStorage("alien://folder=/alice/data/2010/OCDB");
+ // would be needed on grid
+ man->SetDefaultStorage(gSystem->Getenv("TRAIN_CDB_PATH"));
+ cdbpath = TString(gSystem->Getenv("TRAIN_CDB_PATH"));
+ }
+
+ man->SetSpecificStorage("ITS/Align/Data",cdbpath);
+ man->SetSpecificStorage("GRP/Geometry/Data",cdbpath);
+ man->SetRun(run);
+
+ AliCDBEntry* obj = man->Get(AliCDBPath("GRP", "Geometry", "Data"));
+ if (!obj) {
+ printf("freidtlog failed loading geometry object\n");
+ return;
+ }
+ AliGeomManager::SetGeometry((TGeoManager*)obj->GetObject());
+ AliGeomManager::ApplyAlignObjsFromCDB("ITS");
+}
+
+//------------------------------------------------------------------------------
+Bool_t AliCDMesonUtils::SPDLoc2Glo(const Int_t id, const Double_t *loc,
+ Double_t *glo)
+{
+ //
+ //SPDLoc2Glo, do not touch
+ //
+
+ static TGeoHMatrix mat;
+ Int_t vid = AliITSAlignMille2Module::GetVolumeIDFromIndex(id);
+ if (vid<0) {
+ printf("freidtlog Did not find module with such ID %d\n",id);
+ return kFALSE;
+ }
+ AliITSAlignMille2Module::SensVolMatrix(vid,&mat);
+ mat.LocalToMaster(loc,glo);
+ return kTRUE;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::CheckChipEta(const Int_t chipKey, const TString scut,
+ const Double_t vtxPos[],
+ TH2 *hitMapSPDinner, TH2 *hitMapSPDouter)
+{
+ //
+ //CheckChipEta
+ //
+
+ // retrieves the position in eta for a given chip and applies the cut
+ // results:
+ // 0 <= out of range
+ // -1 <= negative pseudo-rapidity position, in range (C-Side)
+ // 1 <= positive pseudo-rapidity position, in range (A-Side)
+ //
+ // scut: "[0.9" or "]0.9", only 3 digits for the value!!
+
+
+ const Bool_t kincl = (scut[0] == '[');
+ const TString cutval = scut(1,3);
+ const Double_t etacut = fabs(cutval.Atof());
+
+ //no eta cut, save time
+ if(kincl && etacut>=2)
+ return kTRUE;
+
+ Int_t etaside = 1;
+ //------------------------------- NOT TO TOUCH ------------------------>>
+ UInt_t module=999, offchip=999;
+ AliSPDUtils::GetOfflineFromOfflineChipKey(chipKey,module,offchip);
+ UInt_t hs = AliSPDUtils::GetOnlineHSFromOffline(module);
+ if(hs<2) offchip = 4 - offchip; // inversion in the inner layer...
+
+ const Int_t col[]={
+ hs<2? 0 : 31,
+ hs<2? 31 : 0,
+ hs<2? 31 : 0,
+ hs<2? 0 : 31};
+ const Int_t aa[]={0, 0, 255, 255};
+ const AliITSsegmentationSPD seg;
+
+ for(Int_t ic=0; ic<4; ic++){
+ Float_t localchip[3]={0.,0.,0.};
+ seg.DetToLocal(aa[ic],col[ic]+32*offchip,localchip[0],localchip[2]);
+ // local coordinate of the chip center
+ //printf("local coordinates %d %d: %f %f \n",chipKey, ic, localchip[0],localchip[2]);
+ const Double_t local[3] = {localchip[0],localchip[1],localchip[2]};
+ Double_t glochip[3]={0.,0.,0.};
+ if(!SPDLoc2Glo(module,local,glochip)){
+ return kFALSE;
+ }
+
+ //-------------------------------------------------------------------<<
+
+ const TVector3 pos(glochip[0]-vtxPos[0], glochip[1]-vtxPos[1],
+ glochip[2]-vtxPos[2]);
+ //pos.Print();
+
+ if (chipKey < 400) { // inner SPD layer
+ if (hitMapSPDinner) {
+ Double_t phi = pos.Phi(); // output in the range -Pi +Pi
+ if (phi < 0.) phi += TMath::TwoPi(); // remap to the interval [0, TwoPi)
+ const Double_t eta = pos.Eta();
+ hitMapSPDinner->Fill(eta, phi);
+ }
+ }
+ else {
+ if (hitMapSPDouter) { // outer SPD layer
+ Double_t phi = pos.Phi(); // output in the range -Pi +Pi
+ if (phi < 0.) phi += TMath::TwoPi(); // remap to the interval [0, TwoPi)
+ const Double_t eta = pos.Eta();
+ hitMapSPDouter->Fill(eta, phi);
+ }
+ }
+
+ if( kincl && fabs(pos.Eta()) > etacut)
+ return kFALSE;
+
+ if(!kincl){
+ if(fabs(pos.Eta()) < etacut)
+ return kFALSE;
+ else if(pos.Eta()<0)
+ etaside = -1;
+ else
+ etaside = 1;
+ }
+ }
+
+ return etaside;
+}
+
+
+//------------------------------------------------------------------------------
+void AliCDMesonUtils::GetNFO(const AliESDEvent *ESDEvent, const TString etacut,
+ Int_t ctr[], TH2 *hitMapSPDinner,
+ TH2 *hitMapSPDouter)
+{
+ //
+ // GetNFO
+ //
+ // analyzes the SPD fastOR for a given eta range and returns
+ // an array with the number of hits in:
+
+ Int_t ninner=0; // inner layer
+ Int_t nouter=0; // outer layer
+ Int_t ipA = 0; // inner layer A side
+ Int_t ipC = 0; // inner layer C side
+ Int_t opA = 0; // outer layer A side
+ Int_t opC = 0; // outer layer C side
+
+ const AliMultiplicity *mult = ESDEvent->GetMultiplicity();
+
+ // position of the primary vertex
+ Double_t tmp[3] = { 0., 0., 0. };
+ ESDEvent->GetPrimaryVertex()->GetXYZ(tmp);
+ Double_t vtxPos[3] = { tmp[0], tmp[1], tmp[2] };
+
+
+ for(Int_t iChipKey=0; iChipKey < 1200; iChipKey++){
+ if(mult->TestFastOrFiredChips(iChipKey)){
+ // here you check if the FastOr bit is 1 or 0
+ const Int_t iseta = CheckChipEta(iChipKey, etacut, vtxPos, hitMapSPDinner,
+ hitMapSPDouter);
+ if(iseta==0)
+ continue;
+
+ if(iChipKey<400) {
+ ninner++; // here you count the FastOr bits in the inner layer
+ if(iseta>0)
+ ipA ++;
+ else
+ ipC ++;
+ }
+ else {
+ nouter++; // here you count the FastOr bits in the outer layer
+ if(iseta>0)
+ opA ++;
+ else
+ opC ++;
+ }
+ }
+ }
+
+ ctr[kInnerPixel]= ninner;
+ ctr[kOuterPixel]= nouter;
+ ctr[kIPA]= ipA;
+ ctr[kIPC]= ipC;
+ ctr[kOPA]= opA;
+ ctr[kOPC]= opC;
+
+ return;
+}
+
+
+//--------------------------------------------------------------------------
+Int_t AliCDMesonUtils::GetFastORmultiplicity(const AliESDEvent* ESDEvent)
+{
+ // determine the number of fired fastOR chips in both layers within
+ // -0.9 < eta < 0.9
+ //
+
+ const AliMultiplicity *mult = ESDEvent->GetMultiplicity();
+
+ // position of the primary vertex
+ Double_t tmp[3] = { 0., 0., 0. };
+ ESDEvent->GetPrimaryVertex()->GetXYZ(tmp);
+ Double_t vtxPos[3] = { tmp[0], tmp[1], tmp[2] };
+
+ Int_t multiplicity = 0;
+
+ for (Int_t iChipKey=0; iChipKey < 1200; iChipKey++) {
+ if(mult->TestFastOrFiredChips(iChipKey)){
+ // here you check if the FastOr bit is 1 or 0
+ const Int_t iseta = CheckChipEta(iChipKey, "[0.9]", vtxPos, 0x0, 0x0);
+ if(iseta==0)
+ continue;
+ else
+ ++multiplicity;
+ }
+ }
+
+ return multiplicity;
+}
+
+
+//--------------------------------------------------------------------------
+void AliCDMesonUtils::FillSPDtrkltMap(const AliVEvent* event,
+ TH2 *hitMapSPDtrklt)
+{
+ //
+ // fill eta phi map of SPD tracklets
+ //
+
+ if (hitMapSPDtrklt) {
+ if (TString(event->ClassName()) == "AliESDEvent") {
+ const AliMultiplicity *mult = ((AliESDEvent*)event)->GetMultiplicity();
+ for (Int_t iTracklet = 0; iTracklet < mult->GetNumberOfTracklets();
+ iTracklet++) {
+ Double_t eta = mult->GetEta(iTracklet);
+ Double_t phi = mult->GetPhi(iTracklet);
+ hitMapSPDtrklt->Fill(eta, phi);
+ }
+ }
+ else if (TString(event->ClassName()) == "AliAODEvent") {
+ const AliAODTracklets* mult = ((AliAODEvent*)event)->GetTracklets();
+ for (Int_t iTracklet = 0; iTracklet < mult->GetNumberOfTracklets();
+ iTracklet++) {
+ Double_t eta = -TMath::Log(TMath::Tan(mult->GetTheta(iTracklet)/2.));
+ Double_t phi = mult->GetPhi(iTracklet);
+ hitMapSPDtrklt->Fill(eta, phi);
+ }
+ }
+ }
+}
+
+
+//==========================================================================
+Int_t AliCDMesonUtils::GetPID(AliPIDResponse *pid, const AliVTrack *trk,
+ const Int_t mode /* = 0 */)
+{
+ // determines PID for ESDs and AODs
+ //
+ //
+
+ if (!pid) return AliCDMesonBase::kBinPIDUnknown; // no pid available
+
+ Double_t tpcpion = -999., tpckaon = -999., tpcproton = -999.,
+ tpcelectron = -999.;
+ Double_t tofpion = -999., tofkaon = -999., tofproton = -999.,
+ tofelectron = -999.;
+ Double_t itspion = -999., itskaon = -999., itsproton = -999.,
+ itselectron = -999.;
+
+
+ // check whether the track was pure ITS standalone track (has not left TPC)
+ const Bool_t kits = !(trk->GetStatus() & AliESDtrack::kTPCout);
+
+ if (kits) { // do ITS pid
+ const Double_t nin = 3.;
+
+ itspion = pid->NumberOfSigmasITS( trk, AliPID::kPion );
+ itskaon = pid->NumberOfSigmasITS( trk, AliPID::kKaon );
+ itsproton = pid->NumberOfSigmasITS( trk, AliPID::kProton );
+ itselectron = pid->NumberOfSigmasITS( trk, AliPID::kElectron );
+
+ if(fabs(itspion) < nin) // inclusion only
+ return AliCDMesonBase::kBinPion;
+ else if(fabs(itskaon) < nin)
+ return AliCDMesonBase::kBinKaon;
+ else if(fabs(itsproton) < nin)
+ return AliCDMesonBase::kBinProton;
+ else if(fabs(itselectron) < nin)
+ return AliCDMesonBase::kBinElectron;
+ else{
+ return AliCDMesonBase::kBinPIDUnknown;
+ }
+ }
+
+ tpcpion = pid->NumberOfSigmasTPC( trk, AliPID::kPion );
+ tpckaon = pid->NumberOfSigmasTPC( trk, AliPID::kKaon );
+ tpcproton = pid->NumberOfSigmasTPC( trk, AliPID::kProton );
+ tpcelectron = pid->NumberOfSigmasTPC( trk, AliPID::kElectron );
+
+ // derive information, whether tof pid is available
+ const Bool_t ka = !(trk->GetStatus() & AliESDtrack::kTOFmismatch);
+ const Bool_t kb = (trk->GetStatus() & AliESDtrack::kTOFpid);
+ const Bool_t ktof = ka && kb;
+
+ // retrieve TOF information if available
+ if(ktof){
+ tofpion = pid->NumberOfSigmasTOF(trk, AliPID::kPion);
+ tofkaon = pid->NumberOfSigmasTOF(trk, AliPID::kKaon);
+ tofproton = pid->NumberOfSigmasTOF(trk, AliPID::kProton);
+ tofelectron = pid->NumberOfSigmasTOF(trk, AliPID::kElectron);
+ }
+
+ if (mode == 0) {
+ // 2010 cuts for PID
+ const Double_t nin = 3.;
+
+ // inclusion + exclusion (either TPC or TOF)
+ if((fabs(tpcpion) < nin && fabs(tpckaon) > nin && fabs(tpcproton) > nin
+ && fabs(tpcelectron) > nin) ||
+ (fabs(tofpion) < nin && fabs(tofkaon) > nin && fabs(tofproton) > nin
+ && fabs(tofelectron) > nin))
+ return AliCDMesonBase::kBinPionE;
+ else if((fabs(tpcpion) > nin && fabs(tpckaon) < nin && fabs(tpcproton) > nin
+ && fabs(tpcelectron) > nin) ||
+ (fabs(tofpion) > nin && fabs(tofkaon) < nin && fabs(tofproton) > nin
+ && fabs(tofelectron) > nin))
+ return AliCDMesonBase::kBinKaonE;
+ else if((fabs(tpcpion) > nin && fabs(tpckaon) > nin && fabs(tpcproton) < nin
+ && fabs(tpcelectron) > nin) ||
+ (fabs(tofpion) > nin && fabs(tofkaon) > nin && fabs(tofproton) < nin
+ && fabs(tofelectron) > nin))
+ return AliCDMesonBase::kBinProtonE;
+ else if((fabs(tpcpion) > nin && fabs(tpckaon) > nin && fabs(tpcproton) > nin
+ && fabs(tpcelectron) < nin) ||
+ (fabs(tofpion) > nin && fabs(tofkaon) > nin && fabs(tofproton) > nin
+ && fabs(tofelectron) < nin))
+ return AliCDMesonBase::kBinElectronE;
+ else if(fabs(tpcpion) < nin && fabs(tofpion) < nin) // inclusion (TPC + TOF)
+ return AliCDMesonBase::kBinPion;
+ else if(fabs(tpckaon) < nin && fabs(tofkaon) < nin)
+ return AliCDMesonBase::kBinKaon;
+ else if(fabs(tpcproton) < nin && fabs(tofproton) < nin)
+ return AliCDMesonBase::kBinProton;
+ else if(fabs(tpcelectron) < nin && fabs(tofelectron) < nin)
+ return AliCDMesonBase::kBinElectron;
+ else{
+ return AliCDMesonBase::kBinPIDUnknown;
+ }
+ }
+ else if (mode == 1) {
+ // 2011 cuts for PID in LHC11f
+ // TPC: [-3,5] sigma (pion)
+ // TOF: 3 sigma for all,
+ // only Pion is tuned!
+ // ONLY INCLUSION CUTS NO EXCLUSION
+ const Double_t nin = 3.;
+
+ if(tpcpion < 4. && tpcpion > -2. && fabs(tofpion) < -3. )
+ return AliCDMesonBase::kBinPion;
+ else if(fabs(tpckaon) < nin && fabs(tofkaon) < nin)
+ return AliCDMesonBase::kBinKaon;
+ else if(fabs(tpcproton) < nin && fabs(tofproton) < nin)
+ return AliCDMesonBase::kBinProton;
+ else if(fabs(tpcelectron) < nin && fabs(tofelectron) < nin)
+ return AliCDMesonBase::kBinElectron;
+ else{
+ return AliCDMesonBase::kBinPIDUnknown;
+ }
+ }
+ return AliCDMesonBase::kBinPIDUnknown;
+}
+
+
+//==========================================================================
+Int_t AliCDMesonUtils::GetPID(const Int_t pdgCode)
+{
+ //
+ // determine particle type based on PDG code
+ //
+
+ if (TMath::Abs(pdgCode) == 211) return AliCDMesonBase::kBinPionE;
+ else if (TMath::Abs(pdgCode) == 321) return AliCDMesonBase::kBinKaonE;
+ else if (TMath::Abs(pdgCode) == 2212) return AliCDMesonBase::kBinProtonE;
+ else if (TMath::Abs(pdgCode) == 11) return AliCDMesonBase::kBinElectronE;
+ else return AliCDMesonBase::kBinPIDUnknown;
+}
+
+
+//------------------------------------------------------------------------------
+Int_t AliCDMesonUtils::CombinePID(const Int_t pid[])
+{
+ //
+ // combine the PID result
+ //
+
+ // determine return value
+ if (pid[0] == pid[1]) { // same result for both tracks
+ return pid[0];
+ }
+ // one track identified with exclusion the other only without
+ else if ((pid[0] == AliCDMesonBase::kBinPionE &&
+ pid[1] == AliCDMesonBase::kBinPion) ||
+ (pid[1] == AliCDMesonBase::kBinPionE &&
+ pid[0] == AliCDMesonBase::kBinPion)) {
+ return AliCDMesonBase::kBinPion;
+ }
+ else if ((pid[0] == AliCDMesonBase::kBinKaonE &&
+ pid[1] == AliCDMesonBase::kBinKaon) ||
+ (pid[1] == AliCDMesonBase::kBinKaonE &&
+ pid[0] == AliCDMesonBase::kBinKaon)) {
+ return AliCDMesonBase::kBinKaon;
+ }
+ else if ((pid[0] == AliCDMesonBase::kBinProtonE &&
+ pid[1] == AliCDMesonBase::kBinProton) ||
+ (pid[1] == AliCDMesonBase::kBinProtonE &&
+ pid[0] == AliCDMesonBase::kBinProton)) {
+ return AliCDMesonBase::kBinProton;
+ }
+ else if ((pid[0] == AliCDMesonBase::kBinElectronE &&
+ pid[1] == AliCDMesonBase::kBinElectron) ||
+ (pid[1] == AliCDMesonBase::kBinElectronE &&
+ pid[0] == AliCDMesonBase::kBinElectron)) {
+ return AliCDMesonBase::kBinElectron;
+ }
+ // one track identified and one not
+ else if (((pid[0] == AliCDMesonBase::kBinPionE ||
+ pid[0] == AliCDMesonBase::kBinPion) &&
+ pid[1] == AliCDMesonBase::kBinPIDUnknown) ||
+ (pid[1] == AliCDMesonBase::kBinPion &&
+ pid[0] == AliCDMesonBase::kBinPIDUnknown)) {
+ return AliCDMesonBase::kBinSinglePion;
+ }
+ else if (((pid[0] == AliCDMesonBase::kBinKaonE ||
+ pid[0] == AliCDMesonBase::kBinKaon)&&
+ pid[1] == AliCDMesonBase::kBinPIDUnknown) ||
+ (pid[1] == AliCDMesonBase::kBinKaon &&
+ pid[0] == AliCDMesonBase::kBinPIDUnknown)) {
+ return AliCDMesonBase::kBinSingleKaon;
+ }
+ else if (((pid[0] == AliCDMesonBase::kBinProtonE ||
+ pid[0] == AliCDMesonBase::kBinProton) &&
+ pid[1] == AliCDMesonBase::kBinPIDUnknown) ||
+ (pid[1] == AliCDMesonBase::kBinProton &&
+ pid[0] == AliCDMesonBase::kBinPIDUnknown)) {
+ return AliCDMesonBase::kBinSingleProton;
+ }
+ else if (((pid[0] == AliCDMesonBase::kBinElectronE ||
+ pid[0] == AliCDMesonBase::kBinElectron) &&
+ pid[1] == AliCDMesonBase::kBinPIDUnknown) ||
+ (pid[1] == AliCDMesonBase::kBinElectron &&
+ pid[0] == AliCDMesonBase::kBinPIDUnknown)) {
+ return AliCDMesonBase::kBinSingleElectron;
+ }
+ else
+ return AliCDMesonBase::kBinPIDUnknown;
+}
+
+
+//------------------------------------------------------------------------------
+TLorentzVector AliCDMesonUtils::GetKinematics(const Double_t *pa,
+ const Double_t *pb,
+ const Double_t ma,
+ const Double_t mb, Double_t& cts)
+{
+ //
+ //get kinematics, cts = cos(theta_{#star})
+ //
+
+ TLorentzVector va, vb;
+ va.SetXYZM(pa[0], pa[1], pa[2], ma);
+ vb.SetXYZM(pb[0], pb[1], pb[2], mb);
+ const TLorentzVector sumv = va+vb;
+
+ const TVector3 bv = -sumv.BoostVector();
+
+ va.Boost(bv);
+ vb.Boost(bv);
+
+ // 3-vectors in the restframe of the mother particle
+ const TVector3 pra = va.Vect();
+ const TVector3 prb = vb.Vect();
+ const TVector3 diff = pra - prb; // their difference
+
+ cts = (diff.Mag2()-pra.Mag2()-prb.Mag2()) / (-2.*pra.Mag()*prb.Mag());
+ // cosine theta star, calculated according to the law-of-cosines
+
+ return sumv;
+}
+
+
+//------------------------------------------------------------------------------
+Double_t AliCDMesonUtils::GetOA(const Double_t *pa, const Double_t *pb)
+{
+ //
+ //cosOpeningAngle
+ //
+
+ TVector3 va, vb;
+ va.SetXYZ(pa[0], pa[1], pa[2]);
+ vb.SetXYZ(pb[0], pb[1], pb[2]);
+
+ const TVector3 ua = va.Unit();
+ const TVector3 ub = vb.Unit();
+
+ return ua.Dot(ub);
+}
--- /dev/null
+/*************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+//
+// AliCDMesonUtils
+// for
+// AliAnalysisTaskCDMeson
+//
+//
+// Author:
+// Xianguo Lu <lu@physi.uni-heidelberg.de>
+// continued by
+// Felix Reidt <Felix.Reidt@cern.ch>
+
+#ifndef ALICDMESONUTILS_H
+#define ALICDMESONUTILS_H
+
+#include <THnSparse.h> // forward declaration not possible
+
+
+#define STD_ALIROOT
+// if this is defined the code should run with a standard aliroot
+
+
+class TH1;
+class TH2;
+class TLorentzVector;
+class TParticle;
+class TVector3;
+class TObjArray;
+
+class AliVEvent;
+class AliESDEvent;
+class AliAODEvent;
+class AliESDtrack;
+class AliESDtrackCuts;
+class AliVTrack;
+class AliPIDResponse;
+
+class AliCDMesonTracks;
+
+class AliCDMesonUtils
+{
+public:
+ enum{
+ kInnerPixel = 0,
+ kOuterPixel,
+ kIPA,
+ kIPC,
+ kOPA,
+ kOPC
+ };
+
+ // ESD only
+ //---------
+
+ // cuts for ESD analysis
+ static Bool_t CutEvent(const AliESDEvent *ESDEvent, TH1 *hspd, TH1 *hpriv,
+ TH1* hpriVtxPos, TH1* hpriVtxDist, TH2 *hfo,
+ TH1* hfochans, Int_t &kfo, Int_t &nip, Int_t &nop,
+ TH1* hpriVtxX, TH1* hpriVtxY, TH1* hpriVtxZ);
+
+ // V0-AND/OR
+ static Bool_t V0AND(const AliESDEvent *ESDEvent) {
+ // were both V0 detectors A/C firing?
+ // (4 = AliCDMesonBase::kBinNG)
+ return (GetV0(ESDEvent) == 4);
+ }
+
+ static Bool_t V0OR(const AliESDEvent *ESDEvent) {
+ // did one of the V0 detectors fire?
+ // (1-4 => single or both detectors active)
+ return (GetV0(ESDEvent) > 0);
+ }
+
+ // VZERO study
+ static void DoVZEROStudy(const AliESDEvent *ESDEvent, TObjArray* hists,
+ const Int_t run);
+
+ // gap determination
+ static Int_t GetGapConfig(const AliESDEvent *ESDEvent, TH2 *hitMapSPDinner,
+ TH2 *hitMapSPDouter, TH2 *hitMapSPDtrklt,
+ TH2 *hitMapFMDa, TH2 *hitMapFMDc,
+ TH1 **fmdSums, TH2 *TPCGapDCAaSide,
+ TH2 *TPCGapDCAcSide);
+
+ static Int_t GetFastORmultiplicity(const AliESDEvent *ESDEvent);
+
+ static void GetSPDTrackletMult(const AliESDEvent *ESDEvent, Int_t& sum,
+ Int_t& forwardA, Int_t& forwardC,
+ Int_t& central);
+
+ static void SPDLoadGeom(const Int_t run); // only needed for ESDs, not in AODs
+
+ // functions needed for the empty event study (ESD only)
+ static Int_t GetEventType(const AliESDEvent *ESDEvent);
+ static void GetMultFMD(const AliESDEvent *ESDEvent, Int_t& fmdA,
+ Int_t& fmdC, Float_t *fmdSums);
+ static void GetMultSPD(const AliESDEvent *ESDEvent, Int_t& spdIA,
+ Int_t& spdIC, Int_t& spdOA, Int_t& spdOC);
+
+
+ // AOD only
+ //---------
+ static Bool_t CutEvent(const AliAODEvent *AODEvent, TH1 *hpriv, TH1* hpriVtxX,
+ TH1* hpriVtxY, TH1* hpriVtxZ, TH1* hpriVtxPos,
+ TH1* hpriVtxDist);
+
+
+ // independent
+ //------------
+ static void FillEtaPhiMap(const AliVEvent *event,
+ const AliCDMesonTracks* tracks, TH2 *map,
+ TH2 *map_c);
+ static void SwapTrack(const AliVTrack *trks[]);
+ static Int_t GetCombCh(const AliVTrack *trks[]);
+ static Int_t GetCombPID(AliPIDResponse *pid, const AliVTrack *trks[],
+ const Int_t mode, TH2 *comb2trkPID = 0x0);
+ static Int_t GetCombPID(const TParticle* particles[], TH2 *comb2trkPID = 0x0);
+ static void GetMassPtCtsOA(const Int_t pid, const TVector3* momenta[],
+ Double_t& mass, Double_t& pt, Double_t& cts,
+ Double_t& oa);
+ static void GetPWAinfo(const Int_t pid, const AliVTrack *trks[],
+ Float_t& theta, Float_t& phi, Float_t& mass,
+ Float_t momentum[]);
+ static void FillSPDtrkltMap(const AliVEvent *event, TH2 *hitMapSPDtrklt);
+
+private:
+ // ESD only
+ //---------
+
+ // Gap determination functions
+ static Int_t GetV0(const AliESDEvent *ESDEvent);
+ static Int_t GetFMD(const AliESDEvent *ESDEvent, TH2 *hitMapFMDa,
+ TH2 *hitMapFMDc, TH1** fmdSums);
+ static Int_t GetSPD(const AliESDEvent *ESDEvent, TH2 *hitMapSPDinner,
+ TH2 *hitMapSPDouter, TH2 *hitMapSPDtrklt);
+ static Int_t GetTPC(const AliESDEvent *ESDEvent, TH2 *TPCGapDCAaSide,
+ TH2 *TPCGapDCAcSide);
+ static Int_t GetZDC(const AliESDEvent *ESDEvent); // not used so far
+
+#ifdef STD_ALIROOT
+ // helpers for the FMD gap determination
+ static Int_t FMDHitCombinations(const AliESDEvent* ESDEvent, Int_t side);
+#endif
+
+ // helpers for the SPD gap determination
+ static Bool_t SPDLoc2Glo(const Int_t id, const Double_t *loc, Double_t *glo);
+ static Int_t CheckChipEta(const Int_t chipKey, const TString scut,
+ const Double_t vtxPos[], TH2 *hitMapSPDinner,
+ TH2 *hitMapSPDouter);
+ static void GetNFO(const AliESDEvent *ESDEvent, const TString etacut,
+ Int_t ctr[], TH2 *hitMapSPDinner, TH2 *hitMapSPDouter);
+
+ // AOD only
+ //---------
+
+
+ // independent
+ //----------
+ static Int_t GetPID(AliPIDResponse *pid, const AliVTrack *trk,
+ const Int_t mode = 0);
+ static Int_t GetPID(const Int_t pdgCode);
+ static Int_t CombinePID(const Int_t pid[]);
+
+ static TLorentzVector GetKinematics(const Double_t *pa, const Double_t *pb,
+ const Double_t ma, const Double_t mb,
+ Double_t & cts);
+ static Double_t GetOA(const Double_t *pa, const Double_t *pb);
+};
+
+#endif
--- /dev/null
+// runMeson.C
+//
+// run macro for central diffractive meson analysis
+//
+// Author: Felix Reidt <felix.reidt@cern.ch>
+//
+// based on:
+// ---
+// Template run macro for AliBasicTask.cxx/.h with example layout of
+// physics selections and options, in macro and task.
+//
+// Author: Arvinder Palaha
+//
+class AliAnalysisGrid;
+
+//______________________________________________________________________________
+void runMeson(
+ const char* runtype = "grid", // local, proof or grid
+ const char *gridmode = "full", // Set the run mode (can be "full", "test", "offline", "submit" or "terminate"). Full & Test work for proof
+ const bool bMCphyssel = 0, // 1 = looking at MC truth or reconstructed, 0 = looking at real data
+ const Long64_t nentries = 2000, // for local and proof mode, ignored in grid mode. Set to 1234567890 for all events.
+ const Long64_t firstentry = 0, // for local and proof mode, ignored in grid mode
+ const char *proofdataset = "/alice/data/LHC10c_000120821_p1", // path to dataset on proof cluster, for proof analysis
+ const char *proofcluster = "alice-caf.cern.ch", // which proof cluster to use in proof mode
+ const char *taskname = "CDMeson" // sets name of grid generated macros
+ )
+{
+ // check run type
+ if(runtype != "local" && runtype != "proof" && runtype != "grid"){
+ Printf("\n\tIncorrect run option, check first argument of run macro");
+ Printf("\tint runtype = local, proof or grid\n");
+ return;
+ }
+ Printf("%s analysis chosen",runtype);
+
+ // load libraries (to be optimized)
+ gSystem->Load("libCore.so");
+ gSystem->Load("libTree.so");
+ gSystem->Load("libPhysics");
+ gSystem->Load("libMinuit");
+ gSystem->Load("libProof");
+ gSystem->Load("libmicrocern");
+ gSystem->Load("liblhapdf");
+ gSystem->Load("libpythia6");
+ gSystem->Load("libEG");
+ gSystem->Load("libGeom");
+ gSystem->Load("libVMC");
+ gSystem->Load("libEGPythia6");
+ gSystem->Load("libSTEERBase");
+ gSystem->Load("libESD");
+ gSystem->Load("libCDB");
+ gSystem->Load("libRAWDatabase");
+ gSystem->Load("libRAWDatarec");
+ gSystem->Load("libAOD");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libANALYSISalice");
+ gSystem->Load("libSTEER");
+ gSystem->Load("libTENDER");
+ gSystem->Load("libRAWDatasim");
+ gSystem->Load("libFASTSIM");
+ gSystem->Load("libEVGEN");
+ gSystem->Load("libAliPythia6");
+ gSystem->Load("libSTAT");
+ gSystem->Load("libhijing");
+ gSystem->Load("libTHijing");
+ gSystem->Load("libSTRUCT");
+ gSystem->Load("libPHOSUtils");
+ gSystem->Load("libPHOSbase");
+ gSystem->Load("libPHOSsim");
+ gSystem->Load("libPHOSrec");
+ gSystem->Load("libMUONcore");
+ gSystem->Load("libMUONmapping");
+ gSystem->Load("libMUONgeometry");
+ gSystem->Load("libMUONcalib");
+ gSystem->Load("libMUONraw");
+ gSystem->Load("libMUONtrigger");
+ gSystem->Load("libMUONbase");
+ gSystem->Load("libMUONsim");
+ gSystem->Load("libMUONrec");
+ gSystem->Load("libMUONevaluation");
+ gSystem->Load("libFMDbase");
+ gSystem->Load("libFMDsim");
+ gSystem->Load("libFMDrec");
+ gSystem->Load("libPMDbase");
+ gSystem->Load("libPMDsim");
+ gSystem->Load("libPMDrec");
+ gSystem->Load("libHMPIDbase");
+ gSystem->Load("libHMPIDsim");
+ gSystem->Load("libHMPIDrec");
+ gSystem->Load("libT0base");
+ gSystem->Load("libT0sim");
+ gSystem->Load("libT0rec");
+ gSystem->Load("libZDCbase");
+ gSystem->Load("libZDCsim");
+ gSystem->Load("libZDCrec");
+ gSystem->Load("libACORDEbase");
+ gSystem->Load("libACORDErec");
+ gSystem->Load("libACORDEsim");
+ gSystem->Load("libVZERObase");
+ gSystem->Load("libVZEROrec");
+ gSystem->Load("libVZEROsim");
+ gSystem->Load("libEMCALraw");
+ gSystem->Load("libEMCALUtils");
+ gSystem->Load("libEMCALbase");
+ gSystem->Load("libEMCALsim");
+ gSystem->Load("libEMCALrec");
+ gSystem->Load("libTPCbase");
+ gSystem->Load("libTPCrec");
+ gSystem->Load("libTPCsim");
+ gSystem->Load("libTPCfast");
+ gSystem->Load("libITSbase");
+ gSystem->Load("libITSsim");
+ gSystem->Load("libITSrec");
+ gSystem->Load("libTRDbase");
+ gSystem->Load("libTRDsim");
+ gSystem->Load("libTRDrec");
+ gSystem->Load("libTOFbase");
+ gSystem->Load("libTOFrec");
+ gSystem->Load("libTOFsim");
+ gSystem->Load("libHLTbase");
+ gSystem->Load("libHLTinterface");
+ gSystem->Load("libHLTsim");
+ gSystem->Load("libHLTrec");
+ gSystem->Load("libPWGPP");
+
+ // add aliroot include path
+ gROOT->ProcessLine(Form(".include %s/include",
+ gSystem->ExpandPathName("$ALICE_ROOT")));
+ gROOT->ProcessLine(Form(".include $ALICE_ROOT/include",
+ gSystem->ExpandPathName("$ALICE_ROOT")));
+ gROOT->ProcessLine(Form(".include $ALICE_ROOT/ITS",
+ gSystem->ExpandPathName("$ALICE_ROOT")));
+ gROOT->ProcessLine(Form(".include $ALICE_ROOT/PWGPP/ITS",
+ gSystem->ExpandPathName("$ALICE_ROOT")));
+ gROOT->ProcessLine(Form(".include $ALICE_ROOT/VZERO",
+ gSystem->ExpandPathName("$ALICE_ROOT")));
+
+ gROOT->SetStyle("Plain");
+
+ // create the alien handler and attach it to the manager
+ AliAnalysisGrid *plugin =
+ CreateAlienHandler(taskname, gridmode, proofcluster, proofdataset);
+
+ // analysis manager
+ AliAnalysisManager* mgr = new AliAnalysisManager("CDMeson-Manager");
+ mgr->SetGridHandler(plugin);
+
+ AliESDInputHandler* esdH = new AliESDInputHandler();
+ mgr->SetInputEventHandler(esdH);
+
+ // === Physics Selection Task ===
+ gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
+ AliPhysicsSelectionTask *physSelTask = AddTaskPhysicsSelection(bMCphyssel);
+ if(!physSelTask) { Printf("no physSelTask"); return; }
+
+ // === Add PID Response Task ===
+ gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
+ AliAnalysisTaskPIDResponse *pidResponseTask = AddTaskPIDResponse(bMCphyssel);
+ if(!pidResponseTask) { Printf("no pidResponseTask"); return; }
+
+
+ // === create user task ===
+ gROOT->LoadMacro("AliCDMesonBase.cxx+g");
+ gROOT->LoadMacro("AliCDMesonTracks.cxx+g");
+ gROOT->LoadMacro("AliCDMesonUtils.cxx+g");
+ gROOT->LoadMacro("AliAnalysisTaskCDMeson.cxx+g");
+
+ Long_t taskConfig = AliCDMesonBase::kBitConfigurationSet; // has to be set
+ taskConfig |= AliCDMesonBase::kBitEtaPhiMaps;
+ taskConfig |= AliCDMesonBase::kBitEtaPhiMapsWithCuts;
+ taskConfig |= AliCDMesonBase::kBitStatsFlow;
+ taskConfig |= AliCDMesonBase::kBitMultPerGapHists;
+ taskConfig |= AliCDMesonBase::kBitRmMultPerGapHists;
+ taskConfig |= AliCDMesonBase::kBitTHnMother;
+ taskConfig |= AliCDMesonBase::kBitFastORStudy;
+ taskConfig |= AliCDMesonBase::kBitHitMapSPD;
+ taskConfig |= AliCDMesonBase::kBitHitMapFMD;
+ taskConfig |= AliCDMesonBase::kBitVtxStudies;
+ taskConfig |= AliCDMesonBase::kBitPIDStudy;
+ taskConfig |= AliCDMesonBase::kBitFMDsum;
+ taskConfig |= AliCDMesonBase::kBitSoftTracks;
+ taskConfig |= AliCDMesonBase::kBitPWAtree;
+ taskConfig |= AliCDMesonBase::kBitMultStudy;
+ //taskConfig |= AliCDMesonBase::kBitReadPreprocessedGap;
+ taskConfig |= AliCDMesonBase::kBitVZEROStudy;
+ taskConfig |= AliCDMesonBase::kBitTPCGapStudy;
+ if(bMCphyssel) taskConfig |= AliCDMesonBase::kBitMCProcess; // only for MC data
+ if(bMCphyssel) taskConfig |= AliCDMesonBase::kBitTHnMC;
+ if(bMCphyssel) taskConfig |= AliCDMesonBase::kBitMultResponseMC;
+ //taskConfig |= AliCDMesonBase::kBitReduceGapEvents;
+ taskConfig |= AliCDMesonBase::kBitConfigurationVersion;
+
+ printf("taskConfig=0x%x\n", taskConfig);
+
+ AliAnalysisTaskSE* task = new AliAnalysisTaskCDMeson(taskname, taskConfig);
+ task->SelectCollisionCandidates(AliVEvent::kMB);
+ mgr->AddTask(task);
+
+ // INPUT ---------------------------------------------------------------------
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+ mgr->ConnectInput(task, 0, cinput );
+
+ // OUTPUT --------------------------------------------------------------------
+ // output filename
+ Char_t foutname[100];
+ Char_t foutnamePWA[100];
+ sprintf(foutname,"freidt_%s.root",taskname);
+ sprintf(foutnamePWA,"freidt_%sPWA.root",taskname);
+
+ // output containers
+ // in AnalysisTaskSE, slot 0 reserved, must start from 1
+ // since the number of output containers is flexible, the slot is assigned
+ // dynamically
+ AliAnalysisDataContainer* outlist[6] = { 0x0, 0x0, 0x0, 0x0, 0x0, 0x0 };
+ Int_t nOutputs = 0;
+
+ outlist[nOutputs] =
+ mgr->CreateContainer("freidt_Hist", TList::Class(),
+ AliAnalysisManager::kOutputContainer,foutname);
+ nOutputs++;
+ mgr->ConnectOutput(task, nOutputs, outlist[nOutputs-1]);
+
+ if (!taskConfig
+ || ((taskConfig & AliCDMesonBase::kBitSoftTracks)
+ && (taskConfig & AliCDMesonBase::kBitTHnMother))) {
+ outlist[nOutputs] =
+ mgr->CreateContainer("freidt_ThnMother", THnSparse::Class(),
+ AliAnalysisManager::kOutputContainer,foutname);
+ nOutputs++;
+ mgr->ConnectOutput(task, nOutputs, outlist[nOutputs-1]);
+ }
+
+ if (!taskConfig || (taskConfig & AliCDMesonBase::kBitSoftTracks)) {
+ outlist[nOutputs] =
+ mgr->CreateContainer("freidt_ThnMotherSoft", THnSparse::Class(),
+ AliAnalysisManager::kOutputContainer,foutname);
+ nOutputs++;
+ mgr->ConnectOutput(task, nOutputs, outlist[nOutputs-1]);
+ }
+
+ if (!taskConfig || (taskConfig & AliCDMesonBase::kBitMultStudy)) {
+ outlist[nOutputs] =
+ mgr->CreateContainer("freidt_ThnMultiplicity", THnSparse::Class(),
+ AliAnalysisManager::kOutputContainer,foutname);
+ nOutputs++;
+ mgr->ConnectOutput(task, nOutputs, outlist[nOutputs-1]);
+ }
+
+ if (!taskConfig || (taskConfig & AliCDMesonBase::kBitTHnMC)) {
+ outlist[nOutputs] =
+ mgr->CreateContainer("freidt_ThnMotherMC", THnSparse::Class(),
+ AliAnalysisManager::kOutputContainer,foutname);
+ nOutputs++;
+ mgr->ConnectOutput(task, nOutputs, outlist[nOutputs-1]);
+ }
+ if (!taskConfig || (taskConfig & AliCDMesonBase::kBitPWAtree)) {
+ outlist[nOutputs] =
+ mgr->CreateContainer("freidt_PWA", TTree::Class(),
+ AliAnalysisManager::kOutputContainer,foutnamePWA);
+ nOutputs++;
+ mgr->ConnectOutput(task, nOutputs, outlist[nOutputs-1]);
+ }
+
+ // enable debug printouts
+ mgr->SetDebugLevel(2);
+ //mgr->SetNSysInfo(100);
+ if (!mgr->InitAnalysis()) return;
+ mgr->PrintStatus();
+
+ // start analysis
+ Printf("Starting Analysis....");
+ mgr->StartAnalysis(runtype); //,nentries,firstentry);
+}
+
+//______________________________________________________________________________
+AliAnalysisGrid* CreateAlienHandler(const char *taskname,
+ const char *gridmode,
+ const char *proofcluster,
+ const char *proofdataset)
+{
+ AliAnalysisAlien *plugin = new AliAnalysisAlien();
+ // Set the run mode (can be "full", "test", "offline", "submit" or "terminate")
+ plugin->SetOverwriteMode();
+ plugin->SetRunMode(gridmode);
+
+ plugin->SetMergeViaJDL(kTRUE);
+
+ // Set versions of used packages
+ plugin->SetAPIVersion("V1.1x");
+ plugin->SetROOTVersion("v5-34-01-1");
+ plugin->SetAliROOTVersion("v5-03-50-AN");
+
+ // Declare input data to be processed.
+ //plugin->SetCheckCopy(kFALSE);
+
+ // Method 1: Create automatically XML collections using alien 'find' command.
+ // Define production directory LFN
+ plugin->SetGridDataDir("/alice/sim/2012/LHC12d4a");
+ // On real reconstructed data:
+ // plugin->SetGridDataDir("/alice/data/2009/LHC09d");
+ // Set data search pattern
+ plugin->SetDataPattern("*AliESDs.root"); // THIS CHOOSES ALL PASSES
+ // Data pattern for reconstructed data
+ //plugin->SetDataPattern("*ESDs/pass2/*ESDs.root"); // CHECK LATEST PASS OF DATA SET IN ALIENSH
+ // plugin->SetDataPattern("ESDs/pass2/AOD038/*AliAOD.root"); // CHECK LATEST PASS OF DATA SET IN ALIENSH
+ plugin->SetRunPrefix(""); //000"); // real data
+ // ...then add run numbers to be considered
+
+ //Int_t runlist[15]={117039, 146859, 146858, 146856, 146824, 146817, 146806, 146805, 146804, 146803, 146802, 146801, 146748, 146747, 146746};
+ //for (Int_t ind=0; ind<1; ind++) {
+ // plugin->AddRunNumber(runlist[ind]);
+ //}
+ //plugin->SetRunRange(114917,115322);
+ //plugin->AddRunNumber(117050);
+ plugin->SetRunRange(115393,126408);
+
+ plugin->SetNrunsPerMaster(1);
+
+ // Define alien work directory where all files will be copied. Relative to alien $HOME.
+ plugin->SetGridWorkingDir(taskname);
+
+ // Declare alien output directory. Relative to working directory.
+ plugin->SetGridOutputDir("out"); // In this case will be $HOME/taskname/out
+
+ // Declare the analysis source files names separated by blancs. To be compiled runtime
+ // using ACLiC on the worker nodes.
+ plugin->SetAnalysisSource("AliCDMesonBase.cxx AliCDMesonTracks.cxx AliCDMesonUtils.cxx AliAnalysisTaskCDMeson.cxx");
+
+ // Declare all libraries (other than the default ones for the framework. These will be
+ // loaded by the generated analysis macro. Add all extra files (task .cxx/.h) here.
+ plugin->SetAdditionalLibs("libGui.so libCore.so libTree.so libPhysics.so libMinuit.so libProof.so libmicrocern.so liblhapdf.so libpythia6.so libEG.so libGeom.so libVMC.so libEGPythia6.so libSTEERBase.so libESD.so libRAWDatabase.so libRAWDatarec.so libAOD.so libANALYSIS.so libANALYSISalice.so libCDB.so libSTEER.so libRAWDatasim.so libFASTSIM.so libEVGEN.so libAliPythia6.so libSTAT.so libhijing.so libTHijing.so libSTRUCT.so libPHOSUtils.so libPHOSbase.so libPHOSsim.so libPHOSrec.so libMUONcore.so libMUONmapping.so libMUONgeometry.so libMUONcalib.so libMUONraw.so libMUONtrigger.so libMUONbase.so libMUONsim.so libMUONrec.so libMUONevaluation.so libFMDbase.so libFMDsim.so libFMDrec.so libPMDbase.so libPMDsim.so libPMDrec.so libHMPIDbase.so libHMPIDsim.so libHMPIDrec.so libT0base.so libT0sim.so libT0rec.so libZDCbase.so libZDCsim.so libZDCrec.so libACORDEbase.so libACORDErec.so libACORDEsim.so libVZERObase.so libVZEROrec.so libVZEROsim.so libEMCALraw.so libEMCALUtils.so libEMCALbase.so libEMCALsim.so libEMCALrec.so libTPCbase.so libTPCrec.so libTPCsim.so libTPCfast.so libITSbase.so libITSsim.so libITSrec.so libTRDbase.so libTRDsim.so libTRDrec.so libTOFbase.so libTOFrec.so libTOFsim.so libHLTbase.so libHLTinterface.so libHLTsim.so libHLTrec.so libTENDER.so libTENDERSupplies.so libPWGPP.so AliCDMesonBase.h AliCDMesonBase.cxx AliCDMesonTracks.h AliCDMesonTracks.cxx AliCDMesonUtils.h AliCDMesonUtils.cxx AliAnalysisTaskCDMeson.h AliAnalysisTaskCDMeson.cxx");
+
+ plugin->AddIncludePath("-I$ALICE_ROOT/ITS -I$ALICE_ROOT/PWGPP/ITS");
+
+ // Declare the output file names separated by blancs.
+ // (can be like: file.root or file.root@ALICE::Niham::File)
+ // To only save certain files, use SetDefaultOutputs(kFALSE), and then
+ // SetOutputFiles("list.root other.filename") to choose which files to save
+ plugin->SetDefaultOutputs();
+
+ // Optionally set a name for the generated analysis macro (default MyAnalysis.C)
+ plugin->SetAnalysisMacro("CDMeson.C");
+
+ // Optionally set maximum number of input files/subjob (default 100, put 0 to ignore)
+ plugin->SetSplitMaxInputFileNumber(100);
+
+ // Optionally modify the executable name (default analysis.sh)
+ plugin->SetExecutable("CDMeson.sh");
+
+ // set number of test files to use in "test" mode
+ plugin->SetNtestFiles(1);
+
+ // Optionally resubmit threshold.
+ //plugin->SetMasterResubmitThreshold(90);
+
+ // Optionally set time to live (default 30000 sec)
+ //plugin->SetTTL(30000);
+
+ // Optionally set input format (default xml-single)
+ //plugin->SetInputFormat("xml-single");
+
+ // Optionally modify the name of the generated JDL (default analysis.jdl)
+ plugin->SetJDLName("CDMeson.jdl");
+
+ // Optionally modify job price (default 1)
+ plugin->SetPrice(1);
+
+ // Optionally modify split mode (default 'se')
+ plugin->SetSplitMode("se");
+
+ //----------------------------------------------------------
+ //--- PROOF MODE SPECIFIC SETTINGS ------------
+ //----------------------------------------------------------
+ // Proof cluster
+ plugin->SetProofCluster(proofcluster);
+ // Dataset to be used
+ plugin->SetProofDataSet(proofdataset);
+ // May need to reset proof. Supported modes: 0-no reset, 1-soft, 2-hard
+ plugin->SetProofReset(0);
+ // May limit number of workers
+ plugin->SetNproofWorkers(0);
+ // May limit the number of workers per slave
+ plugin->SetNproofWorkersPerSlave(1);
+ // May use a specific version of root installed in proof
+ plugin->SetRootVersionForProof("current");
+ // May set the aliroot mode. Check http://aaf.cern.ch/node/83
+ plugin->SetAliRootMode("default"); // Loads AF libs by default
+ // May request ClearPackages (individual ClearPackage not supported)
+ plugin->SetClearPackages(kFALSE);
+ // Plugin test mode works only providing a file containing test file locations, used in "local" mode also
+ plugin->SetFileForTestMode("files.txt"); // file should contain path name to a local directory containg *ESDs.root etc
+ // Request connection to alien upon connection to grid
+ plugin->SetProofConnectGrid(kFALSE);
+ // Other PROOF specific parameters
+ plugin->SetProofParameter("PROOF_UseMergers","-1");
+ printf("Using: PROOF_UseMergers : %s\n", plugin->GetProofParameter("PROOF_UseMergers"));
+
+ return plugin;
+}
#pragma link off all classes;
#pragma link off all functions;
+#pragma link C++ class AliAnalysisTaskCDMeson+;
+#pragma link C++ class AliAnalysisTaskCDskimESD+;
#pragma link C++ class AliAnalysisTaskCDex+;
#pragma link C++ class AliCDMesonBaseStripped+;
+#pragma link C++ class AliCDMesonBase+;
#pragma link C++ class AliCDMesonTracks+;
#pragma link C++ class AliCDMesonUtilsStripped+;
+#pragma link C++ class AliCDMesonUtils+;
+
#endif