-#include <vector>\r
-#include "TChain.h"\r
-#include "TList.h"\r
-#include "TCanvas.h"\r
-#include "TLorentzVector.h"\r
-#include "TGraphErrors.h"\r
-#include "TH1F.h"\r
-#include "TH2F.h"\r
-#include "TArrayF.h"\r
-#include "TF1.h"\r
-#include "TRandom.h"\r
-\r
-#include "AliLog.h"\r
-\r
-#include "AliAnalysisTaskSE.h"\r
-#include "AliAnalysisManager.h"\r
-\r
-#include "AliESDVertex.h"\r
-#include "AliESDEvent.h"\r
-#include "AliESDInputHandler.h"\r
-#include "AliAODEvent.h"\r
-#include "AliAODTrack.h"\r
-#include "AliAODInputHandler.h"\r
-#include "AliGenEventHeader.h"\r
-#include "AliGenHijingEventHeader.h"\r
-#include "AliMCEventHandler.h"\r
-#include "AliMCEvent.h"\r
-#include "AliMixInputEventHandler.h"\r
-#include "AliStack.h"\r
-\r
-#include "TH2D.h" \r
-#include "AliTHn.h" \r
-\r
-#include "AliEventPoolManager.h" \r
-\r
-#include "AliAnalysisTaskTriggeredBF.h"\r
-#include "AliBalanceTriggered.h"\r
-\r
-\r
-// Analysis task for the TriggeredBF code\r
-// Authors: Panos.Christakoglou@nikhef.nl, m.weber@cern.ch\r
-\r
-// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+\r
-//\r
-// For the V0 part:\r
-// --> AliAnalysisTaskExtractV0AOD (by david.chinellato@gmail.com)\r
-//\r
-// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+\r
-\r
-using std::cout;\r
-using std::endl;\r
-using std::vector;\r
-\r
-ClassImp(AliAnalysisTaskTriggeredBF)\r
-\r
-//________________________________________________________________________\r
-AliAnalysisTaskTriggeredBF::AliAnalysisTaskTriggeredBF(const char *name) \r
-: AliAnalysisTaskSE(name), \r
- fBalance(0),\r
- fRunShuffling(kFALSE),\r
- fShuffledBalance(0),\r
- fRunMixing(kFALSE),\r
- fMixingTracks(50000),\r
- fMixedBalance(0),\r
- fPoolMgr(0),\r
- fRunV0(kFALSE),\r
- fPIDResponse(0),\r
- fPIDCombined(0),\r
- fList(0),\r
- fListTriggeredBF(0),\r
- fListTriggeredBFS(0),\r
- fListTriggeredBFM(0),\r
- fHistListPIDQA(0),\r
- fHistListV0(0),\r
- fHistEventStats(0),\r
- fHistCentStats(0),\r
- fHistTriggerStats(0),\r
- fHistTrackStats(0),\r
- fHistVx(0),\r
- fHistVy(0),\r
- fHistVz(0),\r
- fHistClus(0),\r
- fHistDCA(0),\r
- fHistChi2(0),\r
- fHistPt(0),\r
- fHistEta(0),\r
- fHistPhi(0),\r
- fHistPhiBefore(0),\r
- fHistPhiAfter(0),\r
- fHistV0M(0),\r
- fHistRefTracks(0),\r
- fHistV0MultiplicityBeforeTrigSel(0),\r
- fHistV0MultiplicityForTrigEvt(0),\r
- fHistV0MultiplicityForSelEvt(0),\r
- fHistV0MultiplicityForSelEvtNoTPCOnly(0),\r
- fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup(0),\r
- fHistMultiplicityBeforeTrigSel(0),\r
- fHistMultiplicityForTrigEvt(0),\r
- fHistMultiplicity(0),\r
- fHistMultiplicityNoTPCOnly(0),\r
- fHistMultiplicityNoTPCOnlyNoPileup(0),\r
- fHistV0InvMassK0(0),\r
- fHistV0InvMassLambda(0),\r
- fHistV0InvMassAntiLambda(0),\r
- fHistV0Armenteros(0),\r
- fHistV0SelInvMassK0(0),\r
- fHistV0SelInvMassLambda(0),\r
- fHistV0SelInvMassAntiLambda(0),\r
- fHistV0SelArmenteros(0),\r
- fCentralityEstimator("V0M"),\r
- fUseCentrality(kFALSE),\r
- fCentralityPercentileMin(0.), \r
- fCentralityPercentileMax(5.),\r
- fImpactParameterMin(0.),\r
- fImpactParameterMax(20.),\r
- fUseMultiplicity(kFALSE),\r
- fNumberOfAcceptedTracksMin(0),\r
- fNumberOfAcceptedTracksMax(10000),\r
- fHistNumberOfAcceptedTracks(0),\r
- fUseOfflineTrigger(kFALSE),\r
- fVxMax(0.3),\r
- fVyMax(0.3),\r
- fVzMax(10.),\r
- nAODtrackCutBit(128),\r
- fPtMin(0.3),\r
- fPtMax(1.5),\r
- fEtaMin(-0.8),\r
- fEtaMax(-0.8),\r
- fDCAxyCut(-1),\r
- fDCAzCut(-1),\r
- fTPCchi2Cut(-1),\r
- fNClustersTPCCut(-1)\r
-{\r
- // Constructor\r
- // Define input and output slots here\r
- // Input slot #0 works with a TChain\r
- DefineInput(0, TChain::Class());\r
- // Output slot #0 writes into a TH1 container\r
- DefineOutput(1, TList::Class());\r
- DefineOutput(2, TList::Class());\r
- DefineOutput(3, TList::Class());\r
- DefineOutput(4, TList::Class());\r
- DefineOutput(5, TList::Class());\r
-}\r
-\r
-//________________________________________________________________________\r
-AliAnalysisTaskTriggeredBF::~AliAnalysisTaskTriggeredBF() {\r
-\r
- // Destructor\r
-\r
-}\r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskTriggeredBF::UserCreateOutputObjects() {\r
- // Create histograms\r
- // Called once\r
-\r
- // global switch disabling the reference \r
- // (to avoid "Replacing existing TH1" if several wagons are created in train)\r
- Bool_t oldStatus = TH1::AddDirectoryStatus();\r
- TH1::AddDirectory(kFALSE);\r
-\r
- if(!fBalance) {\r
- fBalance = new AliBalanceTriggered();\r
- fBalance->SetAnalysisLevel("AOD");\r
- }\r
- if(fRunShuffling) {\r
- if(!fShuffledBalance) {\r
- fShuffledBalance = new AliBalanceTriggered();\r
- fShuffledBalance->SetAnalysisLevel("AOD");\r
- }\r
- }\r
- if(fRunMixing) {\r
- if(!fMixedBalance) {\r
- fMixedBalance = new AliBalanceTriggered();\r
- fMixedBalance->SetAnalysisLevel("AOD");\r
- }\r
- }\r
-\r
- //QA list\r
- fList = new TList();\r
- fList->SetName("listQA");\r
- fList->SetOwner();\r
-\r
- //Balance Function list\r
- fListTriggeredBF = new TList();\r
- fListTriggeredBF->SetName("listTriggeredBF");\r
- fListTriggeredBF->SetOwner();\r
-\r
- if(fRunShuffling) {\r
- fListTriggeredBFS = new TList();\r
- fListTriggeredBFS->SetName("listTriggeredBFShuffled");\r
- fListTriggeredBFS->SetOwner();\r
- }\r
- if(fRunMixing) {\r
- fListTriggeredBFM = new TList();\r
- fListTriggeredBFM->SetName("listTriggeredBFMixed");\r
- fListTriggeredBFM->SetOwner();\r
- }\r
- \r
- \r
- //Event stats.\r
- TString gCutName[4] = {"Total","Offline trigger",\r
- "Vertex","Analyzed"};\r
- fHistEventStats = new TH1F("fHistEventStats",\r
- "Event statistics;;N_{events}",\r
- 4,0.5,4.5);\r
- for(Int_t i = 1; i <= 4; i++)\r
- fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());\r
- fList->Add(fHistEventStats);\r
- \r
- TString gCentName[9] = {"V0M","FMD","TRK","TKL","CL0","CL1","V0MvsFMD","TKLvsV0M","ZEMvsZDC"};\r
- fHistCentStats = new TH2F("fHistCentStats",\r
- "Centrality statistics;;Cent percentile",\r
- 9,-0.5,8.5,220,-5,105);\r
- for(Int_t i = 1; i <= 9; i++)\r
- fHistCentStats->GetXaxis()->SetBinLabel(i,gCentName[i-1].Data());\r
- fList->Add(fHistCentStats);\r
- \r
- fHistTriggerStats = new TH1F("fHistTriggerStats","Trigger statistics;TriggerBit;N_{events}",130,0,130);\r
- fList->Add(fHistTriggerStats);\r
- \r
- fHistTrackStats = new TH1F("fHistTrackStats","Event statistics;TrackFilterBit;N_{events}",130,0,130);\r
- fList->Add(fHistTrackStats);\r
-\r
- fHistNumberOfAcceptedTracks = new TH1F("fHistNumberOfAcceptedTracks",";N_{acc.};Entries",4001,-0.5,4000.5);\r
- fList->Add(fHistNumberOfAcceptedTracks);\r
-\r
- // Vertex distributions\r
- fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);\r
- fList->Add(fHistVx);\r
- fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);\r
- fList->Add(fHistVy);\r
- fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);\r
- fList->Add(fHistVz);\r
-\r
- // QA histograms\r
- fHistClus = new TH2F("fHistClus","# Cluster (TPC vs. ITS)",10,0,10,200,0,200);\r
- fList->Add(fHistClus);\r
- fHistChi2 = new TH1F("fHistChi2","Chi2/NDF distribution",200,0,10);\r
- fList->Add(fHistChi2);\r
- fHistDCA = new TH2F("fHistDCA","DCA (xy vs. z)",400,-5,5,400,-5,5); \r
- fList->Add(fHistDCA);\r
- fHistPt = new TH1F("fHistPt","p_{T} distribution",200,0,10);\r
- fList->Add(fHistPt);\r
- fHistEta = new TH1F("fHistEta","#eta distribution",200,-2,2);\r
- fList->Add(fHistEta);\r
- fHistPhi = new TH1F("fHistPhi","#phi distribution",200,-20,380);\r
- fList->Add(fHistPhi);\r
- fHistPhiBefore = new TH1F("fHistPhiBefore","#phi distribution",200,0.,2*TMath::Pi());\r
- fList->Add(fHistPhiBefore);\r
- fHistPhiAfter = new TH1F("fHistPhiAfter","#phi distribution",200,0.,2*TMath::Pi());\r
- fList->Add(fHistPhiAfter);\r
- fHistV0M = new TH2F("fHistV0M","V0 Multiplicity C vs. A",500, 0, 20000, 500, 0, 20000);\r
- fList->Add(fHistV0M);\r
- TString gRefTrackName[6] = {"tracks","tracksPos","tracksNeg","tracksTPConly","clusITS0","clusITS1"};\r
- fHistRefTracks = new TH2F("fHistRefTracks","Nr of Ref tracks/event vs. ref track estimator;;Nr of tracks",6, 0, 6, 400, 0, 20000);\r
- for(Int_t i = 1; i <= 6; i++)\r
- fHistRefTracks->GetXaxis()->SetBinLabel(i,gRefTrackName[i-1].Data());\r
- fList->Add(fHistRefTracks);\r
-\r
- //------------------------------------------------\r
- // V0 Multiplicity Histograms\r
- //------------------------------------------------\r
- if(fRunV0){\r
- fHistListV0 = new TList();\r
- fHistListV0->SetOwner(); // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner\r
- \r
- if(! fHistV0MultiplicityBeforeTrigSel) {\r
- fHistV0MultiplicityBeforeTrigSel = new TH1F("fHistV0MultiplicityBeforeTrigSel", \r
- "V0s per event (before Trig. Sel.);Nbr of V0s/Evt;Events", \r
- 25, 0, 25);\r
- fHistListV0->Add(fHistV0MultiplicityBeforeTrigSel);\r
- }\r
- \r
- if(! fHistV0MultiplicityForTrigEvt) {\r
- fHistV0MultiplicityForTrigEvt = new TH1F("fHistV0MultiplicityForTrigEvt", \r
- "V0s per event (for triggered evt);Nbr of V0s/Evt;Events", \r
- 25, 0, 25);\r
- fHistListV0->Add(fHistV0MultiplicityForTrigEvt);\r
- }\r
- \r
- if(! fHistV0MultiplicityForSelEvt) {\r
- fHistV0MultiplicityForSelEvt = new TH1F("fHistV0MultiplicityForSelEvt", \r
- "V0s per event;Nbr of V0s/Evt;Events", \r
- 25, 0, 25);\r
- fHistListV0->Add(fHistV0MultiplicityForSelEvt);\r
- }\r
- \r
- if(! fHistV0MultiplicityForSelEvtNoTPCOnly) {\r
- fHistV0MultiplicityForSelEvtNoTPCOnly = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnly", \r
- "V0s per event;Nbr of V0s/Evt;Events", \r
- 25, 0, 25);\r
- fHistListV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnly);\r
- }\r
- \r
- if(! fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup) {\r
- fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup", \r
- "V0s per event;Nbr of V0s/Evt;Events", \r
- 25, 0, 25);\r
- fHistListV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup);\r
- }\r
- \r
- //------------------------------------------------\r
- // Track Multiplicity Histograms\r
- //------------------------------------------------\r
- \r
- if(! fHistMultiplicityBeforeTrigSel) {\r
- fHistMultiplicityBeforeTrigSel = new TH1F("fHistMultiplicityBeforeTrigSel", \r
- "Tracks per event;Nbr of Tracks;Events", \r
- 200, 0, 200); \r
- fHistListV0->Add(fHistMultiplicityBeforeTrigSel);\r
- }\r
- if(! fHistMultiplicityForTrigEvt) {\r
- fHistMultiplicityForTrigEvt = new TH1F("fHistMultiplicityForTrigEvt", \r
- "Tracks per event;Nbr of Tracks;Events", \r
- 200, 0, 200); \r
- fHistListV0->Add(fHistMultiplicityForTrigEvt);\r
- }\r
- if(! fHistMultiplicity) {\r
- fHistMultiplicity = new TH1F("fHistMultiplicity", \r
- "Tracks per event;Nbr of Tracks;Events", \r
- 200, 0, 200); \r
- fHistListV0->Add(fHistMultiplicity);\r
- }\r
- if(! fHistMultiplicityNoTPCOnly) {\r
- fHistMultiplicityNoTPCOnly = new TH1F("fHistMultiplicityNoTPCOnly", \r
- "Tracks per event;Nbr of Tracks;Events", \r
- 200, 0, 200); \r
- fHistListV0->Add(fHistMultiplicityNoTPCOnly);\r
- }\r
- if(! fHistMultiplicityNoTPCOnlyNoPileup) {\r
- fHistMultiplicityNoTPCOnlyNoPileup = new TH1F("fHistMultiplicityNoTPCOnlyNoPileup", \r
- "Tracks per event;Nbr of Tracks;Events", \r
- 200, 0, 200); \r
- fHistListV0->Add(fHistMultiplicityNoTPCOnlyNoPileup);\r
- }\r
-\r
- //------------------------------------------------\r
- // V0 selection Histograms (before)\r
- //------------------------------------------------\r
- if(!fHistV0InvMassK0) {\r
- fHistV0InvMassK0 = new TH1F("fHistV0InvMassK0",\r
- "Invariant Mass for K0;Mass (GeV/c^{2});Events",\r
- 200,0,2);\r
- fHistListV0->Add(fHistV0InvMassK0);\r
- }\r
- if(!fHistV0InvMassLambda) {\r
- fHistV0InvMassLambda = new TH1F("fHistV0InvMassLambda",\r
- "Invariant Mass for Lambda;Mass (GeV/c^{2});Events",\r
- 200,0,2);\r
- fHistListV0->Add(fHistV0InvMassLambda);\r
- }\r
- if(!fHistV0InvMassAntiLambda) {\r
- fHistV0InvMassAntiLambda = new TH1F("fHistV0InvMassAntiLambda",\r
- "Invariant Mass for AntiLambda;Mass (GeV/c^{2});Events",\r
- 200,0,2);\r
- fHistListV0->Add(fHistV0InvMassAntiLambda);\r
- }\r
- if(!fHistV0Armenteros) {\r
- fHistV0Armenteros = new TH2F("fHistV0Armenteros",\r
- "Armenteros plot;#alpha;q_{t}",\r
- 200,-1,1,200,0,0.5);\r
- fHistListV0->Add(fHistV0Armenteros);\r
- }\r
- \r
- //------------------------------------------------\r
- // V0 selection Histograms (after)\r
- //------------------------------------------------\r
- if(!fHistV0SelInvMassK0) {\r
- fHistV0SelInvMassK0 = new TH1F("fHistV0SelInvMassK0",\r
- "Invariant Mass for K0;Mass (GeV/c^{2});Events",\r
- 200,0,2);\r
- fHistListV0->Add(fHistV0SelInvMassK0);\r
- }\r
- if(!fHistV0SelInvMassLambda) {\r
- fHistV0SelInvMassLambda = new TH1F("fHistV0SelInvMassLambda",\r
- "Invariant Mass for Lambda;Mass (GeV/c^{2});Events",\r
- 200,0,2);\r
- fHistListV0->Add(fHistV0SelInvMassLambda);\r
- }\r
- if(!fHistV0SelInvMassAntiLambda) {\r
- fHistV0SelInvMassAntiLambda = new TH1F("fHistV0SelInvMassAntiLambda",\r
- "Invariant Mass for AntiLambda;Mass (GeV/c^{2});Events",\r
- 200,0,2);\r
- fHistListV0->Add(fHistV0SelInvMassAntiLambda);\r
- }\r
- if(!fHistV0SelArmenteros) {\r
- fHistV0SelArmenteros = new TH2F("fHistV0SelArmenteros",\r
- "Armenteros plot;#alpha;q_{t}",\r
- 200,-1,1,200,0,0.5);\r
- fHistListV0->Add(fHistV0SelArmenteros);\r
- }\r
- }//V0\r
- \r
- // Balance function histograms\r
- // Initialize histograms if not done yet\r
- if(!fBalance->GetHistNp()){\r
- AliWarning("Histograms not yet initialized! --> Will be done now");\r
- AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
- fBalance->InitHistograms();\r
- }\r
-\r
- if(fRunShuffling) {\r
- if(!fShuffledBalance->GetHistNp()) {\r
- AliWarning("Histograms (shuffling) not yet initialized! --> Will be done now");\r
- AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
- fShuffledBalance->InitHistograms();\r
- }\r
- }\r
-\r
- if(fRunMixing) {\r
- if(!fMixedBalance->GetHistNp()) {\r
- AliWarning("Histograms (mixing) not yet initialized! --> Will be done now");\r
- AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
- fMixedBalance->InitHistograms();\r
- }\r
- }\r
-\r
- fListTriggeredBF->Add(fBalance->GetHistNp());\r
- fListTriggeredBF->Add(fBalance->GetHistNn());\r
- fListTriggeredBF->Add(fBalance->GetHistNpn());\r
- fListTriggeredBF->Add(fBalance->GetHistNnn());\r
- fListTriggeredBF->Add(fBalance->GetHistNpp());\r
- fListTriggeredBF->Add(fBalance->GetHistNnp());\r
- \r
- if(fRunShuffling) {\r
- fListTriggeredBFS->Add(fShuffledBalance->GetHistNp());\r
- fListTriggeredBFS->Add(fShuffledBalance->GetHistNn());\r
- fListTriggeredBFS->Add(fShuffledBalance->GetHistNpn());\r
- fListTriggeredBFS->Add(fShuffledBalance->GetHistNnn());\r
- fListTriggeredBFS->Add(fShuffledBalance->GetHistNpp());\r
- fListTriggeredBFS->Add(fShuffledBalance->GetHistNnp());\r
- } \r
-\r
- if(fRunMixing) {\r
- fListTriggeredBFM->Add(fMixedBalance->GetHistNp());\r
- fListTriggeredBFM->Add(fMixedBalance->GetHistNn());\r
- fListTriggeredBFM->Add(fMixedBalance->GetHistNpn());\r
- fListTriggeredBFM->Add(fMixedBalance->GetHistNnn());\r
- fListTriggeredBFM->Add(fMixedBalance->GetHistNpp());\r
- fListTriggeredBFM->Add(fMixedBalance->GetHistNnp());\r
- } \r
-\r
- // PID Response task active?\r
- if(fRunV0) {\r
- fPIDResponse = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetPIDResponse();\r
- if (!fPIDResponse) AliFatal("This Task needs the PID response attached to the inputHandler");\r
- }\r
-\r
- // Event Mixing\r
- Int_t trackDepth = fMixingTracks; \r
- Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemented of AliEventPoolManager\r
- \r
- Double_t centralityBins[] = {0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,15.,20.,25.,30.,35.,40.,45.,50.,55.,60.,65.,70.,75.,80.,90.,100.}; // SHOULD BE DEDUCED FROM CREATED ALITHN!!!\r
- Double_t* centbins = centralityBins;\r
- Int_t nCentralityBins = sizeof(centralityBins) / sizeof(Double_t) - 1;\r
- \r
- // bins for second buffer are shifted by 100 cm\r
- Double_t vertexBins[] = {-10., -7., -5., -3., -1., 1., 3., 5., 7., 10.}; // SHOULD BE DEDUCED FROM CREATED ALITHN!!!\r
- Double_t* vtxbins = vertexBins;\r
- Int_t nVertexBins = sizeof(vertexBins) / sizeof(Double_t) - 1;\r
-\r
- fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, centbins, nVertexBins, vtxbins);\r
-\r
-\r
- // Post output data.\r
- PostData(1, fList);\r
- PostData(2, fListTriggeredBF);\r
- if(fRunShuffling) PostData(3, fListTriggeredBFS);\r
- if(fRunMixing) PostData(4, fListTriggeredBFM);\r
- if(fRunV0) PostData(5,fHistListV0);\r
-\r
- TH1::AddDirectory(oldStatus);\r
-\r
-}\r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskTriggeredBF::UserExec(Option_t *) {\r
- // Main loop\r
- // Called for each event\r
-\r
- TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
- Float_t fCentrality = -1.; \r
-\r
- // ------------------------------------------------------------- \r
- // AOD analysis (vertex and track cuts also here!!!!)\r
- if(gAnalysisLevel == "AOD") {\r
- AliVEvent* eventMain = dynamic_cast<AliVEvent*>(InputEvent()); \r
- if(!eventMain) {\r
- AliError("eventMain not available");\r
- return;\r
- }\r
-\r
- // check event cuts and fill event histograms\r
- if((fCentrality = IsEventAccepted(eventMain)) < 0){\r
- return;\r
- }\r
- \r
- // get the accepted tracks in main event\r
- TObjArray *tracksMain = NULL;\r
- if(fRunV0) tracksMain = GetAcceptedV0s(eventMain);\r
- else tracksMain = GetAcceptedTracks(eventMain);\r
-\r
- // store charges of all accepted tracks, shuffle and reassign (two extra loops!)\r
- TObjArray* tracksShuffled = NULL;\r
- if(fRunShuffling){\r
- tracksShuffled = GetShuffledTracks(tracksMain);\r
- }\r
- \r
- // Event mixing --> UPDATE POOL IS MISSING!!!\r
- if (fRunMixing)\r
- {\r
- // 1. First get an event pool corresponding in mult (cent) and\r
- // zvertex to the current event. Once initialized, the pool\r
- // should contain nMix (reduced) events. This routine does not\r
- // pre-scan the chain. The first several events of every chain\r
- // will be skipped until the needed pools are filled to the\r
- // specified depth. If the pool categories are not too rare, this\r
- // should not be a problem. If they are rare, you could lose`\r
- // statistics.\r
- \r
- // 2. Collect the whole pool's content of tracks into one TObjArray\r
- // (bgTracks), which is effectively a single background super-event.\r
- \r
- // 3. The reduced and bgTracks arrays must both be passed into\r
- // FillCorrelations(). Also nMix should be passed in, so a weight\r
- // of 1./nMix can be applied.\r
- \r
- AliEventPool* pool = fPoolMgr->GetEventPool(fCentrality, eventMain->GetPrimaryVertex()->GetZ());\r
- \r
- if (!pool){\r
- AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCentrality, eventMain->GetPrimaryVertex()->GetZ()));\r
- }\r
- else{\r
-\r
- //pool->SetDebug(1);\r
- \r
- if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5){ \r
- \r
- \r
- Int_t nMix = pool->GetCurrentNEvents();\r
- //cout << "nMix = " << nMix << " tracks in pool = " << pool->NTracksInPool() << endl;\r
- \r
- //((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(2);\r
- //((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());\r
- //if (pool->IsReady())\r
- //((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(3);\r
- \r
- // Fill mixed-event histos here \r
- for (Int_t jMix=0; jMix<nMix; jMix++) \r
- {\r
- TObjArray* tracksMixed = pool->GetEvent(jMix);\r
- fMixedBalance->FillBalance(fCentrality,tracksMain,tracksMixed); \r
- }\r
- }\r
- \r
- // Update the Event pool\r
- pool->UpdatePool(tracksMain);\r
- //pool->PrintInfo();\r
- \r
- }//pool NULL check \r
- }//run mixing\r
- \r
- // calculate balance function\r
- fBalance->FillBalance(fCentrality,tracksMain,NULL);\r
- \r
- // calculate shuffled balance function\r
- if(fRunShuffling && tracksShuffled != NULL) {\r
- fShuffledBalance->FillBalance(fCentrality,tracksShuffled,NULL);\r
- }\r
- \r
- }//AOD analysis\r
- else{\r
- AliError("Triggered Balance Function analysis only for AODs!");\r
- }\r
-} \r
-\r
-//________________________________________________________________________\r
-Float_t AliAnalysisTaskTriggeredBF::IsEventAccepted(AliVEvent *event){\r
- // Checks the Event cuts\r
- // Fills Event statistics histograms\r
- \r
- // event selection done in AliAnalysisTaskSE::Exec() --> this is not used\r
- fHistEventStats->Fill(1); //all events\r
-\r
- Bool_t isSelectedMain = kTRUE;\r
- Float_t fCentrality = -1.;\r
- Int_t nV0s = event->GetNumberOfV0s();\r
- TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
- \r
- if(fUseOfflineTrigger)\r
- isSelectedMain = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();\r
- \r
- //V0 QA histograms (before trigger selection)\r
- if(fRunV0){\r
- fHistMultiplicityBeforeTrigSel->Fill ( -1 );\r
- fHistV0MultiplicityBeforeTrigSel->Fill ( nV0s );\r
- }\r
- \r
- if(isSelectedMain) {\r
- fHistEventStats->Fill(2); //triggered events\r
- \r
- //Centrality stuff \r
- if(fUseCentrality) {\r
- if(gAnalysisLevel == "AOD") { //centrality in AOD header\r
- AliAODHeader *header = (AliAODHeader*) event->GetHeader();\r
- fCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
-\r
- // QA for centrality estimators\r
- fHistCentStats->Fill(0.,header->GetCentralityP()->GetCentralityPercentile("V0M"));\r
- fHistCentStats->Fill(1.,header->GetCentralityP()->GetCentralityPercentile("FMD"));\r
- fHistCentStats->Fill(2.,header->GetCentralityP()->GetCentralityPercentile("TRK"));\r
- fHistCentStats->Fill(3.,header->GetCentralityP()->GetCentralityPercentile("TKL"));\r
- fHistCentStats->Fill(4.,header->GetCentralityP()->GetCentralityPercentile("CL0"));\r
- fHistCentStats->Fill(5.,header->GetCentralityP()->GetCentralityPercentile("CL1"));\r
- fHistCentStats->Fill(6.,header->GetCentralityP()->GetCentralityPercentile("V0MvsFMD"));\r
- fHistCentStats->Fill(7.,header->GetCentralityP()->GetCentralityPercentile("TKLvsV0M"));\r
- fHistCentStats->Fill(8.,header->GetCentralityP()->GetCentralityPercentile("ZEMvsZDC"));\r
- \r
- // centrality QA (V0M)\r
- fHistV0M->Fill(event->GetVZEROData()->GetMTotV0A(), event->GetVZEROData()->GetMTotV0C());\r
- \r
- // centrality QA (reference tracks)\r
- fHistRefTracks->Fill(0.,header->GetRefMultiplicity());\r
- fHistRefTracks->Fill(1.,header->GetRefMultiplicityPos());\r
- fHistRefTracks->Fill(2.,header->GetRefMultiplicityNeg());\r
- fHistRefTracks->Fill(3.,header->GetTPConlyRefMultiplicity());\r
- fHistRefTracks->Fill(4.,header->GetNumberOfITSClusters(0));\r
- fHistRefTracks->Fill(5.,header->GetNumberOfITSClusters(1));\r
- fHistRefTracks->Fill(6.,header->GetNumberOfITSClusters(2));\r
- fHistRefTracks->Fill(7.,header->GetNumberOfITSClusters(3));\r
- fHistRefTracks->Fill(8.,header->GetNumberOfITSClusters(4));\r
-\r
- //V0 QA histograms (after trigger selection)\r
- if(fRunV0){\r
- fHistMultiplicityForTrigEvt->Fill ( fCentrality );\r
- fHistV0MultiplicityForTrigEvt->Fill ( nV0s );\r
- }\r
- }\r
- }\r
- \r
- \r
- const AliVVertex *vertex = event->GetPrimaryVertex();\r
- \r
- if(vertex) {\r
- Double32_t fCov[6];\r
- vertex->GetCovarianceMatrix(fCov);\r
- if(vertex->GetNContributors() > 0) {\r
- if(fCov[5] != 0) {\r
- fHistEventStats->Fill(3); //events with a proper vertex\r
- if(TMath::Abs(vertex->GetX()) < fVxMax) {\r
- if(TMath::Abs(vertex->GetY()) < fVyMax) {\r
- if(TMath::Abs(vertex->GetZ()) < fVzMax) {\r
- fHistEventStats->Fill(4); //analyzed events\r
- fHistVx->Fill(vertex->GetX());\r
- fHistVy->Fill(vertex->GetY());\r
- fHistVz->Fill(vertex->GetZ());\r
-\r
-\r
-\r
- //V0 QA histograms (vertex Z check)\r
- if(fRunV0){\r
- fHistV0MultiplicityForSelEvt ->Fill( nV0s );\r
- fHistMultiplicity->Fill(fCentrality);\r
-\r
- //V0 QA histograms (Only look at events with well-established PV)\r
- const AliAODVertex *lPrimarySPDVtx = ((AliAODEvent*)event)->GetPrimaryVertexSPD();\r
- if(lPrimarySPDVtx){\r
- fHistMultiplicityNoTPCOnly->Fill ( fCentrality );\r
- fHistV0MultiplicityForSelEvtNoTPCOnly->Fill ( nV0s );\r
- \r
- //V0 QA histograms (Pileup Rejection)\r
- // FIXME : quality selection regarding pile-up rejection \r
- fHistMultiplicityNoTPCOnlyNoPileup->Fill(fCentrality);\r
- fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup ->Fill( nV0s );\r
-\r
- }\r
- else{\r
- return -1;\r
- }\r
- }\r
- \r
- \r
- // take only events inside centrality class\r
- if((fCentrality > fCentralityPercentileMin) && (fCentrality < fCentralityPercentileMax)){\r
- return fCentrality; \r
- }//centrality class\r
- }//Vz cut\r
- }//Vy cut\r
- }//Vx cut\r
- }//proper vertex resolution\r
- }//proper number of contributors\r
- }//vertex object valid\r
- }//triggered event \r
- \r
- // in all other cases return -1 (event not accepted)\r
- return -1;\r
-}\r
-\r
-//________________________________________________________________________\r
-TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedTracks(AliVEvent *event){\r
- // Returns TObjArray with tracks after all track cuts (only for AOD!)\r
- // Fills QA histograms\r
-\r
- //output TObjArray holding all good tracks\r
- TObjArray* tracksAccepted = new TObjArray;\r
- tracksAccepted->SetOwner(kTRUE);\r
-\r
- Short_t vCharge = 0;\r
- Double_t vEta = 0.;\r
- Double_t vPhi = 0.;\r
- Double_t vPt = 0.;\r
- \r
- // Loop over tracks in event\r
- for (Int_t iTracks = 0; iTracks < event->GetNumberOfTracks(); iTracks++) {\r
- AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(event->GetTrack(iTracks));\r
- if (!aodTrack) {\r
- AliError(Form("Could not receive track %d", iTracks));\r
- continue;\r
- }\r
- \r
- // AOD track cuts\r
- \r
- // For ESD Filter Information: ANALYSIS/macros/AddTaskESDfilter.C\r
- // take only TPC only tracks \r
- fHistTrackStats->Fill(aodTrack->GetFilterMap());\r
- if(!aodTrack->TestFilterBit(nAODtrackCutBit)) continue;\r
- \r
- vCharge = aodTrack->Charge();\r
- vEta = aodTrack->Eta();\r
- vPhi = aodTrack->Phi() * TMath::RadToDeg();\r
- vPt = aodTrack->Pt();\r
- \r
- Float_t dcaXY = aodTrack->DCA(); // this is the DCA from global track (not exactly what is cut on)\r
- Float_t dcaZ = aodTrack->ZAtDCA(); // this is the DCA from global track (not exactly what is cut on)\r
- \r
- \r
- // Kinematics cuts from ESD track cuts\r
- if( vPt < fPtMin || vPt > fPtMax) continue;\r
- if( vEta < fEtaMin || vEta > fEtaMax) continue;\r
- \r
- // Extra DCA cuts (for systematic studies [!= -1])\r
- if( fDCAxyCut != -1 && fDCAzCut != -1){\r
- if(TMath::Sqrt((dcaXY*dcaXY)/(fDCAxyCut*fDCAxyCut)+(dcaZ*dcaZ)/(fDCAzCut*fDCAzCut)) > 1 ){\r
- continue; // 2D cut\r
- }\r
- }\r
- \r
- // Extra TPC cuts (for systematic studies [!= -1])\r
- if( fTPCchi2Cut != -1 && aodTrack->Chi2perNDF() > fTPCchi2Cut){\r
- continue;\r
- }\r
- if( fNClustersTPCCut != -1 && aodTrack->GetTPCNcls() < fNClustersTPCCut){\r
- continue;\r
- }\r
- \r
- // fill QA histograms\r
- fHistClus->Fill(aodTrack->GetITSNcls(),aodTrack->GetTPCNcls());\r
- fHistDCA->Fill(dcaZ,dcaXY);\r
- fHistChi2->Fill(aodTrack->Chi2perNDF());\r
- fHistPt->Fill(vPt);\r
- fHistEta->Fill(vEta);\r
- fHistPhi->Fill(vPhi);\r
- \r
- // add the track to the TObjArray\r
- tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge,1.));\r
- }\r
-\r
- return tracksAccepted;\r
-}\r
-\r
-//________________________________________________________________________\r
-TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedV0s(AliVEvent *event){\r
- // Returns TObjArray with tracks after all track cuts (only for AOD!)\r
- // Fills QA histograms\r
-\r
- //output TObjArray holding all good tracks\r
- TObjArray* tracksAccepted = new TObjArray;\r
- tracksAccepted->SetOwner(kTRUE);\r
-\r
- Short_t vCharge = 0;\r
- Double_t vEta = 0.;\r
- Double_t vPhi = 0.;\r
- Double_t vPt = 0.;\r
- \r
- //------------------------------------------------\r
- // MAIN LAMBDA LOOP STARTS HERE (basically a copy of AliAnalysisTaskExtractV0AOD)\r
- //------------------------------------------------\r
-\r
- // parameters (for the time being hard coded here) --> from David for EbyE Lambdas\r
- Bool_t fkUseOnTheFly = kFALSE;\r
- Double_t fRapidityBoundary = 0.5; \r
- Double_t fCutDaughterEta = 0.8;\r
- Double_t fCutV0Radius = 0.9;\r
- Double_t fCutDCANegToPV = 0.1;\r
- Double_t fCutDCAPosToPV = 0.1;\r
- Double_t fCutDCAV0Daughters = 1.0;\r
- Double_t fCutV0CosPA = 0.9995;\r
- Double_t fMassLambda = 1.115683;\r
- Double_t fCutMassLambda = 0.007;\r
- Double_t fCutProperLifetime = 3*7.9;\r
- Double_t fCutLeastNumberOfCrossedRows = 70;\r
- Double_t fCutLeastNumberOfCrossedRowsOverFindable = 0.8;\r
- Double_t fCutTPCPIDNSigmasProton = 3.0;\r
- Double_t fCutTPCPIDNSigmasPion = 5.0;\r
-\r
-\r
- //Variable definition\r
- Int_t lOnFlyStatus = 0;// nv0sOn = 0, nv0sOff = 0;\r
- Double_t lChi2V0 = 0;\r
- Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;\r
- Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;\r
- Double_t lV0CosineOfPointingAngle = 0;\r
- Double_t lV0Radius = 0, lPt = 0;\r
- Double_t lEta = 0, lPhi = 0;\r
- Double_t lRap = 0, lRapK0Short = 0, lRapLambda = 0;\r
- Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;\r
- Double_t lAlphaV0 = 0, lPtArmV0 = 0;\r
- \r
- Double_t fMinV0Pt = 0; \r
- Double_t fMaxV0Pt = 100; \r
- \r
-\r
- \r
- // some event observables\r
- Int_t nv0s = event->GetNumberOfV0s();\r
- Double_t tPrimaryVtxPosition[3];\r
- const AliVVertex *primaryVtx = event->GetPrimaryVertex();\r
- tPrimaryVtxPosition[0] = primaryVtx->GetX();\r
- tPrimaryVtxPosition[1] = primaryVtx->GetY();\r
- tPrimaryVtxPosition[2] = primaryVtx->GetZ();\r
-\r
-\r
- //loop over V0s \r
- for (Int_t iV0 = 0; iV0 < nv0s; iV0++) \r
- {// This is the begining of the V0 loop\r
- AliAODv0 *v0 = ((AliAODEvent*)event)->GetV0(iV0);\r
- if (!v0) continue;\r
-\r
- //Obsolete at AOD level... \r
- //---> Fix On-the-Fly candidates, count how many swapped\r
- //if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ){\r
- // fHistSwappedV0Counter -> Fill( 1 );\r
- //}else{\r
- // fHistSwappedV0Counter -> Fill( 0 ); \r
- //}\r
- //if ( fkUseOnTheFly ) CheckChargeV0(v0); \r
- \r
- Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0); \r
- Double_t tV0mom[3];\r
- v0->GetPxPyPz( tV0mom ); \r
- Double_t lV0TotalMomentum = TMath::Sqrt(\r
- tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );\r
- \r
- lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);\r
- lPt = v0->Pt();\r
- lEta = v0->Eta();\r
- lPhi = v0->Phi()*TMath::RadToDeg();\r
- lRapK0Short = v0->RapK0Short();\r
- lRapLambda = v0->RapLambda();\r
- lRap = lRapLambda;//v0->Y(); //FIXME!!!\r
- if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;\r
- \r
- //UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPosID());\r
- //UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetPosID());\r
-\r
- Double_t lMomPos[3]; //v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);\r
- Double_t lMomNeg[3]; //v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);\r
- lMomPos[0] = v0->MomPosX();\r
- lMomPos[1] = v0->MomPosY();\r
- lMomPos[2] = v0->MomPosZ();\r
- lMomNeg[0] = v0->MomNegX();\r
- lMomNeg[1] = v0->MomNegY();\r
- lMomNeg[2] = v0->MomNegZ();\r
- \r
- AliAODTrack *pTrack=(AliAODTrack *)v0->GetDaughter(0); //0->Positive Daughter\r
- AliAODTrack *nTrack=(AliAODTrack *)v0->GetDaughter(1); //1->Negative Daughter\r
- if (!pTrack || !nTrack) {\r
- AliError("ERROR: Could not retreive one of the daughter track");\r
- continue;\r
- }\r
-\r
- //Daughter Eta for Eta selection, afterwards\r
- Double_t lNegEta = nTrack->Eta();\r
- Double_t lPosEta = pTrack->Eta();\r
- \r
- // Filter like-sign V0 (next: add counter and distribution)\r
- if ( pTrack->Charge() == nTrack->Charge()){\r
- continue;\r
- } \r
- \r
- //Quick test this far! \r
- \r
-\r
- //________________________________________________________________________\r
- // Track quality cuts \r
- Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);\r
- Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);\r
- Float_t lLeastNbrCrossedRows = (lPosTrackCrossedRows>lNegTrackCrossedRows) ? lNegTrackCrossedRows : lPosTrackCrossedRows;\r
-\r
- // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)\r
- if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;\r
- if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;\r
- \r
- if ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) ) continue;\r
- \r
- //Findable clusters > 0 condition\r
- if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;\r
- \r
- //Compute ratio Crossed Rows / Findable clusters\r
- //Note: above test avoids division by zero! \r
- Float_t lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF())); \r
- Float_t lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF())); \r
- Float_t lLeastNbrCrossedRowsOverFindable = (lPosTrackCrossedRowsOverFindable>lNegTrackCrossedRowsOverFindable) ? lNegTrackCrossedRowsOverFindable : lPosTrackCrossedRowsOverFindable;\r
-\r
- //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here\r
- if ( lLeastNbrCrossedRowsOverFindable < 0.8) continue;\r
- \r
- //End track Quality Cuts\r
- //________________________________________________________________________\r
- \r
- \r
- lDcaPosToPrimVertex = v0->DcaPosToPrimVertex();\r
- lDcaNegToPrimVertex = v0->DcaNegToPrimVertex();\r
- \r
- lOnFlyStatus = v0->GetOnFlyStatus();\r
- lChi2V0 = v0->Chi2V0();\r
- lDcaV0Daughters = v0->DcaV0Daughters();\r
- lDcaV0ToPrimVertex = v0->DcaV0ToPrimVertex();\r
- lV0CosineOfPointingAngle = v0->CosPointingAngle(tPrimaryVtxPosition);\r
- \r
- // Distance over total momentum\r
- Double_t lDistOverTotMom = TMath::Sqrt(\r
- TMath::Power( tDecayVertexV0[0] - tPrimaryVtxPosition[0] , 2) +\r
- TMath::Power( tDecayVertexV0[1] - tPrimaryVtxPosition[1] , 2) +\r
- TMath::Power( tDecayVertexV0[2] - tPrimaryVtxPosition[2] , 2)\r
- );\r
- lDistOverTotMom /= (lV0TotalMomentum+1e-10); //avoid division by zero, to be sure\r
- \r
- \r
- // Getting invariant mass infos directly from ESD\r
- lInvMassK0s = v0->MassK0Short();\r
- lInvMassLambda = v0->MassLambda();\r
- lInvMassAntiLambda = v0->MassAntiLambda();\r
- lAlphaV0 = v0->AlphaV0();\r
- lPtArmV0 = v0->PtArmV0();\r
-\r
- //Official means of acquiring N-sigmas \r
- Double_t lNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );\r
- Double_t lNSigmasPosPion = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );\r
- Double_t lNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );\r
- Double_t lNSigmasNegPion = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );\r
-\r
- //V0 QA histograms (before V0 selection)\r
- fHistV0InvMassK0->Fill(lInvMassK0s);\r
- fHistV0InvMassLambda->Fill(lInvMassLambda);\r
- fHistV0InvMassAntiLambda->Fill(lInvMassAntiLambda);\r
- fHistV0Armenteros->Fill(lAlphaV0,lPtArmV0);\r
- \r
- \r
- //First Selection: Reject OnFly\r
- if( (lOnFlyStatus == 0 && fkUseOnTheFly == kFALSE) || (lOnFlyStatus != 0 && fkUseOnTheFly == kTRUE ) ){\r
- \r
-\r
- //Second Selection: rough 20-sigma band, parametric. \r
- //K0Short: Enough to parametrize peak broadening with linear function. \r
- Double_t lUpperLimitK0Short = (5.63707e-01) + (1.14979e-02)*lPt; \r
- Double_t lLowerLimitK0Short = (4.30006e-01) - (1.10029e-02)*lPt;\r
- \r
- //Lambda: Linear (for higher pt) plus exponential (for low-pt broadening)\r
- //[0]+[1]*x+[2]*TMath::Exp(-[3]*x)\r
- Double_t lUpperLimitLambda = (1.13688e+00) + (5.27838e-03)*lPt + (8.42220e-02)*TMath::Exp(-(3.80595e+00)*lPt); \r
- Double_t lLowerLimitLambda = (1.09501e+00) - (5.23272e-03)*lPt - (7.52690e-02)*TMath::Exp(-(3.46339e+00)*lPt);\r
- \r
- //Do Selection \r
- if( (lInvMassLambda < lUpperLimitLambda && lInvMassLambda > lLowerLimitLambda ) || \r
- (lInvMassAntiLambda < lUpperLimitLambda && lInvMassAntiLambda > lLowerLimitLambda ) || \r
- (lInvMassK0s < lUpperLimitK0Short && lInvMassK0s > lLowerLimitK0Short ) ){\r
-\r
-\r
- // //Pre-selection in case this is AA...\r
- // //if( fkIsNuclear == kFALSE ) fTree->Fill();\r
- // //if( fkIsNuclear == kTRUE){ \r
- // //If this is a nuclear collision___________________\r
- // // ... pre-filter with TPC, daughter eta selection\r
-\r
- \r
- if( (lInvMassLambda < lUpperLimitLambda && lInvMassLambda > lLowerLimitLambda \r
- && TMath::Abs(lNSigmasPosProton) < 6.0 && TMath::Abs(lNSigmasNegPion) < 6.0 ) || \r
- (lInvMassAntiLambda < lUpperLimitLambda && lInvMassAntiLambda > lLowerLimitLambda \r
- && TMath::Abs(lNSigmasNegProton) < 6.0 && TMath::Abs(lNSigmasPosPion) < 6.0 ) || \r
- (lInvMassK0s < lUpperLimitK0Short && lInvMassK0s > lLowerLimitK0Short \r
- && TMath::Abs(lNSigmasNegPion) < 6.0 && TMath::Abs(lNSigmasPosPion) < 6.0 ) ){\r
- \r
- //insane test\r
- if ( TMath::Abs(lNegEta)<0.8 && TMath::Abs(lPosEta)<0.8 ){\r
-\r
- // start the fine selection (usually done in post processing, but we don't have time to waste) --> Lambdas!\r
- if(\r
- TMath::Abs(lRap)<fRapidityBoundary &&\r
- TMath::Abs(lNegEta) <= fCutDaughterEta && \r
- TMath::Abs(lPosEta) <= fCutDaughterEta &&\r
- lV0Radius >= fCutV0Radius &&\r
- lDcaNegToPrimVertex >= fCutDCANegToPV &&\r
- lDcaPosToPrimVertex >= fCutDCAPosToPV &&\r
- lDcaV0Daughters <= fCutDCAV0Daughters &&\r
- lV0CosineOfPointingAngle >= fCutV0CosPA && \r
- fMassLambda*lDistOverTotMom <= fCutProperLifetime &&\r
- lLeastNbrCrossedRows >= fCutLeastNumberOfCrossedRows &&\r
- lLeastNbrCrossedRowsOverFindable >= fCutLeastNumberOfCrossedRowsOverFindable &&\r
- lPtArmV0 * 5 < TMath::Abs(lAlphaV0) && \r
- ((TMath::Abs(lNSigmasNegPion) <= fCutTPCPIDNSigmasPion &&\r
- TMath::Abs(lNSigmasPosProton) <= fCutTPCPIDNSigmasProton) ||\r
- (TMath::Abs(lNSigmasPosPion) <= fCutTPCPIDNSigmasPion &&\r
- TMath::Abs(lNSigmasNegProton) <= fCutTPCPIDNSigmasProton)) \r
- )\r
- {\r
-\r
- //V0 QA histograms (after V0 selection)\r
- fHistV0SelInvMassK0->Fill(lInvMassK0s);\r
- fHistV0SelInvMassLambda->Fill(lInvMassLambda);\r
- fHistV0SelInvMassAntiLambda->Fill(lInvMassAntiLambda);\r
-\r
- // this means a V0 candidate is found\r
- if(TMath::Abs(lInvMassLambda-fMassLambda) < fCutMassLambda ||\r
- TMath::Abs(lInvMassAntiLambda-fMassLambda) < fCutMassLambda){\r
-\r
- fHistV0SelArmenteros->Fill(lAlphaV0,lPtArmV0); \r
-\r
- vEta = lEta;\r
- vPhi = lPhi;\r
- vPt = lPt;\r
- if(lAlphaV0 > 0) vCharge = 1;\r
- if(lAlphaV0 < 0) vCharge = -1;\r
-\r
- // fill QA histograms\r
- fHistPt->Fill(vPt);\r
- fHistEta->Fill(vEta);\r
- fHistPhi->Fill(vPhi);\r
- \r
- // add the track to the TObjArray\r
- tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge,1.));\r
- }\r
- }\r
- }\r
- }\r
- //}//end nuclear_____________________________________\r
- }\r
- }\r
- }//V0 loop\r
- \r
- return tracksAccepted;\r
-}\r
-\r
-//________________________________________________________________________\r
-TObjArray* AliAnalysisTaskTriggeredBF::GetShuffledTracks(TObjArray *tracks){\r
- // Clones TObjArray and returns it with tracks after shuffling the charges\r
-\r
- TObjArray* tracksShuffled = new TObjArray;\r
- tracksShuffled->SetOwner(kTRUE);\r
-\r
- vector<Short_t> *chargeVector = new vector<Short_t>; //original charge of accepted tracks \r
-\r
- for (Int_t i=0; i<tracks->GetEntriesFast(); i++)\r
- {\r
- AliVParticle* track = (AliVParticle*) tracks->At(i);\r
- chargeVector->push_back(track->Charge());\r
- } \r
- \r
- random_shuffle(chargeVector->begin(), chargeVector->end());\r
- \r
- for(Int_t i = 0; i < tracks->GetEntriesFast(); i++){\r
- AliVParticle* track = (AliVParticle*) tracks->At(i);\r
- tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i),1.));\r
- }\r
-\r
- delete chargeVector;\r
- \r
- return tracksShuffled;\r
-}\r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskTriggeredBF::FinishTaskOutput(){\r
- //checks if Balance Function objects are there (needed to write the histograms)\r
- if (!fBalance) {\r
- AliError("fBalance not available");\r
- return;\r
- } \r
- if(fRunShuffling) {\r
- if (!fShuffledBalance) {\r
- AliError("fShuffledBalance not available");\r
- return;\r
- }\r
- }\r
-\r
-}\r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskTriggeredBF::Terminate(Option_t *) {\r
- // Called once at the end of the query\r
-\r
- // not implemented ...\r
-\r
-}\r
-\r
-void AliAnalysisTaskTriggeredBF::UserExecMix(Option_t *)\r
-{\r
-\r
- // not yet done for event mixing!\r
- return;\r
-\r
-}\r
-\r
+#include <vector>
+#include "TChain.h"
+#include "TList.h"
+#include "TCanvas.h"
+#include "TLorentzVector.h"
+#include "TGraphErrors.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TArrayF.h"
+#include "TF1.h"
+#include "TRandom.h"
+
+#include "AliLog.h"
+
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisManager.h"
+
+#include "AliESDVertex.h"
+#include "AliESDEvent.h"
+#include "AliESDInputHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODTrack.h"
+#include "AliAODInputHandler.h"
+#include "AliGenEventHeader.h"
+#include "AliGenHijingEventHeader.h"
+#include "AliMCEventHandler.h"
+#include "AliMCEvent.h"
+#include "AliMixInputEventHandler.h"
+#include "AliStack.h"
+
+#include "TH2D.h"
+#include "AliTHn.h"
+
+#include "AliEventPoolManager.h"
+
+#include "AliAnalysisTaskTriggeredBF.h"
+#include "AliBalanceTriggered.h"
+
+
+// Analysis task for the TriggeredBF code
+// Authors: Panos.Christakoglou@nikhef.nl, m.weber@cern.ch
+
+// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
+//
+// For the V0 part:
+// --> AliAnalysisTaskExtractV0AOD (by david.chinellato@gmail.com)
+//
+// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
+
+using std::cout;
+using std::endl;
+using std::vector;
+
+ClassImp(AliAnalysisTaskTriggeredBF)
+
+//________________________________________________________________________
+AliAnalysisTaskTriggeredBF::AliAnalysisTaskTriggeredBF(const char *name)
+: AliAnalysisTaskSE(name),
+ fBalance(0),
+ fRunShuffling(kFALSE),
+ fShuffledBalance(0),
+ fRunMixing(kFALSE),
+ fMixingTracks(50000),
+ fMixedBalance(0),
+ fPoolMgr(0),
+ fRunV0(kFALSE),
+ fPIDResponse(0),
+ fPIDCombined(0),
+ fList(0),
+ fListTriggeredBF(0),
+ fListTriggeredBFS(0),
+ fListTriggeredBFM(0),
+ fHistListPIDQA(0),
+ fHistListV0(0),
+ fHistEventStats(0),
+ fHistCentStats(0),
+ fHistTriggerStats(0),
+ fHistTrackStats(0),
+ fHistVx(0),
+ fHistVy(0),
+ fHistVz(0),
+ fHistClus(0),
+ fHistDCA(0),
+ fHistChi2(0),
+ fHistPt(0),
+ fHistEta(0),
+ fHistPhi(0),
+ fHistPhiBefore(0),
+ fHistPhiAfter(0),
+ fHistV0M(0),
+ fHistRefTracks(0),
+ fHistV0MultiplicityBeforeTrigSel(0),
+ fHistV0MultiplicityForTrigEvt(0),
+ fHistV0MultiplicityForSelEvt(0),
+ fHistV0MultiplicityForSelEvtNoTPCOnly(0),
+ fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup(0),
+ fHistMultiplicityBeforeTrigSel(0),
+ fHistMultiplicityForTrigEvt(0),
+ fHistMultiplicity(0),
+ fHistMultiplicityNoTPCOnly(0),
+ fHistMultiplicityNoTPCOnlyNoPileup(0),
+ fHistV0InvMassK0(0),
+ fHistV0InvMassLambda(0),
+ fHistV0InvMassAntiLambda(0),
+ fHistV0Armenteros(0),
+ fHistV0SelInvMassK0(0),
+ fHistV0SelInvMassLambda(0),
+ fHistV0SelInvMassAntiLambda(0),
+ fHistV0SelArmenteros(0),
+ fCentralityEstimator("V0M"),
+ fUseCentrality(kFALSE),
+ fCentralityPercentileMin(0.),
+ fCentralityPercentileMax(5.),
+ fImpactParameterMin(0.),
+ fImpactParameterMax(20.),
+ fUseMultiplicity(kFALSE),
+ fNumberOfAcceptedTracksMin(0),
+ fNumberOfAcceptedTracksMax(10000),
+ fHistNumberOfAcceptedTracks(0),
+ fUseOfflineTrigger(kFALSE),
+ fVxMax(0.3),
+ fVyMax(0.3),
+ fVzMax(10.),
+ nAODtrackCutBit(128),
+ fPtMin(0.3),
+ fPtMax(1.5),
+ fEtaMin(-0.8),
+ fEtaMax(-0.8),
+ fDCAxyCut(-1),
+ fDCAzCut(-1),
+ fTPCchi2Cut(-1),
+ fNClustersTPCCut(-1)
+{
+ // Constructor
+ // Define input and output slots here
+ // Input slot #0 works with a TChain
+ DefineInput(0, TChain::Class());
+ // Output slot #0 writes into a TH1 container
+ DefineOutput(1, TList::Class());
+ DefineOutput(2, TList::Class());
+ DefineOutput(3, TList::Class());
+ DefineOutput(4, TList::Class());
+ DefineOutput(5, TList::Class());
+}
+
+//________________________________________________________________________
+AliAnalysisTaskTriggeredBF::~AliAnalysisTaskTriggeredBF() {
+
+ // Destructor
+
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskTriggeredBF::UserCreateOutputObjects() {
+ // Create histograms
+ // Called once
+
+ // global switch disabling the reference
+ // (to avoid "Replacing existing TH1" if several wagons are created in train)
+ Bool_t oldStatus = TH1::AddDirectoryStatus();
+ TH1::AddDirectory(kFALSE);
+
+ if(!fBalance) {
+ fBalance = new AliBalanceTriggered();
+ fBalance->SetAnalysisLevel("AOD");
+ }
+ if(fRunShuffling) {
+ if(!fShuffledBalance) {
+ fShuffledBalance = new AliBalanceTriggered();
+ fShuffledBalance->SetAnalysisLevel("AOD");
+ }
+ }
+ if(fRunMixing) {
+ if(!fMixedBalance) {
+ fMixedBalance = new AliBalanceTriggered();
+ fMixedBalance->SetAnalysisLevel("AOD");
+ }
+ }
+
+ //QA list
+ fList = new TList();
+ fList->SetName("listQA");
+ fList->SetOwner();
+
+ //Balance Function list
+ fListTriggeredBF = new TList();
+ fListTriggeredBF->SetName("listTriggeredBF");
+ fListTriggeredBF->SetOwner();
+
+ if(fRunShuffling) {
+ fListTriggeredBFS = new TList();
+ fListTriggeredBFS->SetName("listTriggeredBFShuffled");
+ fListTriggeredBFS->SetOwner();
+ }
+ if(fRunMixing) {
+ fListTriggeredBFM = new TList();
+ fListTriggeredBFM->SetName("listTriggeredBFMixed");
+ fListTriggeredBFM->SetOwner();
+ }
+
+
+ //Event stats.
+ TString gCutName[4] = {"Total","Offline trigger",
+ "Vertex","Analyzed"};
+ fHistEventStats = new TH1F("fHistEventStats",
+ "Event statistics;;N_{events}",
+ 4,0.5,4.5);
+ for(Int_t i = 1; i <= 4; i++)
+ fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());
+ fList->Add(fHistEventStats);
+
+ TString gCentName[9] = {"V0M","FMD","TRK","TKL","CL0","CL1","V0MvsFMD","TKLvsV0M","ZEMvsZDC"};
+ fHistCentStats = new TH2F("fHistCentStats",
+ "Centrality statistics;;Cent percentile",
+ 9,-0.5,8.5,220,-5,105);
+ for(Int_t i = 1; i <= 9; i++)
+ fHistCentStats->GetXaxis()->SetBinLabel(i,gCentName[i-1].Data());
+ fList->Add(fHistCentStats);
+
+ fHistTriggerStats = new TH1F("fHistTriggerStats","Trigger statistics;TriggerBit;N_{events}",130,0,130);
+ fList->Add(fHistTriggerStats);
+
+ fHistTrackStats = new TH1F("fHistTrackStats","Event statistics;TrackFilterBit;N_{events}",130,0,130);
+ fList->Add(fHistTrackStats);
+
+ fHistNumberOfAcceptedTracks = new TH1F("fHistNumberOfAcceptedTracks",";N_{acc.};Entries",4001,-0.5,4000.5);
+ fList->Add(fHistNumberOfAcceptedTracks);
+
+ // Vertex distributions
+ fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);
+ fList->Add(fHistVx);
+ fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);
+ fList->Add(fHistVy);
+ fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);
+ fList->Add(fHistVz);
+
+ // QA histograms
+ fHistClus = new TH2F("fHistClus","# Cluster (TPC vs. ITS)",10,0,10,200,0,200);
+ fList->Add(fHistClus);
+ fHistChi2 = new TH1F("fHistChi2","Chi2/NDF distribution",200,0,10);
+ fList->Add(fHistChi2);
+ fHistDCA = new TH2F("fHistDCA","DCA (xy vs. z)",400,-5,5,400,-5,5);
+ fList->Add(fHistDCA);
+ fHistPt = new TH1F("fHistPt","p_{T} distribution",200,0,10);
+ fList->Add(fHistPt);
+ fHistEta = new TH1F("fHistEta","#eta distribution",200,-2,2);
+ fList->Add(fHistEta);
+ fHistPhi = new TH1F("fHistPhi","#phi distribution",200,-20,380);
+ fList->Add(fHistPhi);
+ fHistPhiBefore = new TH1F("fHistPhiBefore","#phi distribution",200,0.,2*TMath::Pi());
+ fList->Add(fHistPhiBefore);
+ fHistPhiAfter = new TH1F("fHistPhiAfter","#phi distribution",200,0.,2*TMath::Pi());
+ fList->Add(fHistPhiAfter);
+ fHistV0M = new TH2F("fHistV0M","V0 Multiplicity C vs. A",500, 0, 20000, 500, 0, 20000);
+ fList->Add(fHistV0M);
+ TString gRefTrackName[6] = {"tracks","tracksPos","tracksNeg","tracksTPConly","clusITS0","clusITS1"};
+ fHistRefTracks = new TH2F("fHistRefTracks","Nr of Ref tracks/event vs. ref track estimator;;Nr of tracks",6, 0, 6, 400, 0, 20000);
+ for(Int_t i = 1; i <= 6; i++)
+ fHistRefTracks->GetXaxis()->SetBinLabel(i,gRefTrackName[i-1].Data());
+ fList->Add(fHistRefTracks);
+
+ //------------------------------------------------
+ // V0 Multiplicity Histograms
+ //------------------------------------------------
+ if(fRunV0){
+ fHistListV0 = new TList();
+ fHistListV0->SetOwner(); // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
+
+ if(! fHistV0MultiplicityBeforeTrigSel) {
+ fHistV0MultiplicityBeforeTrigSel = new TH1F("fHistV0MultiplicityBeforeTrigSel",
+ "V0s per event (before Trig. Sel.);Nbr of V0s/Evt;Events",
+ 25, 0, 25);
+ fHistListV0->Add(fHistV0MultiplicityBeforeTrigSel);
+ }
+
+ if(! fHistV0MultiplicityForTrigEvt) {
+ fHistV0MultiplicityForTrigEvt = new TH1F("fHistV0MultiplicityForTrigEvt",
+ "V0s per event (for triggered evt);Nbr of V0s/Evt;Events",
+ 25, 0, 25);
+ fHistListV0->Add(fHistV0MultiplicityForTrigEvt);
+ }
+
+ if(! fHistV0MultiplicityForSelEvt) {
+ fHistV0MultiplicityForSelEvt = new TH1F("fHistV0MultiplicityForSelEvt",
+ "V0s per event;Nbr of V0s/Evt;Events",
+ 25, 0, 25);
+ fHistListV0->Add(fHistV0MultiplicityForSelEvt);
+ }
+
+ if(! fHistV0MultiplicityForSelEvtNoTPCOnly) {
+ fHistV0MultiplicityForSelEvtNoTPCOnly = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnly",
+ "V0s per event;Nbr of V0s/Evt;Events",
+ 25, 0, 25);
+ fHistListV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnly);
+ }
+
+ if(! fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup) {
+ fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup",
+ "V0s per event;Nbr of V0s/Evt;Events",
+ 25, 0, 25);
+ fHistListV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup);
+ }
+
+ //------------------------------------------------
+ // Track Multiplicity Histograms
+ //------------------------------------------------
+
+ if(! fHistMultiplicityBeforeTrigSel) {
+ fHistMultiplicityBeforeTrigSel = new TH1F("fHistMultiplicityBeforeTrigSel",
+ "Tracks per event;Nbr of Tracks;Events",
+ 200, 0, 200);
+ fHistListV0->Add(fHistMultiplicityBeforeTrigSel);
+ }
+ if(! fHistMultiplicityForTrigEvt) {
+ fHistMultiplicityForTrigEvt = new TH1F("fHistMultiplicityForTrigEvt",
+ "Tracks per event;Nbr of Tracks;Events",
+ 200, 0, 200);
+ fHistListV0->Add(fHistMultiplicityForTrigEvt);
+ }
+ if(! fHistMultiplicity) {
+ fHistMultiplicity = new TH1F("fHistMultiplicity",
+ "Tracks per event;Nbr of Tracks;Events",
+ 200, 0, 200);
+ fHistListV0->Add(fHistMultiplicity);
+ }
+ if(! fHistMultiplicityNoTPCOnly) {
+ fHistMultiplicityNoTPCOnly = new TH1F("fHistMultiplicityNoTPCOnly",
+ "Tracks per event;Nbr of Tracks;Events",
+ 200, 0, 200);
+ fHistListV0->Add(fHistMultiplicityNoTPCOnly);
+ }
+ if(! fHistMultiplicityNoTPCOnlyNoPileup) {
+ fHistMultiplicityNoTPCOnlyNoPileup = new TH1F("fHistMultiplicityNoTPCOnlyNoPileup",
+ "Tracks per event;Nbr of Tracks;Events",
+ 200, 0, 200);
+ fHistListV0->Add(fHistMultiplicityNoTPCOnlyNoPileup);
+ }
+
+ //------------------------------------------------
+ // V0 selection Histograms (before)
+ //------------------------------------------------
+ if(!fHistV0InvMassK0) {
+ fHistV0InvMassK0 = new TH1F("fHistV0InvMassK0",
+ "Invariant Mass for K0;Mass (GeV/c^{2});Events",
+ 200,0,2);
+ fHistListV0->Add(fHistV0InvMassK0);
+ }
+ if(!fHistV0InvMassLambda) {
+ fHistV0InvMassLambda = new TH1F("fHistV0InvMassLambda",
+ "Invariant Mass for Lambda;Mass (GeV/c^{2});Events",
+ 200,0,2);
+ fHistListV0->Add(fHistV0InvMassLambda);
+ }
+ if(!fHistV0InvMassAntiLambda) {
+ fHistV0InvMassAntiLambda = new TH1F("fHistV0InvMassAntiLambda",
+ "Invariant Mass for AntiLambda;Mass (GeV/c^{2});Events",
+ 200,0,2);
+ fHistListV0->Add(fHistV0InvMassAntiLambda);
+ }
+ if(!fHistV0Armenteros) {
+ fHistV0Armenteros = new TH2F("fHistV0Armenteros",
+ "Armenteros plot;#alpha;q_{t}",
+ 200,-1,1,200,0,0.5);
+ fHistListV0->Add(fHistV0Armenteros);
+ }
+
+ //------------------------------------------------
+ // V0 selection Histograms (after)
+ //------------------------------------------------
+ if(!fHistV0SelInvMassK0) {
+ fHistV0SelInvMassK0 = new TH1F("fHistV0SelInvMassK0",
+ "Invariant Mass for K0;Mass (GeV/c^{2});Events",
+ 200,0,2);
+ fHistListV0->Add(fHistV0SelInvMassK0);
+ }
+ if(!fHistV0SelInvMassLambda) {
+ fHistV0SelInvMassLambda = new TH1F("fHistV0SelInvMassLambda",
+ "Invariant Mass for Lambda;Mass (GeV/c^{2});Events",
+ 200,0,2);
+ fHistListV0->Add(fHistV0SelInvMassLambda);
+ }
+ if(!fHistV0SelInvMassAntiLambda) {
+ fHistV0SelInvMassAntiLambda = new TH1F("fHistV0SelInvMassAntiLambda",
+ "Invariant Mass for AntiLambda;Mass (GeV/c^{2});Events",
+ 200,0,2);
+ fHistListV0->Add(fHistV0SelInvMassAntiLambda);
+ }
+ if(!fHistV0SelArmenteros) {
+ fHistV0SelArmenteros = new TH2F("fHistV0SelArmenteros",
+ "Armenteros plot;#alpha;q_{t}",
+ 200,-1,1,200,0,0.5);
+ fHistListV0->Add(fHistV0SelArmenteros);
+ }
+ }//V0
+
+ // Balance function histograms
+ // Initialize histograms if not done yet
+ if(!fBalance->GetHistNp()){
+ AliWarning("Histograms not yet initialized! --> Will be done now");
+ AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
+ fBalance->InitHistograms();
+ }
+
+ if(fRunShuffling) {
+ if(!fShuffledBalance->GetHistNp()) {
+ AliWarning("Histograms (shuffling) not yet initialized! --> Will be done now");
+ AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
+ fShuffledBalance->InitHistograms();
+ }
+ }
+
+ if(fRunMixing) {
+ if(!fMixedBalance->GetHistNp()) {
+ AliWarning("Histograms (mixing) not yet initialized! --> Will be done now");
+ AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
+ fMixedBalance->InitHistograms();
+ }
+ }
+
+ fListTriggeredBF->Add(fBalance->GetHistNp());
+ fListTriggeredBF->Add(fBalance->GetHistNn());
+ fListTriggeredBF->Add(fBalance->GetHistNpn());
+ fListTriggeredBF->Add(fBalance->GetHistNnn());
+ fListTriggeredBF->Add(fBalance->GetHistNpp());
+ fListTriggeredBF->Add(fBalance->GetHistNnp());
+
+ if(fRunShuffling) {
+ fListTriggeredBFS->Add(fShuffledBalance->GetHistNp());
+ fListTriggeredBFS->Add(fShuffledBalance->GetHistNn());
+ fListTriggeredBFS->Add(fShuffledBalance->GetHistNpn());
+ fListTriggeredBFS->Add(fShuffledBalance->GetHistNnn());
+ fListTriggeredBFS->Add(fShuffledBalance->GetHistNpp());
+ fListTriggeredBFS->Add(fShuffledBalance->GetHistNnp());
+ }
+
+ if(fRunMixing) {
+ fListTriggeredBFM->Add(fMixedBalance->GetHistNp());
+ fListTriggeredBFM->Add(fMixedBalance->GetHistNn());
+ fListTriggeredBFM->Add(fMixedBalance->GetHistNpn());
+ fListTriggeredBFM->Add(fMixedBalance->GetHistNnn());
+ fListTriggeredBFM->Add(fMixedBalance->GetHistNpp());
+ fListTriggeredBFM->Add(fMixedBalance->GetHistNnp());
+ }
+
+ // PID Response task active?
+ if(fRunV0) {
+ fPIDResponse = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetPIDResponse();
+ if (!fPIDResponse) AliFatal("This Task needs the PID response attached to the inputHandler");
+ }
+
+ // Event Mixing
+ Int_t trackDepth = fMixingTracks;
+ Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemented of AliEventPoolManager
+
+ Double_t centralityBins[] = {0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,15.,20.,25.,30.,35.,40.,45.,50.,55.,60.,65.,70.,75.,80.,90.,100.}; // SHOULD BE DEDUCED FROM CREATED ALITHN!!!
+ Double_t* centbins = centralityBins;
+ Int_t nCentralityBins = sizeof(centralityBins) / sizeof(Double_t) - 1;
+
+ // bins for second buffer are shifted by 100 cm
+ Double_t vertexBins[] = {-10., -7., -5., -3., -1., 1., 3., 5., 7., 10.}; // SHOULD BE DEDUCED FROM CREATED ALITHN!!!
+ Double_t* vtxbins = vertexBins;
+ Int_t nVertexBins = sizeof(vertexBins) / sizeof(Double_t) - 1;
+
+ fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, centbins, nVertexBins, vtxbins);
+
+
+ // Post output data.
+ PostData(1, fList);
+ PostData(2, fListTriggeredBF);
+ if(fRunShuffling) PostData(3, fListTriggeredBFS);
+ if(fRunMixing) PostData(4, fListTriggeredBFM);
+ if(fRunV0) PostData(5,fHistListV0);
+
+ TH1::AddDirectory(oldStatus);
+
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskTriggeredBF::UserExec(Option_t *) {
+ // Main loop
+ // Called for each event
+
+ TString gAnalysisLevel = fBalance->GetAnalysisLevel();
+ Float_t fCentrality = -1.;
+
+ // -------------------------------------------------------------
+ // AOD analysis (vertex and track cuts also here!!!!)
+ if(gAnalysisLevel == "AOD") {
+ AliVEvent* eventMain = dynamic_cast<AliVEvent*>(InputEvent());
+ if(!eventMain) {
+ AliError("eventMain not available");
+ return;
+ }
+
+ // check event cuts and fill event histograms
+ if((fCentrality = IsEventAccepted(eventMain)) < 0){
+ return;
+ }
+
+ // get the accepted tracks in main event
+ TObjArray *tracksMain = NULL;
+ if(fRunV0) tracksMain = GetAcceptedV0s(eventMain);
+ else tracksMain = GetAcceptedTracks(eventMain);
+
+ // store charges of all accepted tracks, shuffle and reassign (two extra loops!)
+ TObjArray* tracksShuffled = NULL;
+ if(fRunShuffling){
+ tracksShuffled = GetShuffledTracks(tracksMain);
+ }
+
+ // Event mixing --> UPDATE POOL IS MISSING!!!
+ if (fRunMixing)
+ {
+ // 1. First get an event pool corresponding in mult (cent) and
+ // zvertex to the current event. Once initialized, the pool
+ // should contain nMix (reduced) events. This routine does not
+ // pre-scan the chain. The first several events of every chain
+ // will be skipped until the needed pools are filled to the
+ // specified depth. If the pool categories are not too rare, this
+ // should not be a problem. If they are rare, you could lose`
+ // statistics.
+
+ // 2. Collect the whole pool's content of tracks into one TObjArray
+ // (bgTracks), which is effectively a single background super-event.
+
+ // 3. The reduced and bgTracks arrays must both be passed into
+ // FillCorrelations(). Also nMix should be passed in, so a weight
+ // of 1./nMix can be applied.
+
+ AliEventPool* pool = fPoolMgr->GetEventPool(fCentrality, eventMain->GetPrimaryVertex()->GetZ());
+
+ if (!pool){
+ AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCentrality, eventMain->GetPrimaryVertex()->GetZ()));
+ }
+ else{
+
+ //pool->SetDebug(1);
+
+ if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5){
+
+
+ Int_t nMix = pool->GetCurrentNEvents();
+ //cout << "nMix = " << nMix << " tracks in pool = " << pool->NTracksInPool() << endl;
+
+ //((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(2);
+ //((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());
+ //if (pool->IsReady())
+ //((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(3);
+
+ // Fill mixed-event histos here
+ for (Int_t jMix=0; jMix<nMix; jMix++)
+ {
+ TObjArray* tracksMixed = pool->GetEvent(jMix);
+ fMixedBalance->FillBalance(fCentrality,tracksMain,tracksMixed);
+ }
+ }
+
+ // Update the Event pool
+ pool->UpdatePool(tracksMain);
+ //pool->PrintInfo();
+
+ }//pool NULL check
+ }//run mixing
+
+ // calculate balance function
+ fBalance->FillBalance(fCentrality,tracksMain,NULL);
+
+ // calculate shuffled balance function
+ if(fRunShuffling && tracksShuffled != NULL) {
+ fShuffledBalance->FillBalance(fCentrality,tracksShuffled,NULL);
+ }
+
+ }//AOD analysis
+ else{
+ AliError("Triggered Balance Function analysis only for AODs!");
+ }
+}
+
+//________________________________________________________________________
+Float_t AliAnalysisTaskTriggeredBF::IsEventAccepted(AliVEvent *event){
+ // Checks the Event cuts
+ // Fills Event statistics histograms
+
+ // event selection done in AliAnalysisTaskSE::Exec() --> this is not used
+ fHistEventStats->Fill(1); //all events
+
+ Bool_t isSelectedMain = kTRUE;
+ Float_t fCentrality = -1.;
+ Int_t nV0s = event->GetNumberOfV0s();
+ TString gAnalysisLevel = fBalance->GetAnalysisLevel();
+
+ if(fUseOfflineTrigger)
+ isSelectedMain = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
+
+ //V0 QA histograms (before trigger selection)
+ if(fRunV0){
+ fHistMultiplicityBeforeTrigSel->Fill ( -1 );
+ fHistV0MultiplicityBeforeTrigSel->Fill ( nV0s );
+ }
+
+ if(isSelectedMain) {
+ fHistEventStats->Fill(2); //triggered events
+
+ //Centrality stuff
+ if(fUseCentrality) {
+ if(gAnalysisLevel == "AOD") { //centrality in AOD header
+ AliAODHeader *header = (AliAODHeader*) event->GetHeader();
+ fCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());
+
+ // QA for centrality estimators
+ fHistCentStats->Fill(0.,header->GetCentralityP()->GetCentralityPercentile("V0M"));
+ fHistCentStats->Fill(1.,header->GetCentralityP()->GetCentralityPercentile("FMD"));
+ fHistCentStats->Fill(2.,header->GetCentralityP()->GetCentralityPercentile("TRK"));
+ fHistCentStats->Fill(3.,header->GetCentralityP()->GetCentralityPercentile("TKL"));
+ fHistCentStats->Fill(4.,header->GetCentralityP()->GetCentralityPercentile("CL0"));
+ fHistCentStats->Fill(5.,header->GetCentralityP()->GetCentralityPercentile("CL1"));
+ fHistCentStats->Fill(6.,header->GetCentralityP()->GetCentralityPercentile("V0MvsFMD"));
+ fHistCentStats->Fill(7.,header->GetCentralityP()->GetCentralityPercentile("TKLvsV0M"));
+ fHistCentStats->Fill(8.,header->GetCentralityP()->GetCentralityPercentile("ZEMvsZDC"));
+
+ // centrality QA (V0M)
+ fHistV0M->Fill(event->GetVZEROData()->GetMTotV0A(), event->GetVZEROData()->GetMTotV0C());
+
+ // centrality QA (reference tracks)
+ fHistRefTracks->Fill(0.,header->GetRefMultiplicity());
+ fHistRefTracks->Fill(1.,header->GetRefMultiplicityPos());
+ fHistRefTracks->Fill(2.,header->GetRefMultiplicityNeg());
+ fHistRefTracks->Fill(3.,header->GetTPConlyRefMultiplicity());
+ fHistRefTracks->Fill(4.,header->GetNumberOfITSClusters(0));
+ fHistRefTracks->Fill(5.,header->GetNumberOfITSClusters(1));
+ fHistRefTracks->Fill(6.,header->GetNumberOfITSClusters(2));
+ fHistRefTracks->Fill(7.,header->GetNumberOfITSClusters(3));
+ fHistRefTracks->Fill(8.,header->GetNumberOfITSClusters(4));
+
+ //V0 QA histograms (after trigger selection)
+ if(fRunV0){
+ fHistMultiplicityForTrigEvt->Fill ( fCentrality );
+ fHistV0MultiplicityForTrigEvt->Fill ( nV0s );
+ }
+ }
+ }
+
+
+ const AliVVertex *vertex = event->GetPrimaryVertex();
+
+ if(vertex) {
+ Double32_t fCov[6];
+ vertex->GetCovarianceMatrix(fCov);
+ if(vertex->GetNContributors() > 0) {
+ if(fCov[5] != 0) {
+ fHistEventStats->Fill(3); //events with a proper vertex
+ if(TMath::Abs(vertex->GetX()) < fVxMax) {
+ if(TMath::Abs(vertex->GetY()) < fVyMax) {
+ if(TMath::Abs(vertex->GetZ()) < fVzMax) {
+ fHistEventStats->Fill(4); //analyzed events
+ fHistVx->Fill(vertex->GetX());
+ fHistVy->Fill(vertex->GetY());
+ fHistVz->Fill(vertex->GetZ());
+
+
+
+ //V0 QA histograms (vertex Z check)
+ if(fRunV0){
+ fHistV0MultiplicityForSelEvt ->Fill( nV0s );
+ fHistMultiplicity->Fill(fCentrality);
+
+ //V0 QA histograms (Only look at events with well-established PV)
+ const AliAODVertex *lPrimarySPDVtx = ((AliAODEvent*)event)->GetPrimaryVertexSPD();
+ if(lPrimarySPDVtx){
+ fHistMultiplicityNoTPCOnly->Fill ( fCentrality );
+ fHistV0MultiplicityForSelEvtNoTPCOnly->Fill ( nV0s );
+
+ //V0 QA histograms (Pileup Rejection)
+ // FIXME : quality selection regarding pile-up rejection
+ fHistMultiplicityNoTPCOnlyNoPileup->Fill(fCentrality);
+ fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup ->Fill( nV0s );
+
+ }
+ else{
+ return -1;
+ }
+ }
+
+
+ // take only events inside centrality class
+ if((fCentrality > fCentralityPercentileMin) && (fCentrality < fCentralityPercentileMax)){
+ return fCentrality;
+ }//centrality class
+ }//Vz cut
+ }//Vy cut
+ }//Vx cut
+ }//proper vertex resolution
+ }//proper number of contributors
+ }//vertex object valid
+ }//triggered event
+
+ // in all other cases return -1 (event not accepted)
+ return -1;
+}
+
+//________________________________________________________________________
+TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedTracks(AliVEvent *event){
+ // Returns TObjArray with tracks after all track cuts (only for AOD!)
+ // Fills QA histograms
+
+ //output TObjArray holding all good tracks
+ TObjArray* tracksAccepted = new TObjArray;
+ tracksAccepted->SetOwner(kTRUE);
+
+ Short_t vCharge = 0;
+ Double_t vEta = 0.;
+ Double_t vPhi = 0.;
+ Double_t vPt = 0.;
+
+ // Loop over tracks in event
+ for (Int_t iTracks = 0; iTracks < event->GetNumberOfTracks(); iTracks++) {
+ AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(event->GetTrack(iTracks));
+ if (!aodTrack) {
+ AliError(Form("Could not receive track %d", iTracks));
+ continue;
+ }
+
+ // AOD track cuts
+
+ // For ESD Filter Information: ANALYSIS/macros/AddTaskESDfilter.C
+ // take only TPC only tracks
+ fHistTrackStats->Fill(aodTrack->GetFilterMap());
+ if(!aodTrack->TestFilterBit(nAODtrackCutBit)) continue;
+
+ vCharge = aodTrack->Charge();
+ vEta = aodTrack->Eta();
+ vPhi = aodTrack->Phi() * TMath::RadToDeg();
+ vPt = aodTrack->Pt();
+
+ Float_t dcaXY = aodTrack->DCA(); // this is the DCA from global track (not exactly what is cut on)
+ Float_t dcaZ = aodTrack->ZAtDCA(); // this is the DCA from global track (not exactly what is cut on)
+
+
+ // Kinematics cuts from ESD track cuts
+ if( vPt < fPtMin || vPt > fPtMax) continue;
+ if( vEta < fEtaMin || vEta > fEtaMax) continue;
+
+ // Extra DCA cuts (for systematic studies [!= -1])
+ if( fDCAxyCut != -1 && fDCAzCut != -1){
+ if(TMath::Sqrt((dcaXY*dcaXY)/(fDCAxyCut*fDCAxyCut)+(dcaZ*dcaZ)/(fDCAzCut*fDCAzCut)) > 1 ){
+ continue; // 2D cut
+ }
+ }
+
+ // Extra TPC cuts (for systematic studies [!= -1])
+ if( fTPCchi2Cut != -1 && aodTrack->Chi2perNDF() > fTPCchi2Cut){
+ continue;
+ }
+ if( fNClustersTPCCut != -1 && aodTrack->GetTPCNcls() < fNClustersTPCCut){
+ continue;
+ }
+
+ // fill QA histograms
+ fHistClus->Fill(aodTrack->GetITSNcls(),aodTrack->GetTPCNcls());
+ fHistDCA->Fill(dcaZ,dcaXY);
+ fHistChi2->Fill(aodTrack->Chi2perNDF());
+ fHistPt->Fill(vPt);
+ fHistEta->Fill(vEta);
+ fHistPhi->Fill(vPhi);
+
+ // add the track to the TObjArray
+ tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge,1.));
+ }
+
+ return tracksAccepted;
+}
+
+//________________________________________________________________________
+TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedV0s(AliVEvent *event){
+ // Returns TObjArray with tracks after all track cuts (only for AOD!)
+ // Fills QA histograms
+
+ //output TObjArray holding all good tracks
+ TObjArray* tracksAccepted = new TObjArray;
+ tracksAccepted->SetOwner(kTRUE);
+
+ Short_t vCharge = 0;
+ Double_t vEta = 0.;
+ Double_t vPhi = 0.;
+ Double_t vPt = 0.;
+
+ //------------------------------------------------
+ // MAIN LAMBDA LOOP STARTS HERE (basically a copy of AliAnalysisTaskExtractV0AOD)
+ //------------------------------------------------
+
+ // parameters (for the time being hard coded here) --> from David for EbyE Lambdas
+ Bool_t fkUseOnTheFly = kFALSE;
+ Double_t fRapidityBoundary = 0.5;
+ Double_t fCutDaughterEta = 0.8;
+ Double_t fCutV0Radius = 0.9;
+ Double_t fCutDCANegToPV = 0.1;
+ Double_t fCutDCAPosToPV = 0.1;
+ Double_t fCutDCAV0Daughters = 1.0;
+ Double_t fCutV0CosPA = 0.9995;
+ Double_t fMassLambda = 1.115683;
+ Double_t fCutMassLambda = 0.007;
+ Double_t fCutProperLifetime = 3*7.9;
+ Double_t fCutLeastNumberOfCrossedRows = 70;
+ Double_t fCutLeastNumberOfCrossedRowsOverFindable = 0.8;
+ Double_t fCutTPCPIDNSigmasProton = 3.0;
+ Double_t fCutTPCPIDNSigmasPion = 5.0;
+
+
+ //Variable definition
+ Int_t lOnFlyStatus = 0;// nv0sOn = 0, nv0sOff = 0;
+ Double_t lChi2V0 = 0;
+ Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
+ Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
+ Double_t lV0CosineOfPointingAngle = 0;
+ Double_t lV0Radius = 0, lPt = 0;
+ Double_t lEta = 0, lPhi = 0;
+ Double_t lRap = 0, lRapK0Short = 0, lRapLambda = 0;
+ Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
+ Double_t lAlphaV0 = 0, lPtArmV0 = 0;
+
+ Double_t fMinV0Pt = 0;
+ Double_t fMaxV0Pt = 100;
+
+
+
+ // some event observables
+ Int_t nv0s = event->GetNumberOfV0s();
+ Double_t tPrimaryVtxPosition[3];
+ const AliVVertex *primaryVtx = event->GetPrimaryVertex();
+ tPrimaryVtxPosition[0] = primaryVtx->GetX();
+ tPrimaryVtxPosition[1] = primaryVtx->GetY();
+ tPrimaryVtxPosition[2] = primaryVtx->GetZ();
+
+
+ //loop over V0s
+ for (Int_t iV0 = 0; iV0 < nv0s; iV0++)
+ {// This is the begining of the V0 loop
+ AliAODv0 *v0 = ((AliAODEvent*)event)->GetV0(iV0);
+ if (!v0) continue;
+
+ //Obsolete at AOD level...
+ //---> Fix On-the-Fly candidates, count how many swapped
+ //if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ){
+ // fHistSwappedV0Counter -> Fill( 1 );
+ //}else{
+ // fHistSwappedV0Counter -> Fill( 0 );
+ //}
+ //if ( fkUseOnTheFly ) CheckChargeV0(v0);
+
+ Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0);
+ Double_t tV0mom[3];
+ v0->GetPxPyPz( tV0mom );
+ Double_t lV0TotalMomentum = TMath::Sqrt(
+ tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );
+
+ lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);
+ lPt = v0->Pt();
+ lEta = v0->Eta();
+ lPhi = v0->Phi()*TMath::RadToDeg();
+ lRapK0Short = v0->RapK0Short();
+ lRapLambda = v0->RapLambda();
+ lRap = lRapLambda;//v0->Y(); //FIXME!!!
+ if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;
+
+ //UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPosID());
+ //UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetPosID());
+
+ Double_t lMomPos[3]; //v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
+ Double_t lMomNeg[3]; //v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
+ lMomPos[0] = v0->MomPosX();
+ lMomPos[1] = v0->MomPosY();
+ lMomPos[2] = v0->MomPosZ();
+ lMomNeg[0] = v0->MomNegX();
+ lMomNeg[1] = v0->MomNegY();
+ lMomNeg[2] = v0->MomNegZ();
+
+ AliAODTrack *pTrack=(AliAODTrack *)v0->GetDaughter(0); //0->Positive Daughter
+ AliAODTrack *nTrack=(AliAODTrack *)v0->GetDaughter(1); //1->Negative Daughter
+ if (!pTrack || !nTrack) {
+ AliError("ERROR: Could not retreive one of the daughter track");
+ continue;
+ }
+
+ //Daughter Eta for Eta selection, afterwards
+ Double_t lNegEta = nTrack->Eta();
+ Double_t lPosEta = pTrack->Eta();
+
+ // Filter like-sign V0 (next: add counter and distribution)
+ if ( pTrack->Charge() == nTrack->Charge()){
+ continue;
+ }
+
+ //Quick test this far!
+
+
+ //________________________________________________________________________
+ // Track quality cuts
+ Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);
+ Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);
+ Float_t lLeastNbrCrossedRows = (lPosTrackCrossedRows>lNegTrackCrossedRows) ? lNegTrackCrossedRows : lPosTrackCrossedRows;
+
+ // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
+ if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
+ if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
+
+ if ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) ) continue;
+
+ //Findable clusters > 0 condition
+ if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;
+
+ //Compute ratio Crossed Rows / Findable clusters
+ //Note: above test avoids division by zero!
+ Float_t lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF()));
+ Float_t lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF()));
+ Float_t lLeastNbrCrossedRowsOverFindable = (lPosTrackCrossedRowsOverFindable>lNegTrackCrossedRowsOverFindable) ? lNegTrackCrossedRowsOverFindable : lPosTrackCrossedRowsOverFindable;
+
+ //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here
+ if ( lLeastNbrCrossedRowsOverFindable < 0.8) continue;
+
+ //End track Quality Cuts
+ //________________________________________________________________________
+
+
+ lDcaPosToPrimVertex = v0->DcaPosToPrimVertex();
+ lDcaNegToPrimVertex = v0->DcaNegToPrimVertex();
+
+ lOnFlyStatus = v0->GetOnFlyStatus();
+ lChi2V0 = v0->Chi2V0();
+ lDcaV0Daughters = v0->DcaV0Daughters();
+ lDcaV0ToPrimVertex = v0->DcaV0ToPrimVertex();
+ lV0CosineOfPointingAngle = v0->CosPointingAngle(tPrimaryVtxPosition);
+
+ // Distance over total momentum
+ Double_t lDistOverTotMom = TMath::Sqrt(
+ TMath::Power( tDecayVertexV0[0] - tPrimaryVtxPosition[0] , 2) +
+ TMath::Power( tDecayVertexV0[1] - tPrimaryVtxPosition[1] , 2) +
+ TMath::Power( tDecayVertexV0[2] - tPrimaryVtxPosition[2] , 2)
+ );
+ lDistOverTotMom /= (lV0TotalMomentum+1e-10); //avoid division by zero, to be sure
+
+
+ // Getting invariant mass infos directly from ESD
+ lInvMassK0s = v0->MassK0Short();
+ lInvMassLambda = v0->MassLambda();
+ lInvMassAntiLambda = v0->MassAntiLambda();
+ lAlphaV0 = v0->AlphaV0();
+ lPtArmV0 = v0->PtArmV0();
+
+ //Official means of acquiring N-sigmas
+ Double_t lNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );
+ Double_t lNSigmasPosPion = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );
+ Double_t lNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );
+ Double_t lNSigmasNegPion = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );
+
+ //V0 QA histograms (before V0 selection)
+ fHistV0InvMassK0->Fill(lInvMassK0s);
+ fHistV0InvMassLambda->Fill(lInvMassLambda);
+ fHistV0InvMassAntiLambda->Fill(lInvMassAntiLambda);
+ fHistV0Armenteros->Fill(lAlphaV0,lPtArmV0);
+
+
+ //First Selection: Reject OnFly
+ if( (lOnFlyStatus == 0 && fkUseOnTheFly == kFALSE) || (lOnFlyStatus != 0 && fkUseOnTheFly == kTRUE ) ){
+
+
+ //Second Selection: rough 20-sigma band, parametric.
+ //K0Short: Enough to parametrize peak broadening with linear function.
+ Double_t lUpperLimitK0Short = (5.63707e-01) + (1.14979e-02)*lPt;
+ Double_t lLowerLimitK0Short = (4.30006e-01) - (1.10029e-02)*lPt;
+
+ //Lambda: Linear (for higher pt) plus exponential (for low-pt broadening)
+ //[0]+[1]*x+[2]*TMath::Exp(-[3]*x)
+ Double_t lUpperLimitLambda = (1.13688e+00) + (5.27838e-03)*lPt + (8.42220e-02)*TMath::Exp(-(3.80595e+00)*lPt);
+ Double_t lLowerLimitLambda = (1.09501e+00) - (5.23272e-03)*lPt - (7.52690e-02)*TMath::Exp(-(3.46339e+00)*lPt);
+
+ //Do Selection
+ if( (lInvMassLambda < lUpperLimitLambda && lInvMassLambda > lLowerLimitLambda ) ||
+ (lInvMassAntiLambda < lUpperLimitLambda && lInvMassAntiLambda > lLowerLimitLambda ) ||
+ (lInvMassK0s < lUpperLimitK0Short && lInvMassK0s > lLowerLimitK0Short ) ){
+
+
+ // //Pre-selection in case this is AA...
+ // //if( fkIsNuclear == kFALSE ) fTree->Fill();
+ // //if( fkIsNuclear == kTRUE){
+ // //If this is a nuclear collision___________________
+ // // ... pre-filter with TPC, daughter eta selection
+
+
+ if( (lInvMassLambda < lUpperLimitLambda && lInvMassLambda > lLowerLimitLambda
+ && TMath::Abs(lNSigmasPosProton) < 6.0 && TMath::Abs(lNSigmasNegPion) < 6.0 ) ||
+ (lInvMassAntiLambda < lUpperLimitLambda && lInvMassAntiLambda > lLowerLimitLambda
+ && TMath::Abs(lNSigmasNegProton) < 6.0 && TMath::Abs(lNSigmasPosPion) < 6.0 ) ||
+ (lInvMassK0s < lUpperLimitK0Short && lInvMassK0s > lLowerLimitK0Short
+ && TMath::Abs(lNSigmasNegPion) < 6.0 && TMath::Abs(lNSigmasPosPion) < 6.0 ) ){
+
+ //insane test
+ if ( TMath::Abs(lNegEta)<0.8 && TMath::Abs(lPosEta)<0.8 ){
+
+ // start the fine selection (usually done in post processing, but we don't have time to waste) --> Lambdas!
+ if(
+ TMath::Abs(lRap)<fRapidityBoundary &&
+ TMath::Abs(lNegEta) <= fCutDaughterEta &&
+ TMath::Abs(lPosEta) <= fCutDaughterEta &&
+ lV0Radius >= fCutV0Radius &&
+ lDcaNegToPrimVertex >= fCutDCANegToPV &&
+ lDcaPosToPrimVertex >= fCutDCAPosToPV &&
+ lDcaV0Daughters <= fCutDCAV0Daughters &&
+ lV0CosineOfPointingAngle >= fCutV0CosPA &&
+ fMassLambda*lDistOverTotMom <= fCutProperLifetime &&
+ lLeastNbrCrossedRows >= fCutLeastNumberOfCrossedRows &&
+ lLeastNbrCrossedRowsOverFindable >= fCutLeastNumberOfCrossedRowsOverFindable &&
+ lPtArmV0 * 5 < TMath::Abs(lAlphaV0) &&
+ ((TMath::Abs(lNSigmasNegPion) <= fCutTPCPIDNSigmasPion &&
+ TMath::Abs(lNSigmasPosProton) <= fCutTPCPIDNSigmasProton) ||
+ (TMath::Abs(lNSigmasPosPion) <= fCutTPCPIDNSigmasPion &&
+ TMath::Abs(lNSigmasNegProton) <= fCutTPCPIDNSigmasProton))
+ )
+ {
+
+ //V0 QA histograms (after V0 selection)
+ fHistV0SelInvMassK0->Fill(lInvMassK0s);
+ fHistV0SelInvMassLambda->Fill(lInvMassLambda);
+ fHistV0SelInvMassAntiLambda->Fill(lInvMassAntiLambda);
+
+ // this means a V0 candidate is found
+ if(TMath::Abs(lInvMassLambda-fMassLambda) < fCutMassLambda ||
+ TMath::Abs(lInvMassAntiLambda-fMassLambda) < fCutMassLambda){
+
+ fHistV0SelArmenteros->Fill(lAlphaV0,lPtArmV0);
+
+ vEta = lEta;
+ vPhi = lPhi;
+ vPt = lPt;
+ if(lAlphaV0 > 0) vCharge = 1;
+ if(lAlphaV0 < 0) vCharge = -1;
+
+ // fill QA histograms
+ fHistPt->Fill(vPt);
+ fHistEta->Fill(vEta);
+ fHistPhi->Fill(vPhi);
+
+ // add the track to the TObjArray
+ tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge,1.));
+ }
+ }
+ }
+ }
+ //}//end nuclear_____________________________________
+ }
+ }
+ }//V0 loop
+
+ return tracksAccepted;
+}
+
+//________________________________________________________________________
+TObjArray* AliAnalysisTaskTriggeredBF::GetShuffledTracks(TObjArray *tracks){
+ // Clones TObjArray and returns it with tracks after shuffling the charges
+
+ TObjArray* tracksShuffled = new TObjArray;
+ tracksShuffled->SetOwner(kTRUE);
+
+ vector<Short_t> *chargeVector = new vector<Short_t>; //original charge of accepted tracks
+
+ for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
+ {
+ AliVParticle* track = (AliVParticle*) tracks->At(i);
+ chargeVector->push_back(track->Charge());
+ }
+
+ random_shuffle(chargeVector->begin(), chargeVector->end());
+
+ for(Int_t i = 0; i < tracks->GetEntriesFast(); i++){
+ AliVParticle* track = (AliVParticle*) tracks->At(i);
+ tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i),1.));
+ }
+
+ delete chargeVector;
+
+ return tracksShuffled;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskTriggeredBF::FinishTaskOutput(){
+ //checks if Balance Function objects are there (needed to write the histograms)
+ if (!fBalance) {
+ AliError("fBalance not available");
+ return;
+ }
+ if(fRunShuffling) {
+ if (!fShuffledBalance) {
+ AliError("fShuffledBalance not available");
+ return;
+ }
+ }
+
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskTriggeredBF::Terminate(Option_t *) {
+ // Called once at the end of the query
+
+ // not implemented ...
+
+}
+
+void AliAnalysisTaskTriggeredBF::UserExecMix(Option_t *)
+{
+
+ // not yet done for event mixing!
+ return;
+
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff