-/**************************************************************************\r
-* Copyright(c) 1998-2008,ALICE Experiment at CERN,All rights reserved. *\r
-* *\r
-* Author: The ALICE Off-line Project. *\r
-* Contributors are mentioned in the code where appropriate. *\r
-* *\r
-* Permission to use,copy,modify and distribute this software and its *\r
-* documentation strictly for non-commercial purposes is hereby granted *\r
-* without fee,provided that the above copyright notice appears in all *\r
-* copies and that both the copyright notice and this permission notice *\r
-* appear in the supporting documentation. The authors make no claims *\r
-* about the suitability of this software for any purpose. It is *\r
-* provided "as is" without express or implied warranty. *\r
-**************************************************************************/\r
-\r
-/////////////////////////////////////////////////////\r
-// AliAnalysisTaskFlowCascade:\r
-// Analysis task to select Xi and Omega candidates for flow analysis.\r
-//\r
-// Author: Zhong-Bao.Yin@cern.ch\r
-//////////////////////////////////////////////////////\r
-\r
-#include "TChain.h"\r
-#include "TList.h"\r
-#include "TH1D.h"\r
-#include "TH2D.h"\r
-#include "TH3D.h"\r
-#include "TProfile.h"\r
-#include "TVector3.h"\r
-\r
-#include "AliAnalysisTaskSE.h"\r
-#include "AliAnalysisManager.h"\r
-\r
-#include "AliVParticle.h"\r
-#include "AliESDEvent.h"\r
-#include "AliESDInputHandler.h"\r
-#include "AliESDv0.h"\r
-#include "AliESDcascade.h"\r
-#include "AliESDtrack.h"\r
-#include "AliESDtrackCuts.h"\r
-#include "AliCentrality.h"\r
-#include "AliVVertex.h"\r
-#include "AliESDVZERO.h"\r
-#include "AliESDUtils.h"\r
-\r
-#include "AliTPCPIDResponse.h"\r
-#include "AliTOFPIDResponse.h"\r
-#include "AliPIDResponse.h"\r
-\r
-#include "AliAODEvent.h"\r
-#include "AliAODInputHandler.h"\r
-#include "AliAODTrack.h"\r
-#include "AliAODVZERO.h"\r
-#include "AliAODcascade.h"\r
-\r
-#include "TMath.h"\r
-#include "TObjArray.h"\r
-#include "AliFlowCandidateTrack.h"\r
-\r
-#include "AliFlowTrackCuts.h"\r
-#include "AliFlowEventCuts.h"\r
-#include "AliFlowEvent.h"\r
-#include "AliFlowCommonConstants.h"\r
-\r
-#include "AliAnalysisTaskFlowCascade.h"\r
-\r
-ClassImp(AliAnalysisTaskFlowCascade)\r
-\r
-//_____________________________________________________________________________\r
- AliAnalysisTaskFlowCascade::AliAnalysisTaskFlowCascade() :\r
- AliAnalysisTaskSE(),\r
- // fMinCent(0), fMaxCent(0),\r
- fSpecie(0),\r
- fMassBins(0),\r
- fMinMass(0.0),\r
- fMaxMass(0.0),\r
- fCutsEvent(NULL),\r
- fCutsRPTPC(NULL),\r
- fCutsRPVZE(NULL),\r
- fCutsPOI(NULL),\r
- fCutsDau(NULL),\r
- fPIDResponse(NULL),\r
- fFlowEventTPC(NULL),\r
- fFlowEventVZE(NULL),\r
- fCandidates(NULL),\r
- fQAList(NULL)\r
-{\r
- //ctor \r
- for (Int_t i=0; i!=8; ++i)\r
- fCascadeCuts[i] = 0;\r
- \r
-}\r
-\r
-//_____________________________________________________________________________\r
-AliAnalysisTaskFlowCascade\r
-::AliAnalysisTaskFlowCascade(const char *name,\r
- AliFlowEventCuts *cutsEvent, \r
- AliFlowTrackCuts *cutsRPTPC,\r
- AliFlowTrackCuts *cutsRPVZE,\r
- /* AliESDtrackCuts */ AliFlowTrackCuts *cutsDau ) :\r
- AliAnalysisTaskSE(name),\r
- //fMinCent(minCent), fMaxCent(maxCent),\r
- fSpecie(0),\r
- fMassBins(0),\r
- fMinMass(0.0),\r
- fMaxMass(0.0),\r
- fCutsEvent(cutsEvent),\r
- fCutsRPTPC(cutsRPTPC),\r
- fCutsRPVZE(cutsRPVZE),\r
- fCutsPOI(NULL),\r
- fCutsDau(cutsDau),\r
- fPIDResponse(NULL),\r
- fFlowEventTPC(NULL),\r
- fFlowEventVZE(NULL),\r
- fCandidates(NULL),\r
- fQAList(NULL)\r
-{\r
- //ctor \r
- for (Int_t i=0; i!=8; ++i)\r
- fCascadeCuts[i] = 0;\r
-\r
- DefineInput( 0,TChain::Class());\r
- DefineOutput(1,AliFlowEventSimple::Class()); // TPC object\r
- DefineOutput(2,AliFlowEventSimple::Class()); // VZE object\r
- DefineOutput(3,TList::Class());\r
-}\r
-\r
-//_____________________________________________________________________________\r
-AliAnalysisTaskFlowCascade::~AliAnalysisTaskFlowCascade()\r
-{\r
- if(fQAList) delete fQAList;\r
- if (fFlowEventTPC) delete fFlowEventTPC;\r
- if (fFlowEventVZE) delete fFlowEventVZE;\r
- if (fCandidates) delete fCandidates;\r
- if (fCutsDau) delete fCutsDau;\r
- if (fCutsPOI) delete fCutsPOI;\r
- if (fCutsRPTPC) delete fCutsRPTPC;\r
- if (fCutsRPVZE) delete fCutsRPVZE;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::UserCreateOutputObjects()\r
-{\r
-\r
- AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();\r
- AliInputEventHandler* inputHandler\r
- = (AliInputEventHandler*) (man->GetInputEventHandler());\r
- fPIDResponse = inputHandler->GetPIDResponse();\r
- \r
- fQAList = new TList();\r
- fQAList->SetOwner();\r
- AddQAEvents();\r
- AddQACandidates();\r
- // PostData(3,fQAList);\r
-\r
- AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();\r
- cc->SetNbinsMult(1);\r
- cc->SetMultMin(0);\r
- cc->SetMultMax(1);\r
-\r
- cc->SetNbinsPt(20);\r
- cc->SetPtMin(0.0);\r
- cc->SetPtMax(10.0);\r
-\r
- cc->SetNbinsPhi(1);\r
- cc->SetPhiMin(0.0);\r
- cc->SetPhiMax(TMath::TwoPi());\r
-\r
- cc->SetNbinsEta(1);\r
- cc->SetEtaMin(-2.0);\r
- cc->SetEtaMax(+2.0);\r
-\r
- cc->SetNbinsQ(3);\r
- cc->SetQMin(0.0);\r
- cc->SetQMax(3.0);\r
-\r
- cc->SetNbinsMass(fMassBins);\r
- cc->SetMassMin(fMinMass);\r
- cc->SetMassMax(fMaxMass);\r
-\r
- fCutsPOI = new AliFlowTrackCuts("null_cuts");\r
- fCutsPOI->SetParamType(fCutsRPTPC->GetParamType());\r
- fCutsPOI->SetPtRange(+1,-1); // select nothing QUICK \r
- fCutsPOI->SetEtaRange(+1,-1); // select nothing VZERO \r
-\r
- fFlowEventTPC = new AliFlowEvent(3000);\r
- fFlowEventVZE = new AliFlowEvent(1000);\r
- fCandidates = new TObjArray(100);\r
- fCandidates->SetOwner();\r
-\r
- PostData(1,fFlowEventTPC);\r
- PostData(2,fFlowEventVZE);\r
- PostData(3,fQAList);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::AddQAEvents()\r
-{\r
- TList *tQAEvents = new TList();\r
- tQAEvents->SetName("Events");\r
- tQAEvents->SetOwner();\r
- TH1I* tEvent = new TH1I("Event","Number of Events", 3,0,3);\r
- tQAEvents->Add(tEvent);\r
- \r
- TH1D *tTPCRFP = new TH1D("RFPTPC",\r
- "TPC Reference Flow Particles;multiplicity",\r
- 100, 0, 3000); \r
- tQAEvents->Add(tTPCRFP);\r
- TH1D *tVZERFP = new TH1D("RFPVZE", \r
- "VZERO Reference Flow Particles;multiplicity",\r
- 100, 0, 30000); \r
- tQAEvents->Add(tVZERFP);\r
-\r
- TProfile *tCuts = new TProfile("Cuts","Analysis Cuts",10,0,10);\r
- tCuts->Fill(0.5,fCascadeCuts[0],1); \r
- tCuts->GetXaxis()->SetBinLabel(1,"dcaXiDau");\r
- tCuts->Fill(1.5,fCascadeCuts[1],1); \r
- tCuts->GetXaxis()->SetBinLabel(2,"XiCPA");\r
- tCuts->Fill(2.5,fCascadeCuts[2],1); \r
- tCuts->GetXaxis()->SetBinLabel(3,"dcaV0Vtx");\r
- tCuts->Fill(3.5,fCascadeCuts[3],1); \r
- tCuts->GetXaxis()->SetBinLabel(4,"dcaBachVtx");\r
- tCuts->Fill(4.5,fCascadeCuts[4],1); \r
- tCuts->GetXaxis()->SetBinLabel(5,"dcaV0Dau");\r
- tCuts->Fill(5.5,fCascadeCuts[5],1); \r
- tCuts->GetXaxis()->SetBinLabel(6,"V0CPA");\r
- tCuts->Fill(6.5,fCascadeCuts[6],1); \r
- tCuts->GetXaxis()->SetBinLabel(7,"dcaV0DauVtx");\r
- tCuts->Fill(7.5,fCascadeCuts[7],1); \r
- tCuts->GetXaxis()->SetBinLabel(8,"V0Mass");\r
- tQAEvents->Add(tCuts);\r
-\r
- fQAList->Add(tQAEvents);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::AddQACandidates()\r
-{\r
- TList *tQACandidates;\r
- \r
- tQACandidates = new TList();\r
- tQACandidates->SetOwner();\r
- tQACandidates->SetName("Candidates");\r
-\r
- TH1F* tChi2Xi = new TH1F("Chi2Xi", \r
- "Cascade #chi^{2}; #chi^{2}; Number of Cascades", \r
- 160, 0, 160);\r
- tQACandidates->Add(tChi2Xi);\r
-\r
- TH1F* tDCAXiDaughters \r
- = new TH1F( "DcaXiDaughters", \r
- "DCA between Xi Daughters; DCA (cm); Number of Cascades", \r
- 100, 0., 0.5);\r
- tQACandidates->Add(tDCAXiDaughters);\r
-\r
- TH1F * tDCABachToPrimVertex\r
- = new TH1F("DcaBachToPrimVertex", \r
- "DCA of Bach. to Prim. Vertex; DCA (cm);Number of Cascades", \r
- 250, 0., 2.5);\r
- tQACandidates->Add(tDCABachToPrimVertex);\r
- \r
- TH1F * tXiCosOfPointingAngle\r
- = new TH1F("XiCosOfPointingAngle",\r
- "Cos of Xi Pointing Angle; Cos (Xi Point.Angl);Number of Xis", \r
- 200, 0.99, 1.0);\r
- tQACandidates->Add(tXiCosOfPointingAngle);\r
-\r
- TH1F * tXiRadius = new TH1F("XiRadius", \r
- "Casc. decay transv. radius; r (cm); Counts" , \r
- 1050, 0., 105.0 );\r
- tQACandidates->Add(tXiRadius);\r
- \r
- TH1F *tMassLambda \r
- = new TH1F("MassLambdaAsCascDghter",\r
- "#Lambda assoc. to Casc. candidates; Eff. Mass (GeV/c^{2}); Counts", \r
- 300,1.00,1.3);\r
- tQACandidates->Add(tMassLambda);\r
-\r
- TH1F *tV0Chi2 = new TH1F("V0Chi2Xi", \r
- "V0 #chi^{2}, in cascade; #chi^{2};Counts", \r
- 160, 0, 40);\r
- tQACandidates->Add(tV0Chi2);\r
- \r
- TH1F * tV0CosOfPointingAngle \r
- = new TH1F("V0CosOfPointingAngleXi", \r
- "Cos of V0 Pointing Angle, in cascade;Cos(V0 Point. Angl); Counts", \r
- 200, 0.98, 1.0);\r
- tQACandidates->Add(tV0CosOfPointingAngle);\r
-\r
- TH1F *tV0Radius = new TH1F("V0RadiusXi", \r
- "V0 decay radius, in cascade; radius (cm); Counts", \r
- 1050, 0., 105.0);\r
- tQACandidates->Add(tV0Radius);\r
- \r
- TH1F * tDcaV0DaughtersXi \r
- = new TH1F("DcaV0DaughtersXi", \r
- "DCA between V0 daughters, in cascade;DCA (cm);Number of V0s", \r
- 120, 0., 0.6);\r
- tQACandidates->Add(tDcaV0DaughtersXi);\r
-\r
- TH1F * tDcaV0ToPrimVertex \r
- = new TH1F("DcaV0ToPrimVertexXi", \r
- "DCA of V0 to Prim. Vertex, in cascade;DCA (cm);Number of Cascades", 200, 0., 1.);\r
- tQACandidates->Add(tDcaV0ToPrimVertex);\r
-\r
- TH1F * tDCAPosToPrimVertex =\r
- new TH1F("DcaPosToPrimVertexXi", \r
- "DCA of V0 pos daughter to Prim. Vertex;DCA (cm);Counts", \r
- 300, 0, 3);\r
- tQACandidates->Add(tDCAPosToPrimVertex);\r
-\r
- TH1F * tDCANegToPrimVertex \r
- = new TH1F("DcaNegToPrimVertexXi", \r
- "DCA of V0 neg daughter to Prim. Vertex;DCA (cm);Counts", \r
- 300, 0, 3);\r
- tQACandidates->Add(tDCANegToPrimVertex);\r
-\r
- TH1F *tV0toXiCosOfPointingAngle \r
- = new TH1F("V0toXiCosOfPointingAngle", \r
- "Cos. of V0 Ptng Angl Xi vtx; Cos(V0 Point. Angl / Xi vtx); Counts", \r
- 100, 0.99, 1.0);\r
- tQACandidates->Add(tV0toXiCosOfPointingAngle);\r
-\r
- TH2F *th2Armenteros \r
- = new TH2F("Armenteros", \r
- "#alpha_{Arm}(casc. cand.) Vs Pt_{Arm}(casc. cand.); #alpha_{Arm} ; Pt_{Arm} (GeV/c)", \r
- 140, -1.2, 1.2, 300, 0., 0.3);\r
- tQACandidates->Add(th2Armenteros);\r
-\r
- TH2F *th2TPCdEdxOfCascDghters\r
- = new TH2F( "TPCdEdxOfCascDghters",\r
- "TPC dE/dx of the cascade daughters; charge x || #vec{p}_{TPC inner wall}(Casc. daughter) || (GeV/c); TPC signal (ADC) ", \r
- 200, -10.0, 10.0, 450, 0., 900.);\r
- tQACandidates->Add(th2TPCdEdxOfCascDghters);\r
-\r
- TH2F *th2MassVsPtAll \r
- = new TH2F("MassVsPtAll", \r
- "M_{candidates} vs Pt; Pt (GeV/c); M (GeV/c^{2})", \r
- 100, 0., 10., fMassBins, fMinMass, fMaxMass);\r
- tQACandidates->Add(th2MassVsPtAll);\r
-\r
- TH1F *tDistToVtxZAfter \r
- = new TH1F("DistToVtxZAfter", \r
- "Distance to vtx z after propagation to vtx; z [cm]", \r
- 100, -5., 5.);\r
- tQACandidates->Add(tDistToVtxZAfter);\r
-\r
- TH1F * tDistToVtxXYAfter \r
- = new TH1F("DistToVtxXYAfter", \r
- "Distance to vtx xy after propagation to vtx",\r
- 500, 0., 50.);\r
- tQACandidates->Add(tDistToVtxXYAfter);\r
-\r
- TH2F *th2DistToVtxZBeforeVsAfter \r
- = new TH2F("DistToVtxZBeforeVsAfter", \r
- "Distance to vtx z before vs after propagation; Distance before [cm]; Distance after [cm]", \r
- 500, -50., 50., 100, -5., 5.);\r
- tQACandidates->Add(th2DistToVtxZBeforeVsAfter);\r
-\r
- TH2F *th2DistToVtxXYBeforeVsAfter\r
- = new TH2F("DistToVtxXYBeforeVsAfter",\r
- "Distance to vtx xy before vs after propagation; Distance before [cm]; Distance after [cm]", \r
- 500, 0., 50, 500, 0., 50);\r
- tQACandidates->Add(th2DistToVtxXYBeforeVsAfter);\r
-\r
- TH2F * th2PxBeforeVsAfter \r
- = new TH2F("PxBeforeVsAfter", \r
- "Px before vs after propagation; Px [GeV/c]; Px' [GeV/c]", \r
- 200, -10., 10, 200, -10., 10.);\r
- tQACandidates->Add(th2PxBeforeVsAfter);\r
-\r
- TH2F * th2PyBeforeVsAfter\r
- = new TH2F("PyBeforeVsAfter",\r
- "Py before vs after propagation; Py [GeV/c]; Py' [GeV/c]",\r
- 200, -10., 10, 200, -10., 10.);\r
- tQACandidates->Add(th2PyBeforeVsAfter);\r
- \r
- TH2F * th2PhiPosBeforeVsAfter\r
- = new TH2F("PhiPosBeforeVsAfter", \r
- "Phi for positively charged candidates before vs after propagation; #phi; #phi'", \r
- 360, 0., 2.0*TMath::Pi(), 360, 0., 2.0*TMath::Pi());\r
- tQACandidates->Add(th2PhiPosBeforeVsAfter);\r
-\r
- TH2F *th2PhiNegBeforeVsAfter\r
- = new TH2F("PhiNegBeforeVsAfter",\r
- "Phi for negatively charged candidates before vs after propagation; #phi; #phi'",\r
- 360, 0., 2.0*TMath::Pi(), 360, 0., 2.0*TMath::Pi());\r
- tQACandidates->Add(th2PhiNegBeforeVsAfter);\r
-\r
- fQAList->Add(tQACandidates);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::NotifyRun()\r
-{\r
-}\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::UserExec(Option_t *)\r
-{\r
- AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());\r
- AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());\r
- Bool_t acceptEvent=kFALSE; \r
- fCandidates->SetLast(-1);\r
-\r
- if(fESD) {\r
- // recorrecting VZERO (for pass 1 only)\r
- /*\r
- Float_t *vChCorr = new Float_t[64];\r
- Float_t dummy;\r
- AliESDUtils::GetCorrV0(fESD,dummy,NULL,vChCorr);\r
- AliESDVZERO *vzero = (AliESDVZERO*) fESD->GetVZEROData();\r
- vzero->SetMultiplicity( vChCorr );\r
- delete [] vChCorr;\r
- */\r
- //\r
-\r
- ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(0);\r
- \r
- const AliVVertex *vtxGlb = fESD->GetPrimaryVertexTracks();\r
- const AliVVertex *vtxSPD = fESD->GetPrimaryVertexSPD();\r
- if(!vtxGlb || !vtxSPD) return;\r
- if( fCutsEvent->IsSelected(fESD, 0) && (TMath::Abs(vtxSPD->GetZ()-vtxGlb->GetZ()) <= 0.5) ) {\r
- acceptEvent = kTRUE;\r
- ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);\r
- ReadFromESDv0(fESD);\r
- }\r
- } else if(fAOD) {\r
- const AliVVertex *vtxGlb = fAOD->GetPrimaryVertex();\r
- const AliVVertex *vtxSPD = fAOD->GetPrimaryVertexSPD();\r
- if(!vtxGlb || !vtxSPD) return;\r
-\r
- ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(0);\r
- \r
- if(fCutsEvent->IsSelected(fAOD, 0) && (TMath::Abs(vtxSPD->GetZ()-vtxGlb->GetZ()) <= 0.5) ) {\r
- acceptEvent = kTRUE;\r
- ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);\r
- ReadFromAODv0(fAOD);\r
- }\r
-\r
- \r
- /*\r
-\r
- AliAODHeader *aodHeader = fAOD->GetHeader();\r
- if(!aodHeader) return;\r
- AliCentrality *centrality = aodHeader->GetCentralityP();\r
- if(!centrality) return;\r
- Double_t cent = centrality->GetCentralityPercentile("V0M" );\r
- Double_t cent1 = centrality->GetCentralityPercentile("TRK" );\r
- if(TMath::Abs(cent-cent1) >= 5.) return;\r
- \r
- if(cent<fMinCent||cent>=fMaxCent) return; //centrality cut\r
- \r
- Double_t zvtx = fAOD->GetPrimaryVertex()->GetZ();\r
- if(TMath::Abs(zvtx)>10.) return; //vertex cut\r
- \r
- ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);\r
- ReadFromAODv0(fAOD);\r
- */ \r
- }\r
- \r
- if(!acceptEvent) return;\r
-\r
- ((TH1D*)((TList*)fQAList->FindObject("Events"))->FindObject("RFPTPC"))\r
- ->Fill(fFlowEventTPC->GetNumberOfRPs() );\r
- Double_t mult=0;\r
- for(Int_t i=0; i != fFlowEventVZE->GetNumberOfRPs(); ++i) {\r
- AliFlowTrackSimple *pTrack = fFlowEventVZE->GetTrack(i);\r
- mult += pTrack->Weight();\r
- }\r
- ((TH1D*)((TList*)fQAList->FindObject("Events"))->FindObject("RFPVZE"))\r
- ->Fill( mult );\r
-\r
- // if(fDebug) printf("TPCevent %d | VZEevent %d\n",\r
- // fFlowEventTPC->NumberOfTracks(),\r
- // fFlowEventVZE->NumberOfTracks() );\r
- AddCandidates();\r
-\r
- PostData(1,fFlowEventTPC);\r
- PostData(2,fFlowEventVZE);\r
- PostData(3,fQAList);\r
-\r
- return;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::AddCandidates(){\r
- \r
- // if(fDebug) printf("I received %d candidates\n",\r
- // fCandidates->GetEntriesFast());\r
- for(int iCand=0; iCand!=fCandidates->GetEntriesFast(); ++iCand ) {\r
- AliFlowCandidateTrack *cand \r
- = dynamic_cast<AliFlowCandidateTrack*>(fCandidates->At(iCand));\r
- if(!cand) continue;\r
- // if(fDebug) \r
- // printf(" >Checking at candidate %d with %d daughters: mass %f\n",\r
- // iCand, cand->GetNDaughters(), cand->Mass());\r
- \r
- // untagging ===> \r
- for(int iDau=0; iDau != cand->GetNDaughters(); ++iDau) {\r
- // if(fDebug) \r
- // printf(" >Daughter %d with fID %d", iDau, cand->GetIDDaughter(iDau));\r
- for(int iRPs=0; iRPs != fFlowEventTPC->NumberOfTracks(); ++iRPs ) {\r
- AliFlowTrack *iRP \r
- = dynamic_cast<AliFlowTrack*>(fFlowEventTPC->GetTrack( iRPs ));\r
- if (!iRP) continue;\r
- if( !iRP->InRPSelection() ) continue;\r
- if( cand->GetIDDaughter(iDau) == iRP->GetID() ) {\r
- //if(fDebug) printf(" was in RP set");\r
- iRP->SetForRPSelection(kFALSE);\r
- fFlowEventTPC->SetNumberOfRPs( fFlowEventTPC->GetNumberOfRPs() -1 );\r
- }\r
- }\r
- //if(fDebug) printf("\n");\r
- }\r
- // <=== untagging\r
- cand->SetForPOISelection(kTRUE);\r
- fFlowEventTPC->InsertTrack( ((AliFlowTrack*) cand) );\r
- fFlowEventVZE->InsertTrack( ((AliFlowTrack*) cand) );\r
- }\r
-\r
- // if(fDebug) printf("TPCevent %d | VZEevent %d\n",\r
- // fFlowEventTPC->NumberOfTracks(),\r
- // fFlowEventVZE->NumberOfTracks() );\r
-\r
-}\r
-\r
-//______________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::ReadFromESDv0(AliESDEvent *fESD)\r
-{\r
- //AliFlowTrackCuts* cutsPOI = new AliFlowTrackCuts("null_cuts");\r
- //cutsPOI->SetParamType( fCutsRP->GetParamType() );\r
- //cutsPOI->SetParamType( AliFlowTrackCuts::kGlobal );\r
- // cutsPOI->SetPtRange(+1,-1); // select nothing\r
- //cutsPOI->SetEtaRange(+1,-1); // select nothing VZERO\r
-\r
- fCutsRPTPC->SetEvent(fESD,MCEvent());\r
- fCutsRPVZE->SetEvent(fESD,MCEvent());\r
-\r
- fCutsPOI->SetEvent(fESD,MCEvent());\r
-\r
- fFlowEventTPC->Fill(fCutsRPTPC,fCutsPOI);\r
- fFlowEventVZE->Fill(fCutsRPVZE,fCutsPOI);\r
- \r
- Double_t trkPrimaryVtxPos[3] = {-100., -100., -100.};\r
- Double_t bestPrimaryVtxPos[3] = {-100., -100., -100.};\r
- int nCascades=fESD->GetNumberOfCascades();\r
- \r
- const AliESDVertex *primaryTrackingESDVtx = fESD->GetPrimaryVertexTracks();\r
- primaryTrackingESDVtx->GetXYZ(trkPrimaryVtxPos);\r
- \r
- const AliESDVertex *primaryBestESDVtx = fESD->GetPrimaryVertex();\r
- primaryBestESDVtx->GetXYZ(bestPrimaryVtxPos);\r
-\r
- Double_t b = fESD->GetMagneticField();\r
-\r
- for(int i = 0; i != nCascades; ++i) {\r
- \r
- // Double_t trkPrimaryVtxRadius3D = -500.;\r
- // Double_t bestPrimaryVtxRadius3D = -500.;\r
- Double_t effMassXi = 0.;\r
- Double_t chi2Xi = -1.;\r
- Double_t dcaXiDaughters = -1.;\r
- Double_t XiCosOfPointingAngle = -1.;\r
- Double_t posXi[3] = {-1000., -1000., -1000.};\r
- Double_t XiRadius = -1000.;\r
-\r
- // Double_t innerWallMomCascDghters[3] = {-100., -100., -100.};\r
- //Double_t tpcSignalCascDghters[3] = {-100., -100., -100.};\r
-\r
- Double_t invMassLambdaAsCascDghter = 0.;\r
- Double_t V0Chi2Xi = -1.;\r
- Double_t dcaV0DaughtersXi = -1.;\r
- \r
- Double_t dcaBachToPrimaryVtxXi = -1.;\r
- Double_t dcaV0ToPrimaryVtxXi = -1.;\r
- Double_t dcaPosToPrimaryVtxXi = -1.;\r
- Double_t dcaNegToPrimaryVtxXi = -1.;\r
- Double_t V0CosOfPointingAngleXi = -1.;\r
- Double_t posV0Xi[3] = {-1000., -1000., -1000.};\r
- Double_t V0RadiusXi = -1000.;\r
- Double_t V0quality = 0.;\r
-\r
- Double_t invMassXiMinus = 0.;\r
- Double_t invMassXiPlus = 0.;\r
- Double_t invMassOmegaMinus = 0.;\r
- Double_t invMassOmegaPlus = 0.;\r
-\r
- /*\r
- Bool_t isPosInXiProton = kFALSE;\r
- Bool_t isPosInXiPion = kFALSE;\r
- Bool_t isPosInOmegaProton = kFALSE;\r
- Bool_t isPosInOmegaPion = kFALSE;\r
- \r
- Bool_t isNegInXiProton = kFALSE;\r
- Bool_t isNegInXiPion = kFALSE;\r
- Bool_t isNegInOmegaProton = kFALSE;\r
- Bool_t isNegInOmegaPion = kFALSE;\r
-\r
- Bool_t isBachelorKaon = kFALSE;\r
- Bool_t isBachelorPion = kFALSE;\r
- */\r
-\r
- Bool_t isBachelorKaonForTPC = kFALSE;\r
- Bool_t isBachelorPionForTPC = kFALSE;\r
- Bool_t isNegPionForTPC = kFALSE;\r
- Bool_t isPosPionForTPC = kFALSE;\r
- Bool_t isNegProtonForTPC = kFALSE;\r
- Bool_t isPosProtonForTPC = kFALSE;\r
-\r
- Double_t XiPx = 0., XiPy = 0., XiPz = 0.;\r
- Double_t XiPt = 0.;\r
- Double_t XiPtot = 0.;\r
- \r
- Double_t bachPx = 0., bachPy = 0., bachPz = 0.;\r
- Double_t bachPt = 0.;\r
- Double_t bachPtot = 0.;\r
- \r
- //Short_t chargeXi = -2;\r
- Double_t V0toXiCosOfPointingAngle = 0.;\r
- \r
- Double_t rapXi = -20.;\r
- Double_t rapOmega = -20.;\r
- Double_t phi = 6.3;\r
- Double_t alphaXi = -200.;\r
- Double_t ptArmXi = -200.;\r
- // TLorentzVector lv1, lv2, lv3, lv12, lvXi;\r
-\r
- Double_t distToVtxZBefore = -999.;\r
- Double_t distToVtxZAfter = -999.;\r
- Double_t distToVtxXYBefore = -999.;\r
- Double_t distToVtxXYAfter = -999.;\r
- Double_t XiPAfter[3] = {-999., -999., -999.};\r
- Double_t phiAfter = -999.;\r
- \r
- AliESDcascade *xi = fESD->GetCascade(i);\r
- if(!xi) continue;\r
- \r
- if(xi->Charge()<0)\r
- xi->ChangeMassHypothesis(V0quality, 3312); // Xi- hypothesis\r
- else if(xi->Charge() > 0)\r
- xi->ChangeMassHypothesis(V0quality, -3312);\r
- else continue;\r
-\r
- effMassXi = xi->GetEffMassXi();\r
- chi2Xi = xi->GetChi2Xi();\r
- dcaXiDaughters = xi->GetDcaXiDaughters();\r
- XiCosOfPointingAngle \r
- = xi->GetCascadeCosineOfPointingAngle(bestPrimaryVtxPos[0],\r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]);\r
- xi->GetXYZcascade(posXi[0], posXi[1], posXi[2]);\r
- XiRadius = TMath::Sqrt(posXi[0]*posXi[0]\r
- +posXi[1]*posXi[1]\r
- +posXi[2]*posXi[2]);\r
- \r
- UInt_t idxPosXi = (UInt_t)TMath::Abs(xi->GetPindex());\r
- UInt_t idxNegXi = (UInt_t)TMath::Abs(xi->GetNindex());\r
- UInt_t idxBach = (UInt_t)TMath::Abs(xi->GetBindex());\r
-\r
- if(idxBach == idxPosXi || idxBach == idxNegXi) continue;\r
-\r
- AliESDtrack *pTrkXi = fESD->GetTrack(idxPosXi);\r
- AliESDtrack *nTrkXi = fESD->GetTrack(idxNegXi);\r
- AliESDtrack *bTrkXi = fESD->GetTrack(idxBach);\r
-\r
- if( !pTrkXi || !nTrkXi || !bTrkXi ) continue;\r
-\r
- if( !fCutsDau->IsSelected(pTrkXi) \r
- || !fCutsDau->IsSelected(nTrkXi)\r
- || !fCutsDau->IsSelected(bTrkXi) ) continue;\r
-\r
- const AliExternalTrackParam *pExtTrk = pTrkXi->GetInnerParam();\r
- const AliExternalTrackParam *nExtTrk = nTrkXi->GetInnerParam();\r
- const AliExternalTrackParam *bExtTrk = bTrkXi->GetInnerParam();\r
- \r
- if(pExtTrk && pTrkXi->IsOn(AliESDtrack::kTPCin)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(pExtTrk->GetP()*pExtTrk->Charge(), pTrkXi->GetTPCsignal());\r
- }\r
- if(nExtTrk && nTrkXi->IsOn(AliESDtrack::kTPCin)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(nExtTrk->GetP()*nExtTrk->Charge(), nTrkXi->GetTPCsignal());\r
- }\r
- if(bExtTrk && bTrkXi->IsOn(AliESDtrack::kTPCin)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(bExtTrk->GetP()*bExtTrk->Charge(), bTrkXi->GetTPCsignal());\r
- }\r
- \r
- invMassLambdaAsCascDghter = xi->GetEffMass(); // from V0\r
- dcaV0DaughtersXi = xi->GetDcaV0Daughters();\r
- V0Chi2Xi = xi->GetChi2V0();\r
- V0CosOfPointingAngleXi \r
- = xi->GetV0CosineOfPointingAngle(bestPrimaryVtxPos[0],\r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]);\r
- dcaV0ToPrimaryVtxXi = xi->GetD(bestPrimaryVtxPos[0], \r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]);\r
- dcaBachToPrimaryVtxXi = TMath::Abs(bTrkXi->GetD(bestPrimaryVtxPos[0],\r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]));\r
- \r
- //V0\r
- xi->GetXYZ(posV0Xi[0], posV0Xi[1], posV0Xi[2]);\r
- V0RadiusXi = TMath::Sqrt(posV0Xi[0]*posV0Xi[0]\r
- +posV0Xi[1]*posV0Xi[1]\r
- +posV0Xi[2]*posV0Xi[2]);\r
- dcaPosToPrimaryVtxXi = TMath::Abs(pTrkXi->GetD(bestPrimaryVtxPos[0],\r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]));\r
- dcaNegToPrimaryVtxXi = TMath::Abs(nTrkXi->GetD(bestPrimaryVtxPos[0],\r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]));\r
-\r
- //apply cuts\r
- //if(XiRadius < 0.9 || XiRadius > 100.) continue;\r
- //if(dcaXiDaughters > 0.2) continue;\r
- //if(XiCosOfPointingAngle < 0.99) continue;\r
- //if(dcaV0ToPrimaryVtxXi < 0.03) continue;\r
- //if(dcaBachToPrimaryVtxXi < 0.01) continue;\r
- \r
- if(dcaXiDaughters > fCascadeCuts[0]) continue;\r
- if(XiCosOfPointingAngle < fCascadeCuts[1]) continue;\r
- if(dcaV0ToPrimaryVtxXi < fCascadeCuts[2]) continue;\r
- if(dcaBachToPrimaryVtxXi < fCascadeCuts[3]) continue;\r
-\r
- //V0 mass cut?\r
- // if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > 0.01) continue;\r
- if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > fCascadeCuts[7]) \r
- continue;\r
- \r
- //if(dcaV0DaughtersXi > 1.) continue;\r
- //if(V0CosOfPointingAngleXi > 0.9999) continue;\r
- //if(dcaPosToPrimaryVtxXi < 0.1) continue;\r
- //if(dcaNegToPrimaryVtxXi < 0.1) continue;\r
-\r
- if(dcaV0DaughtersXi > fCascadeCuts[4]) continue;\r
- if(V0CosOfPointingAngleXi > fCascadeCuts[5]) continue;\r
- if(dcaPosToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
- if(dcaNegToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
-\r
- //if(V0RadiusXi < 1.0 || V0RadiusXi > 100) continue;\r
- \r
- //other cuts?\r
- // change mass hypothesis to cover all the possibilities\r
- if(bTrkXi->Charge()<0){\r
- V0quality = 0.;\r
- xi->ChangeMassHypothesis(V0quality, 3312); //Xi- hyp.\r
- invMassXiMinus = xi->GetEffMassXi();\r
-\r
- V0quality = 0.;\r
- xi->ChangeMassHypothesis(V0quality, 3334); //Omega- hyp.\r
- invMassOmegaMinus = xi->GetEffMassXi();\r
-\r
- V0quality = 0.;\r
- xi->ChangeMassHypothesis(V0quality, 3312); //back to default hyp.\r
- }\r
-\r
- if(bTrkXi->Charge() > 0){\r
- V0quality = 0.;\r
- xi->ChangeMassHypothesis(V0quality, -3312); //anti-Xi- hyp.\r
- invMassXiPlus = xi->GetEffMassXi();\r
-\r
- V0quality = 0.;\r
- xi->ChangeMassHypothesis(V0quality, -3334); //anti-Omega- hyp.\r
- invMassOmegaPlus = xi->GetEffMassXi();\r
-\r
- V0quality = 0.;\r
- xi->ChangeMassHypothesis(V0quality, -3312); //back to default hyp.\r
- }\r
-\r
- //PID on the daughter tracks\r
- /*\r
- //A - Combined PID\r
- //Resonable guess the priors for the cascade track sample\r
- //(e, mu, pi, K, p)\r
- Double_t priorsGuessXi[5] = {0, 0, 2, 0, 1};\r
- Double_t priorsGuessOmega[5] = {0, 0, 1, 1, 1};\r
-\r
- //Combined bachelor-daughter PID\r
- AliPID pidXi;\r
- pidXi.SetPriors(priorsGuessXi);\r
- AliPID pidOmega;\r
- pidOmega.SetPriors(priorsGuessOmega);\r
-\r
- if(pTrkXi->IsOn(AliESDtrack::kESDpid)){// combined PID exists\r
- Double_t r[10] = {0.};\r
- pTrkXi->GetESDpid(r);\r
- pidXi.SetProbabilities(r);\r
- pidOmega.SetProbabilities(r);\r
-\r
- //Check if the V0 postive track is proton (case for Xi-)\r
- Double_t pProton = pidXi.GetProbability(AliPID::kProton);\r
- if(pProton > pidXi.GetProbability(AliPID::kElectron)\r
- && pProton > pidXi.GetProbability(AliPID::kMuon)\r
- && pProton > pidXi.GetProbability(AliPID::kPion)\r
- && pProton > pidXi.GetProbability(AliPID::kKaon))\r
- isPosInXiProton = kTRUE;\r
- \r
- //Check if the V0 postive track is a pi+ (case for Xi+)\r
- Double_t pPion = pidXi.GetProbability(AliPID::kPion);\r
- if(pPion > pidXi.GetProbability(AliPID::kElectron)\r
- && pPion > pidXi.GetProbability(AliPID::kMuon)\r
- && pPion > pidXi.GetProbability(AliPID::kKaon)\r
- && pPion > pidXi.GetProbability(AliPID::kProton))\r
- isPosInXiPion = kTRUE;\r
- // Check if the V0 positive track is a proton (case for Omega-)\r
- pProton = pidOmega.GetProbability(AliPID::kProton);\r
- if(pProton > pidOmega.GetProbability(AliPID::kElectron)\r
- && pProton > pidOmega.GetProbability(AliPID::kMuon)\r
- && pProton > pidOmega.GetProbability(AliPID::kPion)\r
- && pProton > pidOmega.GetProbability(AliPID::kKaon))\r
- isPosInOmegaProton = kTRUE;\r
- \r
- // Check if the V0 positive track is a pi+ (case for Omega+)\r
- pPion = pidOmega.GetProbability(AliPID::kPion);\r
- if(pPion > pidOmega.GetProbability(AliPID::kElectron)\r
- && pPion > pidOmega.GetProbability(AliPID::kMuon)\r
- && pPion > pidOmega.GetProbability(AliPID::kKaon)\r
- && pPion > pidOmega.GetProbability(AliPID::kProton))\r
- isPosInOmegaPion = kTRUE;\r
- }\r
-\r
- //Combined V0-negative-daughter PID\r
- pidXi.SetPriors(priorsGuessXi);\r
- pidOmega.SetPriors(priorsGuessOmega);\r
- if(nTrkXi->IsOn(AliESDtrack::kESDpid)){\r
- Double_t r[10] = {0.};\r
- nTrkXi->GetESDpid(r);\r
- pidXi.SetProbabilities(r);\r
- pidOmega.SetProbabilities(r);\r
- \r
- // Check if the V0 negative track is a pi- (case for Xi-)\r
- Double_t pPion = pidXi.GetProbability(AliPID::kPion);\r
- if(pPion > pidXi.GetProbability(AliPID::kElectron)\r
- && pPion > pidXi.GetProbability(AliPID::kMuon)\r
- && pPion > pidXi.GetProbability(AliPID::kKaon)\r
- && pPion > pidXi.GetProbability(AliPID::kProton))\r
- isNegInXiPion = kTRUE;\r
-\r
- // Check if the V0 negative track is an anti-proton (case for Xi+)\r
- Double_t pProton = pidXi.GetProbability(AliPID::kProton);\r
- if(pProton > pidXi.GetProbability(AliPID::kElectron)\r
- && pProton > pidXi.GetProbability(AliPID::kMuon)\r
- && pProton > pidXi.GetProbability(AliPID::kPion) \r
- && pProton > pidXi.GetProbability(AliPID::kKaon))\r
- isNegInXiProton = kTRUE;\r
- \r
- // Check if the V0 negative track is a pi- (case for Omega-)\r
- pPion = pidOmega.GetProbability(AliPID::kPion);\r
- if(pPion > pidOmega.GetProbability(AliPID::kElectron)\r
- && pPion > pidOmega.GetProbability(AliPID::kMuon)\r
- && pPion > pidOmega.GetProbability(AliPID::kKaon)\r
- && pPion > pidOmega.GetProbability(AliPID::kProton))\r
- isNegInOmegaPion = kTRUE;\r
- \r
- // Check if the V0 negative track is an anti-proton (case for Omega+) \r
- pProton = pidOmega.GetProbability(AliPID::kProton);\r
- if(pProton > pidOmega.GetProbability(AliPID::kElectron)\r
- && pProton > pidOmega.GetProbability(AliPID::kMuon)\r
- && pProton > pidOmega.GetProbability(AliPID::kPion)\r
- && pProton > pidOmega.GetProbability(AliPID::kKaon))\r
- isNegInOmegaProton = kTRUE;\r
- \r
- }\r
-\r
- // Combined bachelor PID\r
- pidXi.SetPriors(priorsGuessXi);\r
- pidOmega.SetPriors(priorsGuessOmega);\r
- if(bTrkXi->IsOn(AliESDtrack::kESDpid)){//Combined PID exists\r
- Double_t r[10] = {0.};\r
- bTrkXi->GetESDpid(r);\r
- pidXi.SetProbabilities(r);\r
- pidOmega.SetProbabilities(r);\r
- \r
- //Check if the bachelor track is a pion\r
- Double_t pPion = pidXi.GetProbability(AliPID::kPion);\r
- if(pPion > pidXi.GetProbability(AliPID::kElectron)\r
- && pPion > pidXi.GetProbability(AliPID::kMuon)\r
- && pPion > pidXi.GetProbability(AliPID::kKaon)\r
- && pPion > pidXi.GetProbability(AliPID::kProton))\r
- isBachelorPion = kTRUE;\r
- \r
- // Check if the bachelor track is a kaon\r
- Double_t pKaon = pidOmega.GetProbability(AliPID::kKaon);\r
- if(pKaon > pidOmega.GetProbability(AliPID::kElectron)\r
- && pKaon > pidOmega.GetProbability(AliPID::kMuon)\r
- && pKaon > pidOmega.GetProbability(AliPID::kPion)\r
- && pKaon > pidOmega.GetProbability(AliPID::kProton))\r
- isBachelorKaon = kTRUE;\r
- }\r
- */\r
-\r
-\r
- //B - TPC PID: 3-sigma bands on Bethe-Bloch curve\r
- //Bachelor\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kKaon))<3.)\r
- isBachelorKaonForTPC = kTRUE;\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kPion))<3.)\r
- isBachelorPionForTPC = kTRUE;\r
-\r
- //Negative V0 daughter\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kPion))<3.)\r
- isNegPionForTPC = kTRUE;\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kProton))<3.)\r
- isNegProtonForTPC = kTRUE;\r
- \r
- //Positive V0 daughter\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kPion))<3.)\r
- isPosPionForTPC = kTRUE;\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kProton))<3.)\r
- isPosProtonForTPC = kTRUE;\r
- \r
- \r
- //Extra QA information\r
- xi->GetPxPyPz(XiPx, XiPy, XiPz);\r
- XiPt = TMath::Sqrt(XiPx*XiPx + XiPy*XiPy);\r
- XiPtot= TMath::Sqrt(XiPx*XiPx + XiPy*XiPy + XiPz*XiPz);\r
- \r
- XiPAfter[0] = XiPx;\r
- XiPAfter[1] = XiPy;\r
- XiPAfter[2] = XiPz;\r
-\r
- xi->GetBPxPyPz(bachPx, bachPy, bachPz);\r
- bachPt = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy);\r
- bachPtot = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy + bachPz*bachPz);\r
-\r
- //chargeXi = xi->Charge();\r
- \r
- V0toXiCosOfPointingAngle \r
- = xi->GetV0CosineOfPointingAngle(posXi[0], posXi[1], posXi[2]);\r
- rapXi = xi->RapXi();\r
- rapOmega = xi->RapOmega();\r
- phi = xi->Phi();\r
- alphaXi = xi->AlphaXi();\r
- ptArmXi = xi->PtArmXi();\r
-\r
- distToVtxZBefore = posXi[2]-bestPrimaryVtxPos[2];\r
- distToVtxXYBefore \r
- = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
- *(posXi[0] - bestPrimaryVtxPos[0])\r
- +(posXi[1] - bestPrimaryVtxPos[1])\r
- *(posXi[1] - bestPrimaryVtxPos[1]));\r
- \r
- //propagation to the best primary vertex to determine the momentum\r
- Propagate(bestPrimaryVtxPos, posXi, XiPAfter, b, xi->Charge());\r
- distToVtxZAfter = posXi[2] - bestPrimaryVtxPos[2];\r
- distToVtxXYAfter = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
- *(posXi[0] - bestPrimaryVtxPos[0])\r
- +(posXi[1] - bestPrimaryVtxPos[1])\r
- *(posXi[1] - bestPrimaryVtxPos[1]));\r
- phiAfter = TMath::Pi() + TMath::ATan2(-XiPAfter[1],-XiPAfter[0]);\r
- \r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZAfter"))->Fill(distToVtxZAfter);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYAfter"))->Fill(distToVtxXYAfter);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZBeforeVsAfter"))->Fill(distToVtxZBefore, distToVtxZAfter);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYBeforeVsAfter"))->Fill(distToVtxXYBefore, distToVtxXYAfter);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PxBeforeVsAfter"))->Fill(XiPx, XiPAfter[0]);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PyBeforeVsAfter"))->Fill(XiPy, XiPAfter[1]);\r
- if(xi->Charge()>0)\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiPosBeforeVsAfter"))->Fill(phi, phiAfter);\r
- else if(xi->Charge()<0)\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiNegBeforeVsAfter"))->Fill(phi, phiAfter);\r
-\r
- \r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Chi2Xi"))->Fill(chi2Xi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaXiDaughters"))->Fill(dcaXiDaughters);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaBachToPrimVertex"))->Fill(dcaBachToPrimaryVtxXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiCosOfPointingAngle"))->Fill(XiCosOfPointingAngle);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiRadius"))->Fill(XiRadius);\r
- \r
- //V0\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassLambdaAsCascDghter"))->Fill(invMassLambdaAsCascDghter);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0Chi2Xi"))->Fill(V0Chi2Xi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0DaughtersXi"))->Fill(dcaV0DaughtersXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0CosOfPointingAngleXi"))->Fill(V0CosOfPointingAngleXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0RadiusXi"))->Fill(V0RadiusXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0ToPrimVertexXi"))->Fill(dcaV0ToPrimaryVtxXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaPosToPrimVertexXi"))->Fill(dcaPosToPrimaryVtxXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaNegToPrimVertexXi"))->Fill(dcaNegToPrimaryVtxXi);\r
- \r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0toXiCosOfPointingAngle"))->Fill(V0toXiCosOfPointingAngle);\r
- \r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Armenteros"))->Fill(alphaXi, ptArmXi);\r
- \r
- //PID cuts with TPC cuts\r
- if(xi->Charge() < 0){\r
- if(isPosProtonForTPC\r
- && isNegPionForTPC){\r
- \r
- switch(fSpecie) {\r
- case 0:\r
- if(isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiMinus);\r
- \r
- //candidate inserting\r
- MakeTrack(invMassXiMinus, XiPt, /*xi->Phi()*/\r
- phiAfter, xi->Eta(), pTrkXi->GetID(),\r
- nTrkXi->GetID(), bTrkXi->GetID());\r
- }\r
- break;\r
- \r
- case 1:\r
- if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaMinus);\r
- MakeTrack(invMassOmegaMinus, XiPt, /*xi->Phi()*/\r
- phiAfter, xi->Eta(),\r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }\r
- break;\r
- }\r
- }\r
- }\r
-\r
- if(xi->Charge() > 0){\r
- if(isNegProtonForTPC\r
- && isPosPionForTPC){\r
- \r
- switch (fSpecie){\r
- case 0:\r
- if(isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiPlus);\r
- \r
- //candidate inserting \r
- MakeTrack(invMassXiPlus, XiPt, /*xi->Phi()*/\r
- phiAfter, xi->Eta(),\r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }\r
- break;\r
-\r
- case 1:\r
- if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaPlus);\r
- MakeTrack(invMassOmegaPlus, XiPt, /*xi->Phi()*/\r
- phiAfter, xi->Eta(),\r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }\r
- break;\r
- }\r
- }\r
- }\r
- \r
- }\r
-\r
- return;\r
-}\r
-\r
-//______________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::ReadFromAODv0(AliAODEvent *fAOD)\r
-{\r
- \r
- fCutsRPTPC->SetEvent(fAOD, MCEvent());\r
- fCutsRPVZE->SetEvent(fAOD, MCEvent());\r
- fCutsPOI->SetEvent(fAOD, MCEvent());\r
- fFlowEventTPC->Fill(fCutsRPTPC, fCutsPOI);\r
- fFlowEventVZE->Fill(fCutsRPVZE, fCutsPOI);\r
-\r
- // Double_t trkPrimaryVtxPos[3] = {-100., -100., -100.};\r
- Double_t bestPrimaryVtxPos[3] = {-100., -100., -100.};\r
-\r
- Double_t b = fAOD->GetMagneticField();\r
-\r
- int nCascades=fAOD->GetNumberOfCascades();\r
- const AliAODVertex *primaryBestAODVtx = fAOD->GetPrimaryVertex();\r
- primaryBestAODVtx->GetXYZ(bestPrimaryVtxPos);\r
- \r
- // calculation part dedicated to Xi vertices\r
- for(Int_t iXi = 0; iXi < nCascades; iXi++){\r
- Double_t effMassXi = 0.;\r
- Double_t chi2Xi = -1.;\r
- Double_t dcaXiDaughters = -1.;\r
- Double_t XiCosOfPointingAngle = -1.;\r
- Double_t posXi[3] = {-1000., -1000., -1000.};\r
- Double_t XiRadius = -1000.;\r
- \r
- Double_t invMassLambdaAsCascDghter = 0.;\r
- Double_t V0Chi2Xi = -1.;\r
- Double_t dcaV0DaughtersXi = -1.;\r
- \r
- Double_t dcaBachToPrimaryVtxXi = -1.;\r
- Double_t dcaV0ToPrimaryVtxXi = -1.;\r
- Double_t dcaPosToPrimaryVtxXi = -1.;\r
- Double_t dcaNegToPrimaryVtxXi = -1.;\r
- Double_t V0CosOfPointingAngleXi = -1.;\r
- Double_t posV0Xi[3] = {-1000., -1000., -1000.};\r
- Double_t V0RadiusXi = -1000.;\r
- // Double_t V0quality = 0.;\r
-\r
- Double_t invMassXiMinus = 0.;\r
- Double_t invMassXiPlus = 0.;\r
- Double_t invMassOmegaMinus = 0.;\r
- Double_t invMassOmegaPlus = 0.;\r
- \r
- /*\r
- Bool_t isPosInXiProton = kFALSE;\r
- Bool_t isPosInXiPion = kFALSE;\r
- Bool_t isPosInOmegaProton = kFALSE;\r
- Bool_t isPosInOmegaPion = kFALSE;\r
- \r
- Bool_t isNegInXiProton = kFALSE;\r
- Bool_t isNegInXiPion = kFALSE;\r
- Bool_t isNegInOmegaProton = kFALSE;\r
- Bool_t isNegInOmegaPion = kFALSE;\r
-\r
- Bool_t isBachelorKaon = kFALSE;\r
- Bool_t isBachelorPion = kFALSE;\r
- */\r
-\r
-\r
- Bool_t isBachelorKaonForTPC = kFALSE;\r
- Bool_t isBachelorPionForTPC = kFALSE;\r
- Bool_t isNegPionForTPC = kFALSE;\r
- Bool_t isPosPionForTPC = kFALSE;\r
- Bool_t isNegProtonForTPC = kFALSE;\r
- Bool_t isPosProtonForTPC = kFALSE;\r
- \r
- Double_t XiPx = 0., XiPy = 0., XiPz = 0.;\r
- Double_t XiPt = 0.;\r
- Double_t XiPtot = 0.;\r
- \r
- Double_t bachPx = 0., bachPy = 0., bachPz = 0.;\r
- Double_t bachPt = 0.;\r
- Double_t bachPtot = 0.;\r
- \r
- //Short_t chargeXi = -2;\r
- Double_t V0toXiCosOfPointingAngle = 0.;\r
- \r
- Double_t rapXi = -20.;\r
- Double_t rapOmega = -20.;\r
- Double_t phi = 6.3;\r
- Double_t alphaXi = -200.;\r
- Double_t ptArmXi = -200.;\r
-\r
- Double_t distToVtxZBefore = -999.;\r
- Double_t distToVtxZAfter = -999.;\r
- Double_t distToVtxXYBefore = -999.;\r
- Double_t distToVtxXYAfter = -999.;\r
- Double_t XiPAfter[3] = {-999., -999., -999.};\r
- Double_t phiAfter = -999.;\r
-\r
- const AliAODcascade *xi = fAOD->GetCascade(iXi);\r
- if (!xi) continue;\r
-\r
- effMassXi = xi->MassXi(); //default working hypothesis: Xi- decay\r
- chi2Xi = xi->Chi2Xi();\r
- dcaXiDaughters = xi->DcaXiDaughters();\r
- XiCosOfPointingAngle = xi->CosPointingAngleXi(bestPrimaryVtxPos[0],\r
- bestPrimaryVtxPos[1],\r
- bestPrimaryVtxPos[2]);\r
- posXi[0] = xi->DecayVertexXiX();\r
- posXi[1] = xi->DecayVertexXiY();\r
- posXi[2] = xi->DecayVertexXiZ();\r
- XiRadius = TMath::Sqrt(posXi[0]*posXi[0]\r
- +posXi[1]*posXi[1]\r
- +posXi[2]*posXi[2]);\r
-\r
- AliAODTrack *pTrkXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(0) );\r
- AliAODTrack *nTrkXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(1) );\r
- AliAODTrack *bTrkXi \r
- = dynamic_cast<AliAODTrack*>( xi->GetDecayVertexXi()->GetDaughter(0) );\r
-\r
- if(!pTrkXi || !nTrkXi || !bTrkXi) continue;\r
-\r
- UInt_t idxPosXi = (UInt_t) TMath::Abs( pTrkXi->GetID() );\r
- UInt_t idxNegXi = (UInt_t) TMath::Abs( nTrkXi->GetID() );\r
- UInt_t idxBach = (UInt_t) TMath::Abs( bTrkXi->GetID() );\r
-\r
- if(idxBach == idxNegXi || idxBach == idxPosXi) continue;\r
-\r
- if( !fCutsDau->IsSelected(pTrkXi) \r
- || !fCutsDau->IsSelected(nTrkXi)\r
- || !fCutsDau->IsSelected(bTrkXi) ) continue;\r
-\r
- \r
- if(pTrkXi->IsOn(AliESDtrack::kTPCin)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(pTrkXi->P()*pTrkXi->Charge(), pTrkXi->GetTPCsignal());\r
- }\r
- if( nTrkXi->IsOn(AliESDtrack::kTPCin) ){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(nTrkXi->P()*nTrkXi->Charge(), nTrkXi->GetTPCsignal());\r
- }\r
- if(bTrkXi->IsOn(AliESDtrack::kTPCin)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(bTrkXi->P()*bTrkXi->Charge(), bTrkXi->GetTPCsignal());\r
- }\r
- \r
- if(xi->ChargeXi() < 0)\r
- invMassLambdaAsCascDghter = xi->MassLambda();\r
- else\r
- invMassLambdaAsCascDghter = xi->MassAntiLambda();\r
- \r
- dcaV0DaughtersXi = xi->DcaV0Daughters();\r
- V0Chi2Xi = xi->Chi2V0();\r
- V0CosOfPointingAngleXi \r
- = xi->CosPointingAngle(bestPrimaryVtxPos);\r
- dcaV0ToPrimaryVtxXi = xi->DcaV0ToPrimVertex();\r
- dcaBachToPrimaryVtxXi = xi->DcaBachToPrimVertex();\r
- \r
- //V0\r
- posV0Xi[0] = xi->DecayVertexV0X();\r
- posV0Xi[1] = xi->DecayVertexV0Y();\r
- posV0Xi[2] = xi->DecayVertexV0Z();\r
- V0RadiusXi = TMath::Sqrt(posV0Xi[0]*posV0Xi[0]\r
- +posV0Xi[1]*posV0Xi[1]\r
- +posV0Xi[2]*posV0Xi[2]);\r
- dcaPosToPrimaryVtxXi = xi->DcaPosToPrimVertex();\r
- dcaNegToPrimaryVtxXi = xi->DcaNegToPrimVertex();\r
-\r
- //apply cuts ?\r
- // if(XiRadius < 1. || XiRadius > 100.) continue;\r
- //if(dcaXiDaughters > 0.1) continue;\r
- //if(XiCosOfPointingAngle < 0.999) continue;\r
- //if(dcaV0ToPrimaryVtxXi < 0.05) continue;\r
- //if(dcaBachToPrimaryVtxXi < 0.03) continue;\r
-\r
- if(dcaXiDaughters > fCascadeCuts[0]) continue;\r
- if(XiCosOfPointingAngle < fCascadeCuts[1]) continue;\r
- if(dcaV0ToPrimaryVtxXi < fCascadeCuts[2]) continue;\r
- if(dcaBachToPrimaryVtxXi < fCascadeCuts[3]) continue;\r
- \r
- //V0 mass cut?\r
- //if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > 0.006) continue;\r
- if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > fCascadeCuts[7]) \r
- continue;\r
-\r
- //if(dcaV0DaughtersXi > 1.) continue;\r
- //if(V0CosOfPointingAngleXi > 0.9999) continue;\r
- //if(dcaPosToPrimaryVtxXi < 0.1) continue;\r
- //if(dcaNegToPrimaryVtxXi < 0.1) continue;\r
- if(dcaV0DaughtersXi > fCascadeCuts[4]) continue;\r
- if(V0CosOfPointingAngleXi > fCascadeCuts[5]) continue;\r
- if(dcaPosToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
- if(dcaNegToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
- \r
- // if(V0RadiusXi < 1.0 || V0RadiusXi > 100) continue;\r
- \r
- //other cuts?\r
-\r
-\r
- //???\r
- if(xi->ChargeXi()<0){\r
- invMassXiMinus = xi->MassXi();\r
- invMassOmegaMinus = xi->MassOmega();\r
- }else{\r
- invMassXiPlus = xi->MassXi();\r
- invMassOmegaPlus = xi->MassOmega();\r
- }\r
-\r
- /*\r
- if(pTrkXi->GetMostProbablePID() == AliAODTrack::kProton) {\r
- isPosInXiProton = kTRUE;\r
- isPosInOmegaProton = kTRUE;\r
- }\r
- if(pTrkXi->GetMostProbablePID() == AliAODTrack::kPion){\r
- isPosInXiPion = kTRUE;\r
- isPosInOmegaPion = kTRUE;\r
- }\r
- \r
- if(nTrkXi->GetMostProbablePID() == AliAODTrack::kPion){\r
- isNegInXiPion = kTRUE;\r
- isNegInOmegaPion = kTRUE;\r
- }\r
- if(nTrkXi->GetMostProbablePID() == AliAODTrack::kProton){\r
- isNegInXiProton = kTRUE;\r
- isNegInOmegaProton = kTRUE;\r
- }\r
-\r
- if(bTrkXi->GetMostProbablePID() == AliAODTrack::kPion)\r
- isBachelorPion = kTRUE;\r
- if(bTrkXi->GetMostProbablePID() == AliAODTrack::kKaon)\r
- isBachelorKaon = kTRUE;\r
- */\r
-\r
- //PID with TPC only: \r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kKaon))<3.)\r
- isBachelorKaonForTPC = kTRUE;\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kPion))<3.)\r
- isBachelorPionForTPC = kTRUE;\r
-\r
- //Negative V0 daughter\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kPion))<3.)\r
- isNegPionForTPC = kTRUE;\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kProton))<3.)\r
- isNegProtonForTPC = kTRUE;\r
- \r
- //Positive V0 daughter\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kPion))<3.)\r
- isPosPionForTPC = kTRUE;\r
- if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kProton))<3.)\r
- isPosProtonForTPC = kTRUE;\r
-\r
- //Extra QA information\r
- XiPx = xi->MomXiX();\r
- XiPy = xi->MomXiY();\r
- XiPz = xi->MomXiZ();\r
- XiPt = TMath::Sqrt(XiPx*XiPx + XiPy*XiPy);\r
- XiPtot= TMath::Sqrt(XiPx*XiPx + XiPy*XiPy + XiPz*XiPz);\r
- \r
- bachPx = xi->MomBachX();\r
- bachPy = xi->MomBachY();\r
- bachPz = xi->MomBachZ();\r
- \r
- bachPt = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy);\r
- bachPtot = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy + bachPz*bachPz);\r
- \r
- V0toXiCosOfPointingAngle = xi->CosPointingAngle( xi->GetDecayVertexXi() );\r
- \r
- rapXi = xi->RapXi();\r
- rapOmega = xi->RapOmega();\r
- phi = xi->Phi();\r
- alphaXi = xi->AlphaXi();\r
- ptArmXi = xi->PtArmXi();\r
-\r
- distToVtxZBefore = posXi[2]-bestPrimaryVtxPos[2];\r
- distToVtxXYBefore\r
- = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
- *(posXi[0] - bestPrimaryVtxPos[0])\r
- +(posXi[1] - bestPrimaryVtxPos[1])\r
- *(posXi[1] - bestPrimaryVtxPos[1]));\r
-\r
-\r
- XiPAfter[0] = XiPx;\r
- XiPAfter[1] = XiPy;\r
- XiPAfter[2] = XiPz;\r
- //propagation to the best primary vertex to determine the momentum\r
- Propagate(bestPrimaryVtxPos, posXi, XiPAfter, b, xi->ChargeXi());\r
- distToVtxZAfter = posXi[2] - bestPrimaryVtxPos[2];\r
- distToVtxXYAfter = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
- *(posXi[0] - bestPrimaryVtxPos[0])\r
- +(posXi[1] - bestPrimaryVtxPos[1])\r
- *(posXi[1] - bestPrimaryVtxPos[1]));\r
- phiAfter = TMath::Pi() + TMath::ATan2(-XiPAfter[1],-XiPAfter[0]);\r
-\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZAfter"))->Fill(distToVtxZAfter);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYAfter"))->Fill(distToVtxXYAfter);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZBeforeVsAfter"))->Fill(distToVtxZBefore, distToVtxZAfter);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYBeforeVsAfter"))->Fill(distToVtxXYBefore, distToVtxXYAfter);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PxBeforeVsAfter"))->Fill(XiPx, XiPAfter[0]);\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PyBeforeVsAfter"))->Fill(XiPy, XiPAfter[1]);\r
- if(xi->ChargeXi()>0)\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiPosBeforeVsAfter"))->Fill(phi, phiAfter);\r
- else if(xi->ChargeXi()<0)\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiNegBeforeVsAfter"))->Fill(phi, phiAfter);\r
- \r
- //for default hypothesis\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Chi2Xi"))->Fill(chi2Xi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaXiDaughters"))->Fill(dcaXiDaughters);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaBachToPrimVertex"))->Fill(dcaBachToPrimaryVtxXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiCosOfPointingAngle"))->Fill(XiCosOfPointingAngle);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiRadius"))->Fill(XiRadius);\r
- \r
- //V0\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassLambdaAsCascDghter"))->Fill(invMassLambdaAsCascDghter);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0Chi2Xi"))->Fill(V0Chi2Xi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0DaughtersXi"))->Fill(dcaV0DaughtersXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0CosOfPointingAngleXi"))->Fill(V0CosOfPointingAngleXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0RadiusXi"))->Fill(V0RadiusXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0ToPrimVertexXi"))->Fill(dcaV0ToPrimaryVtxXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaPosToPrimVertexXi"))->Fill(dcaPosToPrimaryVtxXi);\r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaNegToPrimVertexXi"))->Fill(dcaNegToPrimaryVtxXi);\r
- \r
- ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0toXiCosOfPointingAngle"))->Fill(V0toXiCosOfPointingAngle);\r
- \r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Armenteros"))->Fill(alphaXi, ptArmXi);\r
- \r
- //with PID cuts\r
- if(xi->ChargeXi()<0){\r
- if(isPosProtonForTPC && isNegPionForTPC){\r
- switch (fSpecie){\r
- case 0:\r
- if( isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiMinus);\r
- MakeTrack(invMassXiMinus, XiPt, /*xi->Phi(),*/\r
- phiAfter, xi->Eta(), \r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }// endif\r
- \r
- break;\r
-\r
- case 1:\r
- if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8\r
- && (invMassXiMinus > 1.32486 || invMassXiMinus < 1.30486)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaMinus); \r
- \r
- MakeTrack(invMassOmegaMinus, XiPt, /* xi->Phi(),*/\r
- phiAfter, xi->Eta(),\r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }//endif\r
- break;\r
- }\r
- }\r
- }//end if ChargeXi()<0\r
- \r
- if(xi->ChargeXi() > 0){\r
- if(isNegProtonForTPC && isPosPionForTPC){ \r
- switch(fSpecie){\r
- case 0:\r
- if (isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiPlus);\r
- \r
- //candidate inserting \r
- MakeTrack(invMassXiPlus, XiPt, /* xi->Phi(),*/\r
- phiAfter, xi->Eta(),\r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }//endif particle id\r
- break;\r
- \r
- case 1:\r
- if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8\r
- && (invMassXiPlus > 1.32486 || invMassXiPlus < 1.30486)){\r
- ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaPlus);\r
- MakeTrack(invMassOmegaPlus, XiPt, /* xi->Phi(),*/\r
- phiAfter, xi->Eta(),\r
- pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
- }//endif particle id\r
- }\r
- }\r
- }//endif ChargeXi()>0 \r
- }//for Xi candidate loop\r
-\r
- return;\r
-\r
-}\r
-\r
-\r
-void AliAnalysisTaskFlowCascade::MakeTrack( double mass, \r
- double pt, \r
- double phi, \r
- double eta, \r
- int iid, \r
- int jid,\r
- int kid) {\r
- // create track for flow tasks \r
- if(fCandidates->GetLast()+1>=fCandidates->GetSize()) {\r
- fCandidates->Expand( 2*fCandidates->GetSize() );\r
- }\r
- Bool_t overwrite = kTRUE;\r
-\r
- AliFlowCandidateTrack *sTrack \r
- = (static_cast<AliFlowCandidateTrack*> (fCandidates->At( fCandidates->GetLast()+1 )));\r
- if( !sTrack ) { // creates new\r
- sTrack = new AliFlowCandidateTrack();\r
- overwrite = kFALSE;\r
- } else { // overwrites\r
- sTrack->ClearMe();\r
- }\r
-\r
-\r
- sTrack->SetMass(mass);\r
- sTrack->SetPt(pt);\r
- sTrack->SetPhi(phi);\r
- sTrack->SetEta(eta);\r
- sTrack->AddDaughter(iid);\r
- sTrack->AddDaughter(jid);\r
- sTrack->AddDaughter(kid);\r
- sTrack->SetForPOISelection(kTRUE);\r
- sTrack->SetForRPSelection(kFALSE);\r
- \r
- if(overwrite) {\r
- fCandidates->SetLast( fCandidates->GetLast()+1 );\r
- } else {\r
- fCandidates->AddLast(sTrack);\r
- }\r
- \r
- return;\r
-}\r
-//_____________________________________________________________________________\r
-void AliAnalysisTaskFlowCascade::Terminate(Option_t *)\r
-{\r
-\r
-}\r
-\r
-void AliAnalysisTaskFlowCascade::Propagate(Double_t vv[3], \r
- Double_t x[3], \r
- Double_t p[3], \r
- Double_t bz, \r
- Short_t sign){\r
- //Propagation to the primary vertex to determine the px and py\r
- //x, p are the position and momentum as input and output\r
- //bz is the magnetic field along z direction\r
- //sign is the charge of particle for propagation\r
-\r
- Double_t pp = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);\r
- Double_t len = (vv[2]-x[2])*pp/p[2];\r
- Double_t a = -kB2C*bz*sign; \r
-\r
- Double_t rho = a/pp;\r
- x[0] += p[0]*TMath::Sin(rho*len)/a - p[1]*(1-TMath::Cos(rho*len))/a;\r
- x[1] += p[1]*TMath::Sin(rho*len)/a + p[0]*(1-TMath::Cos(rho*len))/a;\r
- x[2] += p[2]*len/pp;\r
-\r
- Double_t p0=p[0];\r
- p[0] = p0 *TMath::Cos(rho*len) - p[1]*TMath::Sin(rho*len);\r
- p[1] = p[1]*TMath::Cos(rho*len) + p0 *TMath::Sin(rho*len);\r
-}\r
-\r
-\r
-//===================================================================== \r
-void AliAnalysisTaskFlowCascade::SetCommonConstants(Int_t massBins, \r
- Double_t minMass, \r
- Double_t maxMass)\r
-{\r
- // setter for mass bins \r
- \r
- fMassBins = massBins;\r
- fMinMass = minMass;\r
- fMaxMass = maxMass;\r
-}\r
-\r
-\r
-//==================================================================== \r
-void AliAnalysisTaskFlowCascade::SetCuts2010(int set) {\r
-\r
- // fCascadeCuts[0]: DcaXiDaughter; fCascadeCuts[1]: XiCosOfPointingAngle\r
- // fCascadeCuts[2]: DcaV0ToPrimaryVtxXi; fCascadeCuts[3]: DcaBachToPrimaryVtxXi\r
- // fCascadeCuts[4]: DcaV0DaughersXi; fCascadeCuts[5]: V0CosOfPointingAngleXi\r
- // fCascadeCuts[6]: DcaV0DaughterToPrimaryVtxXi; fCascadeCuts[7]: V0MassWidth\r
- \r
- switch(set){\r
- \r
- case 0: //tighter\r
- fCascadeCuts[0] = 0.2; fCascadeCuts[1] = 0.999;\r
- fCascadeCuts[2] = 0.03; fCascadeCuts[3] = 0.05;\r
- fCascadeCuts[4] = .5; fCascadeCuts[5] = 0.9998;\r
- fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.006;\r
- break;\r
- \r
- case 1: //middle\r
- fCascadeCuts[0] = 0.3; fCascadeCuts[1] = 0.99;\r
- fCascadeCuts[2] = 0.01; fCascadeCuts[3] = 0.03;\r
- fCascadeCuts[4] = .6; fCascadeCuts[5] = 0.9999;\r
- fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.008;\r
- break;\r
-\r
- case 2: //looser\r
- fCascadeCuts[0] = 0.3; fCascadeCuts[1] = 0.99;\r
- fCascadeCuts[2] = 0.01; fCascadeCuts[3] = 0.03;\r
- fCascadeCuts[4] = 1.; fCascadeCuts[5] = 1.;\r
- fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.01;\r
- break;\r
- }\r
- \r
-}\r
+/**************************************************************************
+* Copyright(c) 1998-2008,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. *
+**************************************************************************/
+
+/////////////////////////////////////////////////////
+// AliAnalysisTaskFlowCascade:
+// Analysis task to select Xi and Omega candidates for flow analysis.
+//
+// Author: Zhong-Bao.Yin@cern.ch
+//////////////////////////////////////////////////////
+
+#include "TChain.h"
+#include "TList.h"
+#include "TH1D.h"
+#include "TH2D.h"
+#include "TH3D.h"
+#include "TProfile.h"
+#include "TVector3.h"
+
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisManager.h"
+
+#include "AliVParticle.h"
+#include "AliESDEvent.h"
+#include "AliESDInputHandler.h"
+#include "AliESDv0.h"
+#include "AliESDcascade.h"
+#include "AliESDtrack.h"
+#include "AliESDtrackCuts.h"
+#include "AliCentrality.h"
+#include "AliVVertex.h"
+#include "AliESDVZERO.h"
+#include "AliESDUtils.h"
+
+#include "AliTPCPIDResponse.h"
+#include "AliTOFPIDResponse.h"
+#include "AliPIDResponse.h"
+
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
+#include "AliAODTrack.h"
+#include "AliAODVZERO.h"
+#include "AliAODcascade.h"
+
+#include "TMath.h"
+#include "TObjArray.h"
+#include "AliFlowCandidateTrack.h"
+
+#include "AliFlowTrackCuts.h"
+#include "AliFlowEventCuts.h"
+#include "AliFlowEvent.h"
+#include "AliFlowCommonConstants.h"
+
+#include "AliAnalysisTaskFlowCascade.h"
+
+ClassImp(AliAnalysisTaskFlowCascade)
+
+//_____________________________________________________________________________
+ AliAnalysisTaskFlowCascade::AliAnalysisTaskFlowCascade() :
+ AliAnalysisTaskSE(),
+ // fMinCent(0), fMaxCent(0),
+ fSpecie(0),
+ fMassBins(0),
+ fMinMass(0.0),
+ fMaxMass(0.0),
+ fCutsEvent(NULL),
+ fCutsRPTPC(NULL),
+ fCutsRPVZE(NULL),
+ fCutsPOI(NULL),
+ fCutsDau(NULL),
+ fPIDResponse(NULL),
+ fFlowEventTPC(NULL),
+ fFlowEventVZE(NULL),
+ fCandidates(NULL),
+ fQAList(NULL)
+{
+ //ctor
+ for (Int_t i=0; i!=8; ++i)
+ fCascadeCuts[i] = 0;
+
+}
+
+//_____________________________________________________________________________
+AliAnalysisTaskFlowCascade
+::AliAnalysisTaskFlowCascade(const char *name,
+ AliFlowEventCuts *cutsEvent,
+ AliFlowTrackCuts *cutsRPTPC,
+ AliFlowTrackCuts *cutsRPVZE,
+ /* AliESDtrackCuts */ AliFlowTrackCuts *cutsDau ) :
+ AliAnalysisTaskSE(name),
+ //fMinCent(minCent), fMaxCent(maxCent),
+ fSpecie(0),
+ fMassBins(0),
+ fMinMass(0.0),
+ fMaxMass(0.0),
+ fCutsEvent(cutsEvent),
+ fCutsRPTPC(cutsRPTPC),
+ fCutsRPVZE(cutsRPVZE),
+ fCutsPOI(NULL),
+ fCutsDau(cutsDau),
+ fPIDResponse(NULL),
+ fFlowEventTPC(NULL),
+ fFlowEventVZE(NULL),
+ fCandidates(NULL),
+ fQAList(NULL)
+{
+ //ctor
+ for (Int_t i=0; i!=8; ++i)
+ fCascadeCuts[i] = 0;
+
+ DefineInput( 0,TChain::Class());
+ DefineOutput(1,AliFlowEventSimple::Class()); // TPC object
+ DefineOutput(2,AliFlowEventSimple::Class()); // VZE object
+ DefineOutput(3,TList::Class());
+}
+
+//_____________________________________________________________________________
+AliAnalysisTaskFlowCascade::~AliAnalysisTaskFlowCascade()
+{
+ if(fQAList) delete fQAList;
+ if (fFlowEventTPC) delete fFlowEventTPC;
+ if (fFlowEventVZE) delete fFlowEventVZE;
+ if (fCandidates) delete fCandidates;
+ if (fCutsDau) delete fCutsDau;
+ if (fCutsPOI) delete fCutsPOI;
+ if (fCutsRPTPC) delete fCutsRPTPC;
+ if (fCutsRPVZE) delete fCutsRPVZE;
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::UserCreateOutputObjects()
+{
+
+ AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+ AliInputEventHandler* inputHandler
+ = (AliInputEventHandler*) (man->GetInputEventHandler());
+ fPIDResponse = inputHandler->GetPIDResponse();
+
+ fQAList = new TList();
+ fQAList->SetOwner();
+ AddQAEvents();
+ AddQACandidates();
+ // PostData(3,fQAList);
+
+ AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();
+ cc->SetNbinsMult(1);
+ cc->SetMultMin(0);
+ cc->SetMultMax(1);
+
+ cc->SetNbinsPt(20);
+ cc->SetPtMin(0.0);
+ cc->SetPtMax(10.0);
+
+ cc->SetNbinsPhi(1);
+ cc->SetPhiMin(0.0);
+ cc->SetPhiMax(TMath::TwoPi());
+
+ cc->SetNbinsEta(1);
+ cc->SetEtaMin(-2.0);
+ cc->SetEtaMax(+2.0);
+
+ cc->SetNbinsQ(3);
+ cc->SetQMin(0.0);
+ cc->SetQMax(3.0);
+
+ cc->SetNbinsMass(fMassBins);
+ cc->SetMassMin(fMinMass);
+ cc->SetMassMax(fMaxMass);
+
+ fCutsPOI = new AliFlowTrackCuts("null_cuts");
+ fCutsPOI->SetParamType(fCutsRPTPC->GetParamType());
+ fCutsPOI->SetPtRange(+1,-1); // select nothing QUICK
+ fCutsPOI->SetEtaRange(+1,-1); // select nothing VZERO
+
+ fFlowEventTPC = new AliFlowEvent(3000);
+ fFlowEventVZE = new AliFlowEvent(1000);
+ fCandidates = new TObjArray(100);
+ fCandidates->SetOwner();
+
+ PostData(1,fFlowEventTPC);
+ PostData(2,fFlowEventVZE);
+ PostData(3,fQAList);
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::AddQAEvents()
+{
+ TList *tQAEvents = new TList();
+ tQAEvents->SetName("Events");
+ tQAEvents->SetOwner();
+ TH1I* tEvent = new TH1I("Event","Number of Events", 3,0,3);
+ tQAEvents->Add(tEvent);
+
+ TH1D *tTPCRFP = new TH1D("RFPTPC",
+ "TPC Reference Flow Particles;multiplicity",
+ 100, 0, 3000);
+ tQAEvents->Add(tTPCRFP);
+ TH1D *tVZERFP = new TH1D("RFPVZE",
+ "VZERO Reference Flow Particles;multiplicity",
+ 100, 0, 30000);
+ tQAEvents->Add(tVZERFP);
+
+ TProfile *tCuts = new TProfile("Cuts","Analysis Cuts",10,0,10);
+ tCuts->Fill(0.5,fCascadeCuts[0],1);
+ tCuts->GetXaxis()->SetBinLabel(1,"dcaXiDau");
+ tCuts->Fill(1.5,fCascadeCuts[1],1);
+ tCuts->GetXaxis()->SetBinLabel(2,"XiCPA");
+ tCuts->Fill(2.5,fCascadeCuts[2],1);
+ tCuts->GetXaxis()->SetBinLabel(3,"dcaV0Vtx");
+ tCuts->Fill(3.5,fCascadeCuts[3],1);
+ tCuts->GetXaxis()->SetBinLabel(4,"dcaBachVtx");
+ tCuts->Fill(4.5,fCascadeCuts[4],1);
+ tCuts->GetXaxis()->SetBinLabel(5,"dcaV0Dau");
+ tCuts->Fill(5.5,fCascadeCuts[5],1);
+ tCuts->GetXaxis()->SetBinLabel(6,"V0CPA");
+ tCuts->Fill(6.5,fCascadeCuts[6],1);
+ tCuts->GetXaxis()->SetBinLabel(7,"dcaV0DauVtx");
+ tCuts->Fill(7.5,fCascadeCuts[7],1);
+ tCuts->GetXaxis()->SetBinLabel(8,"V0Mass");
+ tQAEvents->Add(tCuts);
+
+ fQAList->Add(tQAEvents);
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::AddQACandidates()
+{
+ TList *tQACandidates;
+
+ tQACandidates = new TList();
+ tQACandidates->SetOwner();
+ tQACandidates->SetName("Candidates");
+
+ TH1F* tChi2Xi = new TH1F("Chi2Xi",
+ "Cascade #chi^{2}; #chi^{2}; Number of Cascades",
+ 160, 0, 160);
+ tQACandidates->Add(tChi2Xi);
+
+ TH1F* tDCAXiDaughters
+ = new TH1F( "DcaXiDaughters",
+ "DCA between Xi Daughters; DCA (cm); Number of Cascades",
+ 100, 0., 0.5);
+ tQACandidates->Add(tDCAXiDaughters);
+
+ TH1F * tDCABachToPrimVertex
+ = new TH1F("DcaBachToPrimVertex",
+ "DCA of Bach. to Prim. Vertex; DCA (cm);Number of Cascades",
+ 250, 0., 2.5);
+ tQACandidates->Add(tDCABachToPrimVertex);
+
+ TH1F * tXiCosOfPointingAngle
+ = new TH1F("XiCosOfPointingAngle",
+ "Cos of Xi Pointing Angle; Cos (Xi Point.Angl);Number of Xis",
+ 200, 0.99, 1.0);
+ tQACandidates->Add(tXiCosOfPointingAngle);
+
+ TH1F * tXiRadius = new TH1F("XiRadius",
+ "Casc. decay transv. radius; r (cm); Counts" ,
+ 1050, 0., 105.0 );
+ tQACandidates->Add(tXiRadius);
+
+ TH1F *tMassLambda
+ = new TH1F("MassLambdaAsCascDghter",
+ "#Lambda assoc. to Casc. candidates; Eff. Mass (GeV/c^{2}); Counts",
+ 300,1.00,1.3);
+ tQACandidates->Add(tMassLambda);
+
+ TH1F *tV0Chi2 = new TH1F("V0Chi2Xi",
+ "V0 #chi^{2}, in cascade; #chi^{2};Counts",
+ 160, 0, 40);
+ tQACandidates->Add(tV0Chi2);
+
+ TH1F * tV0CosOfPointingAngle
+ = new TH1F("V0CosOfPointingAngleXi",
+ "Cos of V0 Pointing Angle, in cascade;Cos(V0 Point. Angl); Counts",
+ 200, 0.98, 1.0);
+ tQACandidates->Add(tV0CosOfPointingAngle);
+
+ TH1F *tV0Radius = new TH1F("V0RadiusXi",
+ "V0 decay radius, in cascade; radius (cm); Counts",
+ 1050, 0., 105.0);
+ tQACandidates->Add(tV0Radius);
+
+ TH1F * tDcaV0DaughtersXi
+ = new TH1F("DcaV0DaughtersXi",
+ "DCA between V0 daughters, in cascade;DCA (cm);Number of V0s",
+ 120, 0., 0.6);
+ tQACandidates->Add(tDcaV0DaughtersXi);
+
+ TH1F * tDcaV0ToPrimVertex
+ = new TH1F("DcaV0ToPrimVertexXi",
+ "DCA of V0 to Prim. Vertex, in cascade;DCA (cm);Number of Cascades", 200, 0., 1.);
+ tQACandidates->Add(tDcaV0ToPrimVertex);
+
+ TH1F * tDCAPosToPrimVertex =
+ new TH1F("DcaPosToPrimVertexXi",
+ "DCA of V0 pos daughter to Prim. Vertex;DCA (cm);Counts",
+ 300, 0, 3);
+ tQACandidates->Add(tDCAPosToPrimVertex);
+
+ TH1F * tDCANegToPrimVertex
+ = new TH1F("DcaNegToPrimVertexXi",
+ "DCA of V0 neg daughter to Prim. Vertex;DCA (cm);Counts",
+ 300, 0, 3);
+ tQACandidates->Add(tDCANegToPrimVertex);
+
+ TH1F *tV0toXiCosOfPointingAngle
+ = new TH1F("V0toXiCosOfPointingAngle",
+ "Cos. of V0 Ptng Angl Xi vtx; Cos(V0 Point. Angl / Xi vtx); Counts",
+ 100, 0.99, 1.0);
+ tQACandidates->Add(tV0toXiCosOfPointingAngle);
+
+ TH2F *th2Armenteros
+ = new TH2F("Armenteros",
+ "#alpha_{Arm}(casc. cand.) Vs Pt_{Arm}(casc. cand.); #alpha_{Arm} ; Pt_{Arm} (GeV/c)",
+ 140, -1.2, 1.2, 300, 0., 0.3);
+ tQACandidates->Add(th2Armenteros);
+
+ TH2F *th2TPCdEdxOfCascDghters
+ = new TH2F( "TPCdEdxOfCascDghters",
+ "TPC dE/dx of the cascade daughters; charge x || #vec{p}_{TPC inner wall}(Casc. daughter) || (GeV/c); TPC signal (ADC) ",
+ 200, -10.0, 10.0, 450, 0., 900.);
+ tQACandidates->Add(th2TPCdEdxOfCascDghters);
+
+ TH2F *th2MassVsPtAll
+ = new TH2F("MassVsPtAll",
+ "M_{candidates} vs Pt; Pt (GeV/c); M (GeV/c^{2})",
+ 100, 0., 10., fMassBins, fMinMass, fMaxMass);
+ tQACandidates->Add(th2MassVsPtAll);
+
+ TH1F *tDistToVtxZAfter
+ = new TH1F("DistToVtxZAfter",
+ "Distance to vtx z after propagation to vtx; z [cm]",
+ 100, -5., 5.);
+ tQACandidates->Add(tDistToVtxZAfter);
+
+ TH1F * tDistToVtxXYAfter
+ = new TH1F("DistToVtxXYAfter",
+ "Distance to vtx xy after propagation to vtx",
+ 500, 0., 50.);
+ tQACandidates->Add(tDistToVtxXYAfter);
+
+ TH2F *th2DistToVtxZBeforeVsAfter
+ = new TH2F("DistToVtxZBeforeVsAfter",
+ "Distance to vtx z before vs after propagation; Distance before [cm]; Distance after [cm]",
+ 500, -50., 50., 100, -5., 5.);
+ tQACandidates->Add(th2DistToVtxZBeforeVsAfter);
+
+ TH2F *th2DistToVtxXYBeforeVsAfter
+ = new TH2F("DistToVtxXYBeforeVsAfter",
+ "Distance to vtx xy before vs after propagation; Distance before [cm]; Distance after [cm]",
+ 500, 0., 50, 500, 0., 50);
+ tQACandidates->Add(th2DistToVtxXYBeforeVsAfter);
+
+ TH2F * th2PxBeforeVsAfter
+ = new TH2F("PxBeforeVsAfter",
+ "Px before vs after propagation; Px [GeV/c]; Px' [GeV/c]",
+ 200, -10., 10, 200, -10., 10.);
+ tQACandidates->Add(th2PxBeforeVsAfter);
+
+ TH2F * th2PyBeforeVsAfter
+ = new TH2F("PyBeforeVsAfter",
+ "Py before vs after propagation; Py [GeV/c]; Py' [GeV/c]",
+ 200, -10., 10, 200, -10., 10.);
+ tQACandidates->Add(th2PyBeforeVsAfter);
+
+ TH2F * th2PhiPosBeforeVsAfter
+ = new TH2F("PhiPosBeforeVsAfter",
+ "Phi for positively charged candidates before vs after propagation; #phi; #phi'",
+ 360, 0., 2.0*TMath::Pi(), 360, 0., 2.0*TMath::Pi());
+ tQACandidates->Add(th2PhiPosBeforeVsAfter);
+
+ TH2F *th2PhiNegBeforeVsAfter
+ = new TH2F("PhiNegBeforeVsAfter",
+ "Phi for negatively charged candidates before vs after propagation; #phi; #phi'",
+ 360, 0., 2.0*TMath::Pi(), 360, 0., 2.0*TMath::Pi());
+ tQACandidates->Add(th2PhiNegBeforeVsAfter);
+
+ fQAList->Add(tQACandidates);
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::NotifyRun()
+{
+}
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::UserExec(Option_t *)
+{
+ AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());
+ AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
+ Bool_t acceptEvent=kFALSE;
+ fCandidates->SetLast(-1);
+
+ if(fESD) {
+ // recorrecting VZERO (for pass 1 only)
+ /*
+ Float_t *vChCorr = new Float_t[64];
+ Float_t dummy;
+ AliESDUtils::GetCorrV0(fESD,dummy,NULL,vChCorr);
+ AliESDVZERO *vzero = (AliESDVZERO*) fESD->GetVZEROData();
+ vzero->SetMultiplicity( vChCorr );
+ delete [] vChCorr;
+ */
+ //
+
+ ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(0);
+
+ const AliVVertex *vtxGlb = fESD->GetPrimaryVertexTracks();
+ const AliVVertex *vtxSPD = fESD->GetPrimaryVertexSPD();
+ if(!vtxGlb || !vtxSPD) return;
+ if( fCutsEvent->IsSelected(fESD, 0) && (TMath::Abs(vtxSPD->GetZ()-vtxGlb->GetZ()) <= 0.5) ) {
+ acceptEvent = kTRUE;
+ ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);
+ ReadFromESDv0(fESD);
+ }
+ } else if(fAOD) {
+ const AliVVertex *vtxGlb = fAOD->GetPrimaryVertex();
+ const AliVVertex *vtxSPD = fAOD->GetPrimaryVertexSPD();
+ if(!vtxGlb || !vtxSPD) return;
+
+ ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(0);
+
+ if(fCutsEvent->IsSelected(fAOD, 0) && (TMath::Abs(vtxSPD->GetZ()-vtxGlb->GetZ()) <= 0.5) ) {
+ acceptEvent = kTRUE;
+ ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);
+ ReadFromAODv0(fAOD);
+ }
+
+
+ /*
+
+ AliAODHeader *aodHeader = fAOD->GetHeader();
+ if(!aodHeader) return;
+ AliCentrality *centrality = aodHeader->GetCentralityP();
+ if(!centrality) return;
+ Double_t cent = centrality->GetCentralityPercentile("V0M" );
+ Double_t cent1 = centrality->GetCentralityPercentile("TRK" );
+ if(TMath::Abs(cent-cent1) >= 5.) return;
+
+ if(cent<fMinCent||cent>=fMaxCent) return; //centrality cut
+
+ Double_t zvtx = fAOD->GetPrimaryVertex()->GetZ();
+ if(TMath::Abs(zvtx)>10.) return; //vertex cut
+
+ ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);
+ ReadFromAODv0(fAOD);
+ */
+ }
+
+ if(!acceptEvent) return;
+
+ ((TH1D*)((TList*)fQAList->FindObject("Events"))->FindObject("RFPTPC"))
+ ->Fill(fFlowEventTPC->GetNumberOfRPs() );
+ Double_t mult=0;
+ for(Int_t i=0; i != fFlowEventVZE->GetNumberOfRPs(); ++i) {
+ AliFlowTrackSimple *pTrack = fFlowEventVZE->GetTrack(i);
+ mult += pTrack->Weight();
+ }
+ ((TH1D*)((TList*)fQAList->FindObject("Events"))->FindObject("RFPVZE"))
+ ->Fill( mult );
+
+ // if(fDebug) printf("TPCevent %d | VZEevent %d\n",
+ // fFlowEventTPC->NumberOfTracks(),
+ // fFlowEventVZE->NumberOfTracks() );
+ AddCandidates();
+
+ PostData(1,fFlowEventTPC);
+ PostData(2,fFlowEventVZE);
+ PostData(3,fQAList);
+
+ return;
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::AddCandidates(){
+
+ // if(fDebug) printf("I received %d candidates\n",
+ // fCandidates->GetEntriesFast());
+ for(int iCand=0; iCand!=fCandidates->GetEntriesFast(); ++iCand ) {
+ AliFlowCandidateTrack *cand
+ = dynamic_cast<AliFlowCandidateTrack*>(fCandidates->At(iCand));
+ if(!cand) continue;
+ // if(fDebug)
+ // printf(" >Checking at candidate %d with %d daughters: mass %f\n",
+ // iCand, cand->GetNDaughters(), cand->Mass());
+
+ // untagging ===>
+ for(int iDau=0; iDau != cand->GetNDaughters(); ++iDau) {
+ // if(fDebug)
+ // printf(" >Daughter %d with fID %d", iDau, cand->GetIDDaughter(iDau));
+ for(int iRPs=0; iRPs != fFlowEventTPC->NumberOfTracks(); ++iRPs ) {
+ AliFlowTrack *iRP
+ = dynamic_cast<AliFlowTrack*>(fFlowEventTPC->GetTrack( iRPs ));
+ if (!iRP) continue;
+ if( !iRP->InRPSelection() ) continue;
+ if( cand->GetIDDaughter(iDau) == iRP->GetID() ) {
+ //if(fDebug) printf(" was in RP set");
+ iRP->SetForRPSelection(kFALSE);
+ fFlowEventTPC->SetNumberOfRPs( fFlowEventTPC->GetNumberOfRPs() -1 );
+ }
+ }
+ //if(fDebug) printf("\n");
+ }
+ // <=== untagging
+ cand->SetForPOISelection(kTRUE);
+ fFlowEventTPC->InsertTrack( ((AliFlowTrack*) cand) );
+ fFlowEventVZE->InsertTrack( ((AliFlowTrack*) cand) );
+ }
+
+ // if(fDebug) printf("TPCevent %d | VZEevent %d\n",
+ // fFlowEventTPC->NumberOfTracks(),
+ // fFlowEventVZE->NumberOfTracks() );
+
+}
+
+//______________________________________________________________________________
+void AliAnalysisTaskFlowCascade::ReadFromESDv0(AliESDEvent *fESD)
+{
+ //AliFlowTrackCuts* cutsPOI = new AliFlowTrackCuts("null_cuts");
+ //cutsPOI->SetParamType( fCutsRP->GetParamType() );
+ //cutsPOI->SetParamType( AliFlowTrackCuts::kGlobal );
+ // cutsPOI->SetPtRange(+1,-1); // select nothing
+ //cutsPOI->SetEtaRange(+1,-1); // select nothing VZERO
+
+ fCutsRPTPC->SetEvent(fESD,MCEvent());
+ fCutsRPVZE->SetEvent(fESD,MCEvent());
+
+ fCutsPOI->SetEvent(fESD,MCEvent());
+
+ fFlowEventTPC->Fill(fCutsRPTPC,fCutsPOI);
+ fFlowEventVZE->Fill(fCutsRPVZE,fCutsPOI);
+
+ Double_t trkPrimaryVtxPos[3] = {-100., -100., -100.};
+ Double_t bestPrimaryVtxPos[3] = {-100., -100., -100.};
+ int nCascades=fESD->GetNumberOfCascades();
+
+ const AliESDVertex *primaryTrackingESDVtx = fESD->GetPrimaryVertexTracks();
+ primaryTrackingESDVtx->GetXYZ(trkPrimaryVtxPos);
+
+ const AliESDVertex *primaryBestESDVtx = fESD->GetPrimaryVertex();
+ primaryBestESDVtx->GetXYZ(bestPrimaryVtxPos);
+
+ Double_t b = fESD->GetMagneticField();
+
+ for(int i = 0; i != nCascades; ++i) {
+
+ // Double_t trkPrimaryVtxRadius3D = -500.;
+ // Double_t bestPrimaryVtxRadius3D = -500.;
+ Double_t effMassXi = 0.;
+ Double_t chi2Xi = -1.;
+ Double_t dcaXiDaughters = -1.;
+ Double_t XiCosOfPointingAngle = -1.;
+ Double_t posXi[3] = {-1000., -1000., -1000.};
+ Double_t XiRadius = -1000.;
+
+ // Double_t innerWallMomCascDghters[3] = {-100., -100., -100.};
+ //Double_t tpcSignalCascDghters[3] = {-100., -100., -100.};
+
+ Double_t invMassLambdaAsCascDghter = 0.;
+ Double_t V0Chi2Xi = -1.;
+ Double_t dcaV0DaughtersXi = -1.;
+
+ Double_t dcaBachToPrimaryVtxXi = -1.;
+ Double_t dcaV0ToPrimaryVtxXi = -1.;
+ Double_t dcaPosToPrimaryVtxXi = -1.;
+ Double_t dcaNegToPrimaryVtxXi = -1.;
+ Double_t V0CosOfPointingAngleXi = -1.;
+ Double_t posV0Xi[3] = {-1000., -1000., -1000.};
+ Double_t V0RadiusXi = -1000.;
+ Double_t V0quality = 0.;
+
+ Double_t invMassXiMinus = 0.;
+ Double_t invMassXiPlus = 0.;
+ Double_t invMassOmegaMinus = 0.;
+ Double_t invMassOmegaPlus = 0.;
+
+ /*
+ Bool_t isPosInXiProton = kFALSE;
+ Bool_t isPosInXiPion = kFALSE;
+ Bool_t isPosInOmegaProton = kFALSE;
+ Bool_t isPosInOmegaPion = kFALSE;
+
+ Bool_t isNegInXiProton = kFALSE;
+ Bool_t isNegInXiPion = kFALSE;
+ Bool_t isNegInOmegaProton = kFALSE;
+ Bool_t isNegInOmegaPion = kFALSE;
+
+ Bool_t isBachelorKaon = kFALSE;
+ Bool_t isBachelorPion = kFALSE;
+ */
+
+ Bool_t isBachelorKaonForTPC = kFALSE;
+ Bool_t isBachelorPionForTPC = kFALSE;
+ Bool_t isNegPionForTPC = kFALSE;
+ Bool_t isPosPionForTPC = kFALSE;
+ Bool_t isNegProtonForTPC = kFALSE;
+ Bool_t isPosProtonForTPC = kFALSE;
+
+ Double_t XiPx = 0., XiPy = 0., XiPz = 0.;
+ Double_t XiPt = 0.;
+ Double_t XiPtot = 0.;
+
+ Double_t bachPx = 0., bachPy = 0., bachPz = 0.;
+ Double_t bachPt = 0.;
+ Double_t bachPtot = 0.;
+
+ //Short_t chargeXi = -2;
+ Double_t V0toXiCosOfPointingAngle = 0.;
+
+ Double_t rapXi = -20.;
+ Double_t rapOmega = -20.;
+ Double_t phi = 6.3;
+ Double_t alphaXi = -200.;
+ Double_t ptArmXi = -200.;
+ // TLorentzVector lv1, lv2, lv3, lv12, lvXi;
+
+ Double_t distToVtxZBefore = -999.;
+ Double_t distToVtxZAfter = -999.;
+ Double_t distToVtxXYBefore = -999.;
+ Double_t distToVtxXYAfter = -999.;
+ Double_t XiPAfter[3] = {-999., -999., -999.};
+ Double_t phiAfter = -999.;
+
+ AliESDcascade *xi = fESD->GetCascade(i);
+ if(!xi) continue;
+
+ if(xi->Charge()<0)
+ xi->ChangeMassHypothesis(V0quality, 3312); // Xi- hypothesis
+ else if(xi->Charge() > 0)
+ xi->ChangeMassHypothesis(V0quality, -3312);
+ else continue;
+
+ effMassXi = xi->GetEffMassXi();
+ chi2Xi = xi->GetChi2Xi();
+ dcaXiDaughters = xi->GetDcaXiDaughters();
+ XiCosOfPointingAngle
+ = xi->GetCascadeCosineOfPointingAngle(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]);
+ xi->GetXYZcascade(posXi[0], posXi[1], posXi[2]);
+ XiRadius = TMath::Sqrt(posXi[0]*posXi[0]
+ +posXi[1]*posXi[1]
+ +posXi[2]*posXi[2]);
+
+ UInt_t idxPosXi = (UInt_t)TMath::Abs(xi->GetPindex());
+ UInt_t idxNegXi = (UInt_t)TMath::Abs(xi->GetNindex());
+ UInt_t idxBach = (UInt_t)TMath::Abs(xi->GetBindex());
+
+ if(idxBach == idxPosXi || idxBach == idxNegXi) continue;
+
+ AliESDtrack *pTrkXi = fESD->GetTrack(idxPosXi);
+ AliESDtrack *nTrkXi = fESD->GetTrack(idxNegXi);
+ AliESDtrack *bTrkXi = fESD->GetTrack(idxBach);
+
+ if( !pTrkXi || !nTrkXi || !bTrkXi ) continue;
+
+ if( !fCutsDau->IsSelected(pTrkXi)
+ || !fCutsDau->IsSelected(nTrkXi)
+ || !fCutsDau->IsSelected(bTrkXi) ) continue;
+
+ const AliExternalTrackParam *pExtTrk = pTrkXi->GetInnerParam();
+ const AliExternalTrackParam *nExtTrk = nTrkXi->GetInnerParam();
+ const AliExternalTrackParam *bExtTrk = bTrkXi->GetInnerParam();
+
+ if(pExtTrk && pTrkXi->IsOn(AliESDtrack::kTPCin)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(pExtTrk->GetP()*pExtTrk->Charge(), pTrkXi->GetTPCsignal());
+ }
+ if(nExtTrk && nTrkXi->IsOn(AliESDtrack::kTPCin)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(nExtTrk->GetP()*nExtTrk->Charge(), nTrkXi->GetTPCsignal());
+ }
+ if(bExtTrk && bTrkXi->IsOn(AliESDtrack::kTPCin)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(bExtTrk->GetP()*bExtTrk->Charge(), bTrkXi->GetTPCsignal());
+ }
+
+ invMassLambdaAsCascDghter = xi->GetEffMass(); // from V0
+ dcaV0DaughtersXi = xi->GetDcaV0Daughters();
+ V0Chi2Xi = xi->GetChi2V0();
+ V0CosOfPointingAngleXi
+ = xi->GetV0CosineOfPointingAngle(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]);
+ dcaV0ToPrimaryVtxXi = xi->GetD(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]);
+ dcaBachToPrimaryVtxXi = TMath::Abs(bTrkXi->GetD(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]));
+
+ //V0
+ xi->GetXYZ(posV0Xi[0], posV0Xi[1], posV0Xi[2]);
+ V0RadiusXi = TMath::Sqrt(posV0Xi[0]*posV0Xi[0]
+ +posV0Xi[1]*posV0Xi[1]
+ +posV0Xi[2]*posV0Xi[2]);
+ dcaPosToPrimaryVtxXi = TMath::Abs(pTrkXi->GetD(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]));
+ dcaNegToPrimaryVtxXi = TMath::Abs(nTrkXi->GetD(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]));
+
+ //apply cuts
+ //if(XiRadius < 0.9 || XiRadius > 100.) continue;
+ //if(dcaXiDaughters > 0.2) continue;
+ //if(XiCosOfPointingAngle < 0.99) continue;
+ //if(dcaV0ToPrimaryVtxXi < 0.03) continue;
+ //if(dcaBachToPrimaryVtxXi < 0.01) continue;
+
+ if(dcaXiDaughters > fCascadeCuts[0]) continue;
+ if(XiCosOfPointingAngle < fCascadeCuts[1]) continue;
+ if(dcaV0ToPrimaryVtxXi < fCascadeCuts[2]) continue;
+ if(dcaBachToPrimaryVtxXi < fCascadeCuts[3]) continue;
+
+ //V0 mass cut?
+ // if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > 0.01) continue;
+ if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > fCascadeCuts[7])
+ continue;
+
+ //if(dcaV0DaughtersXi > 1.) continue;
+ //if(V0CosOfPointingAngleXi > 0.9999) continue;
+ //if(dcaPosToPrimaryVtxXi < 0.1) continue;
+ //if(dcaNegToPrimaryVtxXi < 0.1) continue;
+
+ if(dcaV0DaughtersXi > fCascadeCuts[4]) continue;
+ if(V0CosOfPointingAngleXi > fCascadeCuts[5]) continue;
+ if(dcaPosToPrimaryVtxXi < fCascadeCuts[6]) continue;
+ if(dcaNegToPrimaryVtxXi < fCascadeCuts[6]) continue;
+
+ //if(V0RadiusXi < 1.0 || V0RadiusXi > 100) continue;
+
+ //other cuts?
+ // change mass hypothesis to cover all the possibilities
+ if(bTrkXi->Charge()<0){
+ V0quality = 0.;
+ xi->ChangeMassHypothesis(V0quality, 3312); //Xi- hyp.
+ invMassXiMinus = xi->GetEffMassXi();
+
+ V0quality = 0.;
+ xi->ChangeMassHypothesis(V0quality, 3334); //Omega- hyp.
+ invMassOmegaMinus = xi->GetEffMassXi();
+
+ V0quality = 0.;
+ xi->ChangeMassHypothesis(V0quality, 3312); //back to default hyp.
+ }
+
+ if(bTrkXi->Charge() > 0){
+ V0quality = 0.;
+ xi->ChangeMassHypothesis(V0quality, -3312); //anti-Xi- hyp.
+ invMassXiPlus = xi->GetEffMassXi();
+
+ V0quality = 0.;
+ xi->ChangeMassHypothesis(V0quality, -3334); //anti-Omega- hyp.
+ invMassOmegaPlus = xi->GetEffMassXi();
+
+ V0quality = 0.;
+ xi->ChangeMassHypothesis(V0quality, -3312); //back to default hyp.
+ }
+
+ //PID on the daughter tracks
+ /*
+ //A - Combined PID
+ //Resonable guess the priors for the cascade track sample
+ //(e, mu, pi, K, p)
+ Double_t priorsGuessXi[5] = {0, 0, 2, 0, 1};
+ Double_t priorsGuessOmega[5] = {0, 0, 1, 1, 1};
+
+ //Combined bachelor-daughter PID
+ AliPID pidXi;
+ pidXi.SetPriors(priorsGuessXi);
+ AliPID pidOmega;
+ pidOmega.SetPriors(priorsGuessOmega);
+
+ if(pTrkXi->IsOn(AliESDtrack::kESDpid)){// combined PID exists
+ Double_t r[10] = {0.};
+ pTrkXi->GetESDpid(r);
+ pidXi.SetProbabilities(r);
+ pidOmega.SetProbabilities(r);
+
+ //Check if the V0 postive track is proton (case for Xi-)
+ Double_t pProton = pidXi.GetProbability(AliPID::kProton);
+ if(pProton > pidXi.GetProbability(AliPID::kElectron)
+ && pProton > pidXi.GetProbability(AliPID::kMuon)
+ && pProton > pidXi.GetProbability(AliPID::kPion)
+ && pProton > pidXi.GetProbability(AliPID::kKaon))
+ isPosInXiProton = kTRUE;
+
+ //Check if the V0 postive track is a pi+ (case for Xi+)
+ Double_t pPion = pidXi.GetProbability(AliPID::kPion);
+ if(pPion > pidXi.GetProbability(AliPID::kElectron)
+ && pPion > pidXi.GetProbability(AliPID::kMuon)
+ && pPion > pidXi.GetProbability(AliPID::kKaon)
+ && pPion > pidXi.GetProbability(AliPID::kProton))
+ isPosInXiPion = kTRUE;
+ // Check if the V0 positive track is a proton (case for Omega-)
+ pProton = pidOmega.GetProbability(AliPID::kProton);
+ if(pProton > pidOmega.GetProbability(AliPID::kElectron)
+ && pProton > pidOmega.GetProbability(AliPID::kMuon)
+ && pProton > pidOmega.GetProbability(AliPID::kPion)
+ && pProton > pidOmega.GetProbability(AliPID::kKaon))
+ isPosInOmegaProton = kTRUE;
+
+ // Check if the V0 positive track is a pi+ (case for Omega+)
+ pPion = pidOmega.GetProbability(AliPID::kPion);
+ if(pPion > pidOmega.GetProbability(AliPID::kElectron)
+ && pPion > pidOmega.GetProbability(AliPID::kMuon)
+ && pPion > pidOmega.GetProbability(AliPID::kKaon)
+ && pPion > pidOmega.GetProbability(AliPID::kProton))
+ isPosInOmegaPion = kTRUE;
+ }
+
+ //Combined V0-negative-daughter PID
+ pidXi.SetPriors(priorsGuessXi);
+ pidOmega.SetPriors(priorsGuessOmega);
+ if(nTrkXi->IsOn(AliESDtrack::kESDpid)){
+ Double_t r[10] = {0.};
+ nTrkXi->GetESDpid(r);
+ pidXi.SetProbabilities(r);
+ pidOmega.SetProbabilities(r);
+
+ // Check if the V0 negative track is a pi- (case for Xi-)
+ Double_t pPion = pidXi.GetProbability(AliPID::kPion);
+ if(pPion > pidXi.GetProbability(AliPID::kElectron)
+ && pPion > pidXi.GetProbability(AliPID::kMuon)
+ && pPion > pidXi.GetProbability(AliPID::kKaon)
+ && pPion > pidXi.GetProbability(AliPID::kProton))
+ isNegInXiPion = kTRUE;
+
+ // Check if the V0 negative track is an anti-proton (case for Xi+)
+ Double_t pProton = pidXi.GetProbability(AliPID::kProton);
+ if(pProton > pidXi.GetProbability(AliPID::kElectron)
+ && pProton > pidXi.GetProbability(AliPID::kMuon)
+ && pProton > pidXi.GetProbability(AliPID::kPion)
+ && pProton > pidXi.GetProbability(AliPID::kKaon))
+ isNegInXiProton = kTRUE;
+
+ // Check if the V0 negative track is a pi- (case for Omega-)
+ pPion = pidOmega.GetProbability(AliPID::kPion);
+ if(pPion > pidOmega.GetProbability(AliPID::kElectron)
+ && pPion > pidOmega.GetProbability(AliPID::kMuon)
+ && pPion > pidOmega.GetProbability(AliPID::kKaon)
+ && pPion > pidOmega.GetProbability(AliPID::kProton))
+ isNegInOmegaPion = kTRUE;
+
+ // Check if the V0 negative track is an anti-proton (case for Omega+)
+ pProton = pidOmega.GetProbability(AliPID::kProton);
+ if(pProton > pidOmega.GetProbability(AliPID::kElectron)
+ && pProton > pidOmega.GetProbability(AliPID::kMuon)
+ && pProton > pidOmega.GetProbability(AliPID::kPion)
+ && pProton > pidOmega.GetProbability(AliPID::kKaon))
+ isNegInOmegaProton = kTRUE;
+
+ }
+
+ // Combined bachelor PID
+ pidXi.SetPriors(priorsGuessXi);
+ pidOmega.SetPriors(priorsGuessOmega);
+ if(bTrkXi->IsOn(AliESDtrack::kESDpid)){//Combined PID exists
+ Double_t r[10] = {0.};
+ bTrkXi->GetESDpid(r);
+ pidXi.SetProbabilities(r);
+ pidOmega.SetProbabilities(r);
+
+ //Check if the bachelor track is a pion
+ Double_t pPion = pidXi.GetProbability(AliPID::kPion);
+ if(pPion > pidXi.GetProbability(AliPID::kElectron)
+ && pPion > pidXi.GetProbability(AliPID::kMuon)
+ && pPion > pidXi.GetProbability(AliPID::kKaon)
+ && pPion > pidXi.GetProbability(AliPID::kProton))
+ isBachelorPion = kTRUE;
+
+ // Check if the bachelor track is a kaon
+ Double_t pKaon = pidOmega.GetProbability(AliPID::kKaon);
+ if(pKaon > pidOmega.GetProbability(AliPID::kElectron)
+ && pKaon > pidOmega.GetProbability(AliPID::kMuon)
+ && pKaon > pidOmega.GetProbability(AliPID::kPion)
+ && pKaon > pidOmega.GetProbability(AliPID::kProton))
+ isBachelorKaon = kTRUE;
+ }
+ */
+
+
+ //B - TPC PID: 3-sigma bands on Bethe-Bloch curve
+ //Bachelor
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kKaon))<3.)
+ isBachelorKaonForTPC = kTRUE;
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kPion))<3.)
+ isBachelorPionForTPC = kTRUE;
+
+ //Negative V0 daughter
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kPion))<3.)
+ isNegPionForTPC = kTRUE;
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kProton))<3.)
+ isNegProtonForTPC = kTRUE;
+
+ //Positive V0 daughter
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kPion))<3.)
+ isPosPionForTPC = kTRUE;
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kProton))<3.)
+ isPosProtonForTPC = kTRUE;
+
+
+ //Extra QA information
+ xi->GetPxPyPz(XiPx, XiPy, XiPz);
+ XiPt = TMath::Sqrt(XiPx*XiPx + XiPy*XiPy);
+ XiPtot= TMath::Sqrt(XiPx*XiPx + XiPy*XiPy + XiPz*XiPz);
+
+ XiPAfter[0] = XiPx;
+ XiPAfter[1] = XiPy;
+ XiPAfter[2] = XiPz;
+
+ xi->GetBPxPyPz(bachPx, bachPy, bachPz);
+ bachPt = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy);
+ bachPtot = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy + bachPz*bachPz);
+
+ //chargeXi = xi->Charge();
+
+ V0toXiCosOfPointingAngle
+ = xi->GetV0CosineOfPointingAngle(posXi[0], posXi[1], posXi[2]);
+ rapXi = xi->RapXi();
+ rapOmega = xi->RapOmega();
+ phi = xi->Phi();
+ alphaXi = xi->AlphaXi();
+ ptArmXi = xi->PtArmXi();
+
+ distToVtxZBefore = posXi[2]-bestPrimaryVtxPos[2];
+ distToVtxXYBefore
+ = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])
+ *(posXi[0] - bestPrimaryVtxPos[0])
+ +(posXi[1] - bestPrimaryVtxPos[1])
+ *(posXi[1] - bestPrimaryVtxPos[1]));
+
+ //propagation to the best primary vertex to determine the momentum
+ Propagate(bestPrimaryVtxPos, posXi, XiPAfter, b, xi->Charge());
+ distToVtxZAfter = posXi[2] - bestPrimaryVtxPos[2];
+ distToVtxXYAfter = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])
+ *(posXi[0] - bestPrimaryVtxPos[0])
+ +(posXi[1] - bestPrimaryVtxPos[1])
+ *(posXi[1] - bestPrimaryVtxPos[1]));
+ phiAfter = TMath::Pi() + TMath::ATan2(-XiPAfter[1],-XiPAfter[0]);
+
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZAfter"))->Fill(distToVtxZAfter);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYAfter"))->Fill(distToVtxXYAfter);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZBeforeVsAfter"))->Fill(distToVtxZBefore, distToVtxZAfter);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYBeforeVsAfter"))->Fill(distToVtxXYBefore, distToVtxXYAfter);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PxBeforeVsAfter"))->Fill(XiPx, XiPAfter[0]);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PyBeforeVsAfter"))->Fill(XiPy, XiPAfter[1]);
+ if(xi->Charge()>0)
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiPosBeforeVsAfter"))->Fill(phi, phiAfter);
+ else if(xi->Charge()<0)
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiNegBeforeVsAfter"))->Fill(phi, phiAfter);
+
+
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Chi2Xi"))->Fill(chi2Xi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaXiDaughters"))->Fill(dcaXiDaughters);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaBachToPrimVertex"))->Fill(dcaBachToPrimaryVtxXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiCosOfPointingAngle"))->Fill(XiCosOfPointingAngle);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiRadius"))->Fill(XiRadius);
+
+ //V0
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassLambdaAsCascDghter"))->Fill(invMassLambdaAsCascDghter);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0Chi2Xi"))->Fill(V0Chi2Xi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0DaughtersXi"))->Fill(dcaV0DaughtersXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0CosOfPointingAngleXi"))->Fill(V0CosOfPointingAngleXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0RadiusXi"))->Fill(V0RadiusXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0ToPrimVertexXi"))->Fill(dcaV0ToPrimaryVtxXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaPosToPrimVertexXi"))->Fill(dcaPosToPrimaryVtxXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaNegToPrimVertexXi"))->Fill(dcaNegToPrimaryVtxXi);
+
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0toXiCosOfPointingAngle"))->Fill(V0toXiCosOfPointingAngle);
+
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Armenteros"))->Fill(alphaXi, ptArmXi);
+
+ //PID cuts with TPC cuts
+ if(xi->Charge() < 0){
+ if(isPosProtonForTPC
+ && isNegPionForTPC){
+
+ switch(fSpecie) {
+ case 0:
+ if(isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiMinus);
+
+ //candidate inserting
+ MakeTrack(invMassXiMinus, XiPt, /*xi->Phi()*/
+ phiAfter, xi->Eta(), pTrkXi->GetID(),
+ nTrkXi->GetID(), bTrkXi->GetID());
+ }
+ break;
+
+ case 1:
+ if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaMinus);
+ MakeTrack(invMassOmegaMinus, XiPt, /*xi->Phi()*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }
+ break;
+ }
+ }
+ }
+
+ if(xi->Charge() > 0){
+ if(isNegProtonForTPC
+ && isPosPionForTPC){
+
+ switch (fSpecie){
+ case 0:
+ if(isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiPlus);
+
+ //candidate inserting
+ MakeTrack(invMassXiPlus, XiPt, /*xi->Phi()*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }
+ break;
+
+ case 1:
+ if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaPlus);
+ MakeTrack(invMassOmegaPlus, XiPt, /*xi->Phi()*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }
+ break;
+ }
+ }
+ }
+
+ }
+
+ return;
+}
+
+//______________________________________________________________________________
+void AliAnalysisTaskFlowCascade::ReadFromAODv0(AliAODEvent *fAOD)
+{
+
+ fCutsRPTPC->SetEvent(fAOD, MCEvent());
+ fCutsRPVZE->SetEvent(fAOD, MCEvent());
+ fCutsPOI->SetEvent(fAOD, MCEvent());
+ fFlowEventTPC->Fill(fCutsRPTPC, fCutsPOI);
+ fFlowEventVZE->Fill(fCutsRPVZE, fCutsPOI);
+
+ // Double_t trkPrimaryVtxPos[3] = {-100., -100., -100.};
+ Double_t bestPrimaryVtxPos[3] = {-100., -100., -100.};
+
+ Double_t b = fAOD->GetMagneticField();
+
+ int nCascades=fAOD->GetNumberOfCascades();
+ const AliAODVertex *primaryBestAODVtx = fAOD->GetPrimaryVertex();
+ primaryBestAODVtx->GetXYZ(bestPrimaryVtxPos);
+
+ // calculation part dedicated to Xi vertices
+ for(Int_t iXi = 0; iXi < nCascades; iXi++){
+ Double_t effMassXi = 0.;
+ Double_t chi2Xi = -1.;
+ Double_t dcaXiDaughters = -1.;
+ Double_t XiCosOfPointingAngle = -1.;
+ Double_t posXi[3] = {-1000., -1000., -1000.};
+ Double_t XiRadius = -1000.;
+
+ Double_t invMassLambdaAsCascDghter = 0.;
+ Double_t V0Chi2Xi = -1.;
+ Double_t dcaV0DaughtersXi = -1.;
+
+ Double_t dcaBachToPrimaryVtxXi = -1.;
+ Double_t dcaV0ToPrimaryVtxXi = -1.;
+ Double_t dcaPosToPrimaryVtxXi = -1.;
+ Double_t dcaNegToPrimaryVtxXi = -1.;
+ Double_t V0CosOfPointingAngleXi = -1.;
+ Double_t posV0Xi[3] = {-1000., -1000., -1000.};
+ Double_t V0RadiusXi = -1000.;
+ // Double_t V0quality = 0.;
+
+ Double_t invMassXiMinus = 0.;
+ Double_t invMassXiPlus = 0.;
+ Double_t invMassOmegaMinus = 0.;
+ Double_t invMassOmegaPlus = 0.;
+
+ /*
+ Bool_t isPosInXiProton = kFALSE;
+ Bool_t isPosInXiPion = kFALSE;
+ Bool_t isPosInOmegaProton = kFALSE;
+ Bool_t isPosInOmegaPion = kFALSE;
+
+ Bool_t isNegInXiProton = kFALSE;
+ Bool_t isNegInXiPion = kFALSE;
+ Bool_t isNegInOmegaProton = kFALSE;
+ Bool_t isNegInOmegaPion = kFALSE;
+
+ Bool_t isBachelorKaon = kFALSE;
+ Bool_t isBachelorPion = kFALSE;
+ */
+
+
+ Bool_t isBachelorKaonForTPC = kFALSE;
+ Bool_t isBachelorPionForTPC = kFALSE;
+ Bool_t isNegPionForTPC = kFALSE;
+ Bool_t isPosPionForTPC = kFALSE;
+ Bool_t isNegProtonForTPC = kFALSE;
+ Bool_t isPosProtonForTPC = kFALSE;
+
+ Double_t XiPx = 0., XiPy = 0., XiPz = 0.;
+ Double_t XiPt = 0.;
+ Double_t XiPtot = 0.;
+
+ Double_t bachPx = 0., bachPy = 0., bachPz = 0.;
+ Double_t bachPt = 0.;
+ Double_t bachPtot = 0.;
+
+ //Short_t chargeXi = -2;
+ Double_t V0toXiCosOfPointingAngle = 0.;
+
+ Double_t rapXi = -20.;
+ Double_t rapOmega = -20.;
+ Double_t phi = 6.3;
+ Double_t alphaXi = -200.;
+ Double_t ptArmXi = -200.;
+
+ Double_t distToVtxZBefore = -999.;
+ Double_t distToVtxZAfter = -999.;
+ Double_t distToVtxXYBefore = -999.;
+ Double_t distToVtxXYAfter = -999.;
+ Double_t XiPAfter[3] = {-999., -999., -999.};
+ Double_t phiAfter = -999.;
+
+ const AliAODcascade *xi = fAOD->GetCascade(iXi);
+ if (!xi) continue;
+
+ effMassXi = xi->MassXi(); //default working hypothesis: Xi- decay
+ chi2Xi = xi->Chi2Xi();
+ dcaXiDaughters = xi->DcaXiDaughters();
+ XiCosOfPointingAngle = xi->CosPointingAngleXi(bestPrimaryVtxPos[0],
+ bestPrimaryVtxPos[1],
+ bestPrimaryVtxPos[2]);
+ posXi[0] = xi->DecayVertexXiX();
+ posXi[1] = xi->DecayVertexXiY();
+ posXi[2] = xi->DecayVertexXiZ();
+ XiRadius = TMath::Sqrt(posXi[0]*posXi[0]
+ +posXi[1]*posXi[1]
+ +posXi[2]*posXi[2]);
+
+ AliAODTrack *pTrkXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(0) );
+ AliAODTrack *nTrkXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(1) );
+ AliAODTrack *bTrkXi
+ = dynamic_cast<AliAODTrack*>( xi->GetDecayVertexXi()->GetDaughter(0) );
+
+ if(!pTrkXi || !nTrkXi || !bTrkXi) continue;
+
+ UInt_t idxPosXi = (UInt_t) TMath::Abs( pTrkXi->GetID() );
+ UInt_t idxNegXi = (UInt_t) TMath::Abs( nTrkXi->GetID() );
+ UInt_t idxBach = (UInt_t) TMath::Abs( bTrkXi->GetID() );
+
+ if(idxBach == idxNegXi || idxBach == idxPosXi) continue;
+
+ if( !fCutsDau->IsSelected(pTrkXi)
+ || !fCutsDau->IsSelected(nTrkXi)
+ || !fCutsDau->IsSelected(bTrkXi) ) continue;
+
+
+ if(pTrkXi->IsOn(AliESDtrack::kTPCin)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(pTrkXi->P()*pTrkXi->Charge(), pTrkXi->GetTPCsignal());
+ }
+ if( nTrkXi->IsOn(AliESDtrack::kTPCin) ){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(nTrkXi->P()*nTrkXi->Charge(), nTrkXi->GetTPCsignal());
+ }
+ if(bTrkXi->IsOn(AliESDtrack::kTPCin)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(bTrkXi->P()*bTrkXi->Charge(), bTrkXi->GetTPCsignal());
+ }
+
+ if(xi->ChargeXi() < 0)
+ invMassLambdaAsCascDghter = xi->MassLambda();
+ else
+ invMassLambdaAsCascDghter = xi->MassAntiLambda();
+
+ dcaV0DaughtersXi = xi->DcaV0Daughters();
+ V0Chi2Xi = xi->Chi2V0();
+ V0CosOfPointingAngleXi
+ = xi->CosPointingAngle(bestPrimaryVtxPos);
+ dcaV0ToPrimaryVtxXi = xi->DcaV0ToPrimVertex();
+ dcaBachToPrimaryVtxXi = xi->DcaBachToPrimVertex();
+
+ //V0
+ posV0Xi[0] = xi->DecayVertexV0X();
+ posV0Xi[1] = xi->DecayVertexV0Y();
+ posV0Xi[2] = xi->DecayVertexV0Z();
+ V0RadiusXi = TMath::Sqrt(posV0Xi[0]*posV0Xi[0]
+ +posV0Xi[1]*posV0Xi[1]
+ +posV0Xi[2]*posV0Xi[2]);
+ dcaPosToPrimaryVtxXi = xi->DcaPosToPrimVertex();
+ dcaNegToPrimaryVtxXi = xi->DcaNegToPrimVertex();
+
+ //apply cuts ?
+ // if(XiRadius < 1. || XiRadius > 100.) continue;
+ //if(dcaXiDaughters > 0.1) continue;
+ //if(XiCosOfPointingAngle < 0.999) continue;
+ //if(dcaV0ToPrimaryVtxXi < 0.05) continue;
+ //if(dcaBachToPrimaryVtxXi < 0.03) continue;
+
+ if(dcaXiDaughters > fCascadeCuts[0]) continue;
+ if(XiCosOfPointingAngle < fCascadeCuts[1]) continue;
+ if(dcaV0ToPrimaryVtxXi < fCascadeCuts[2]) continue;
+ if(dcaBachToPrimaryVtxXi < fCascadeCuts[3]) continue;
+
+ //V0 mass cut?
+ //if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > 0.006) continue;
+ if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > fCascadeCuts[7])
+ continue;
+
+ //if(dcaV0DaughtersXi > 1.) continue;
+ //if(V0CosOfPointingAngleXi > 0.9999) continue;
+ //if(dcaPosToPrimaryVtxXi < 0.1) continue;
+ //if(dcaNegToPrimaryVtxXi < 0.1) continue;
+ if(dcaV0DaughtersXi > fCascadeCuts[4]) continue;
+ if(V0CosOfPointingAngleXi > fCascadeCuts[5]) continue;
+ if(dcaPosToPrimaryVtxXi < fCascadeCuts[6]) continue;
+ if(dcaNegToPrimaryVtxXi < fCascadeCuts[6]) continue;
+
+ // if(V0RadiusXi < 1.0 || V0RadiusXi > 100) continue;
+
+ //other cuts?
+
+
+ //???
+ if(xi->ChargeXi()<0){
+ invMassXiMinus = xi->MassXi();
+ invMassOmegaMinus = xi->MassOmega();
+ }else{
+ invMassXiPlus = xi->MassXi();
+ invMassOmegaPlus = xi->MassOmega();
+ }
+
+ /*
+ if(pTrkXi->GetMostProbablePID() == AliAODTrack::kProton) {
+ isPosInXiProton = kTRUE;
+ isPosInOmegaProton = kTRUE;
+ }
+ if(pTrkXi->GetMostProbablePID() == AliAODTrack::kPion){
+ isPosInXiPion = kTRUE;
+ isPosInOmegaPion = kTRUE;
+ }
+
+ if(nTrkXi->GetMostProbablePID() == AliAODTrack::kPion){
+ isNegInXiPion = kTRUE;
+ isNegInOmegaPion = kTRUE;
+ }
+ if(nTrkXi->GetMostProbablePID() == AliAODTrack::kProton){
+ isNegInXiProton = kTRUE;
+ isNegInOmegaProton = kTRUE;
+ }
+
+ if(bTrkXi->GetMostProbablePID() == AliAODTrack::kPion)
+ isBachelorPion = kTRUE;
+ if(bTrkXi->GetMostProbablePID() == AliAODTrack::kKaon)
+ isBachelorKaon = kTRUE;
+ */
+
+ //PID with TPC only:
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kKaon))<3.)
+ isBachelorKaonForTPC = kTRUE;
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kPion))<3.)
+ isBachelorPionForTPC = kTRUE;
+
+ //Negative V0 daughter
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kPion))<3.)
+ isNegPionForTPC = kTRUE;
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kProton))<3.)
+ isNegProtonForTPC = kTRUE;
+
+ //Positive V0 daughter
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kPion))<3.)
+ isPosPionForTPC = kTRUE;
+ if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kProton))<3.)
+ isPosProtonForTPC = kTRUE;
+
+ //Extra QA information
+ XiPx = xi->MomXiX();
+ XiPy = xi->MomXiY();
+ XiPz = xi->MomXiZ();
+ XiPt = TMath::Sqrt(XiPx*XiPx + XiPy*XiPy);
+ XiPtot= TMath::Sqrt(XiPx*XiPx + XiPy*XiPy + XiPz*XiPz);
+
+ bachPx = xi->MomBachX();
+ bachPy = xi->MomBachY();
+ bachPz = xi->MomBachZ();
+
+ bachPt = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy);
+ bachPtot = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy + bachPz*bachPz);
+
+ V0toXiCosOfPointingAngle = xi->CosPointingAngle( xi->GetDecayVertexXi() );
+
+ rapXi = xi->RapXi();
+ rapOmega = xi->RapOmega();
+ phi = xi->Phi();
+ alphaXi = xi->AlphaXi();
+ ptArmXi = xi->PtArmXi();
+
+ distToVtxZBefore = posXi[2]-bestPrimaryVtxPos[2];
+ distToVtxXYBefore
+ = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])
+ *(posXi[0] - bestPrimaryVtxPos[0])
+ +(posXi[1] - bestPrimaryVtxPos[1])
+ *(posXi[1] - bestPrimaryVtxPos[1]));
+
+
+ XiPAfter[0] = XiPx;
+ XiPAfter[1] = XiPy;
+ XiPAfter[2] = XiPz;
+ //propagation to the best primary vertex to determine the momentum
+ Propagate(bestPrimaryVtxPos, posXi, XiPAfter, b, xi->ChargeXi());
+ distToVtxZAfter = posXi[2] - bestPrimaryVtxPos[2];
+ distToVtxXYAfter = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])
+ *(posXi[0] - bestPrimaryVtxPos[0])
+ +(posXi[1] - bestPrimaryVtxPos[1])
+ *(posXi[1] - bestPrimaryVtxPos[1]));
+ phiAfter = TMath::Pi() + TMath::ATan2(-XiPAfter[1],-XiPAfter[0]);
+
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZAfter"))->Fill(distToVtxZAfter);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYAfter"))->Fill(distToVtxXYAfter);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZBeforeVsAfter"))->Fill(distToVtxZBefore, distToVtxZAfter);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYBeforeVsAfter"))->Fill(distToVtxXYBefore, distToVtxXYAfter);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PxBeforeVsAfter"))->Fill(XiPx, XiPAfter[0]);
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PyBeforeVsAfter"))->Fill(XiPy, XiPAfter[1]);
+ if(xi->ChargeXi()>0)
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiPosBeforeVsAfter"))->Fill(phi, phiAfter);
+ else if(xi->ChargeXi()<0)
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiNegBeforeVsAfter"))->Fill(phi, phiAfter);
+
+ //for default hypothesis
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Chi2Xi"))->Fill(chi2Xi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaXiDaughters"))->Fill(dcaXiDaughters);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaBachToPrimVertex"))->Fill(dcaBachToPrimaryVtxXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiCosOfPointingAngle"))->Fill(XiCosOfPointingAngle);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiRadius"))->Fill(XiRadius);
+
+ //V0
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassLambdaAsCascDghter"))->Fill(invMassLambdaAsCascDghter);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0Chi2Xi"))->Fill(V0Chi2Xi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0DaughtersXi"))->Fill(dcaV0DaughtersXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0CosOfPointingAngleXi"))->Fill(V0CosOfPointingAngleXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0RadiusXi"))->Fill(V0RadiusXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0ToPrimVertexXi"))->Fill(dcaV0ToPrimaryVtxXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaPosToPrimVertexXi"))->Fill(dcaPosToPrimaryVtxXi);
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaNegToPrimVertexXi"))->Fill(dcaNegToPrimaryVtxXi);
+
+ ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0toXiCosOfPointingAngle"))->Fill(V0toXiCosOfPointingAngle);
+
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Armenteros"))->Fill(alphaXi, ptArmXi);
+
+ //with PID cuts
+ if(xi->ChargeXi()<0){
+ if(isPosProtonForTPC && isNegPionForTPC){
+ switch (fSpecie){
+ case 0:
+ if( isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiMinus);
+ MakeTrack(invMassXiMinus, XiPt, /*xi->Phi(),*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }// endif
+
+ break;
+
+ case 1:
+ if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8
+ && (invMassXiMinus > 1.32486 || invMassXiMinus < 1.30486)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaMinus);
+
+ MakeTrack(invMassOmegaMinus, XiPt, /* xi->Phi(),*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }//endif
+ break;
+ }
+ }
+ }//end if ChargeXi()<0
+
+ if(xi->ChargeXi() > 0){
+ if(isNegProtonForTPC && isPosPionForTPC){
+ switch(fSpecie){
+ case 0:
+ if (isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiPlus);
+
+ //candidate inserting
+ MakeTrack(invMassXiPlus, XiPt, /* xi->Phi(),*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }//endif particle id
+ break;
+
+ case 1:
+ if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8
+ && (invMassXiPlus > 1.32486 || invMassXiPlus < 1.30486)){
+ ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaPlus);
+ MakeTrack(invMassOmegaPlus, XiPt, /* xi->Phi(),*/
+ phiAfter, xi->Eta(),
+ pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());
+ }//endif particle id
+ }
+ }
+ }//endif ChargeXi()>0
+ }//for Xi candidate loop
+
+ return;
+
+}
+
+
+void AliAnalysisTaskFlowCascade::MakeTrack( double mass,
+ double pt,
+ double phi,
+ double eta,
+ int iid,
+ int jid,
+ int kid) {
+ // create track for flow tasks
+ if(fCandidates->GetLast()+1>=fCandidates->GetSize()) {
+ fCandidates->Expand( 2*fCandidates->GetSize() );
+ }
+ Bool_t overwrite = kTRUE;
+
+ AliFlowCandidateTrack *sTrack
+ = (static_cast<AliFlowCandidateTrack*> (fCandidates->At( fCandidates->GetLast()+1 )));
+ if( !sTrack ) { // creates new
+ sTrack = new AliFlowCandidateTrack();
+ overwrite = kFALSE;
+ } else { // overwrites
+ sTrack->ClearMe();
+ }
+
+
+ sTrack->SetMass(mass);
+ sTrack->SetPt(pt);
+ sTrack->SetPhi(phi);
+ sTrack->SetEta(eta);
+ sTrack->AddDaughter(iid);
+ sTrack->AddDaughter(jid);
+ sTrack->AddDaughter(kid);
+ sTrack->SetForPOISelection(kTRUE);
+ sTrack->SetForRPSelection(kFALSE);
+
+ if(overwrite) {
+ fCandidates->SetLast( fCandidates->GetLast()+1 );
+ } else {
+ fCandidates->AddLast(sTrack);
+ }
+
+ return;
+}
+//_____________________________________________________________________________
+void AliAnalysisTaskFlowCascade::Terminate(Option_t *)
+{
+
+}
+
+void AliAnalysisTaskFlowCascade::Propagate(Double_t vv[3],
+ Double_t x[3],
+ Double_t p[3],
+ Double_t bz,
+ Short_t sign){
+ //Propagation to the primary vertex to determine the px and py
+ //x, p are the position and momentum as input and output
+ //bz is the magnetic field along z direction
+ //sign is the charge of particle for propagation
+
+ Double_t pp = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
+ Double_t len = (vv[2]-x[2])*pp/p[2];
+ Double_t a = -kB2C*bz*sign;
+
+ Double_t rho = a/pp;
+ x[0] += p[0]*TMath::Sin(rho*len)/a - p[1]*(1-TMath::Cos(rho*len))/a;
+ x[1] += p[1]*TMath::Sin(rho*len)/a + p[0]*(1-TMath::Cos(rho*len))/a;
+ x[2] += p[2]*len/pp;
+
+ Double_t p0=p[0];
+ p[0] = p0 *TMath::Cos(rho*len) - p[1]*TMath::Sin(rho*len);
+ p[1] = p[1]*TMath::Cos(rho*len) + p0 *TMath::Sin(rho*len);
+}
+
+
+//=====================================================================
+void AliAnalysisTaskFlowCascade::SetCommonConstants(Int_t massBins,
+ Double_t minMass,
+ Double_t maxMass)
+{
+ // setter for mass bins
+
+ fMassBins = massBins;
+ fMinMass = minMass;
+ fMaxMass = maxMass;
+}
+
+
+//====================================================================
+void AliAnalysisTaskFlowCascade::SetCuts2010(int set) {
+
+ // fCascadeCuts[0]: DcaXiDaughter; fCascadeCuts[1]: XiCosOfPointingAngle
+ // fCascadeCuts[2]: DcaV0ToPrimaryVtxXi; fCascadeCuts[3]: DcaBachToPrimaryVtxXi
+ // fCascadeCuts[4]: DcaV0DaughersXi; fCascadeCuts[5]: V0CosOfPointingAngleXi
+ // fCascadeCuts[6]: DcaV0DaughterToPrimaryVtxXi; fCascadeCuts[7]: V0MassWidth
+
+ switch(set){
+
+ case 0: //tighter
+ fCascadeCuts[0] = 0.2; fCascadeCuts[1] = 0.999;
+ fCascadeCuts[2] = 0.03; fCascadeCuts[3] = 0.05;
+ fCascadeCuts[4] = .5; fCascadeCuts[5] = 0.9998;
+ fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.006;
+ break;
+
+ case 1: //middle
+ fCascadeCuts[0] = 0.3; fCascadeCuts[1] = 0.99;
+ fCascadeCuts[2] = 0.01; fCascadeCuts[3] = 0.03;
+ fCascadeCuts[4] = .6; fCascadeCuts[5] = 0.9999;
+ fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.008;
+ break;
+
+ case 2: //looser
+ fCascadeCuts[0] = 0.3; fCascadeCuts[1] = 0.99;
+ fCascadeCuts[2] = 0.01; fCascadeCuts[3] = 0.03;
+ fCascadeCuts[4] = 1.; fCascadeCuts[5] = 1.;
+ fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.01;
+ break;
+ }
+
+}