PWG/FLOW/Tasks/AliAnalysisTaskFlowCascade.cxx

   1 /**************************************************************************\r
   2 * Copyright(c) 1998-2008,ALICE Experiment at CERN,All rights reserved.    *\r
   3 *                                                                         *\r
   4 * Author: The ALICE Off-line Project.                                     *\r
   5 * Contributors are mentioned in the code where appropriate.               *\r
   6 *                                                                         *\r
   7 * Permission to use,copy,modify and distribute this software and its      *\r
   8 * documentation strictly for non-commercial purposes is hereby granted    *\r
   9 * without fee,provided that the above copyright notice appears in all     *\r
  10 * copies and that both the copyright notice and this permission notice    *\r
  11 * appear in the supporting documentation. The authors make no claims      *\r
  12 * about the suitability of this software for any purpose. It is           *\r
  13 * provided "as is" without express or implied warranty.                   *\r
  14 **************************************************************************/\r
  15 \r
  16 /////////////////////////////////////////////////////\r
  17 // AliAnalysisTaskFlowCascade:\r
  18 // Analysis task to select K0/Lambda candidates for flow analysis.\r
  19 //\r
  20 // Author: \r
  21 //////////////////////////////////////////////////////\r
  22 \r
  23 #include "TChain.h"\r
  24 #include "TList.h"\r
  25 #include "TH1D.h"\r
  26 #include "TH2D.h"\r
  27 #include "TH3D.h"\r
  28 #include "TProfile.h"\r
  29 #include "TVector3.h"\r
  30 \r
  31 #include "AliAnalysisTaskSE.h"\r
  32 #include "AliAnalysisManager.h"\r
  33 \r
  34 #include "AliVParticle.h"\r
  35 #include "AliESDEvent.h"\r
  36 #include "AliESDInputHandler.h"\r
  37 #include "AliESDv0.h"\r
  38 #include "AliESDcascade.h"\r
  39 #include "AliESDtrack.h"\r
  40 #include "AliESDtrackCuts.h"\r
  41 #include "AliCentrality.h"\r
  42 #include "AliVVertex.h"\r
  43 #include "AliESDVZERO.h"\r
  44 #include "AliESDUtils.h"\r
  45 \r
  46 #include "AliTPCPIDResponse.h"\r
  47 #include "AliTOFPIDResponse.h"\r
  48 #include "AliPIDResponse.h"\r
  49 \r
  50 #include "AliAODEvent.h"\r
  51 #include "AliAODInputHandler.h"\r
  52 #include "AliAODTrack.h"\r
  53 #include "AliAODVZERO.h"\r
  54 #include "AliAODcascade.h"\r
  55 \r
  56 #include "TMath.h"\r
  57 #include "TObjArray.h"\r
  58 #include "AliFlowCandidateTrack.h"\r
  59 \r
  60 #include "AliFlowTrackCuts.h"\r
  61 #include "AliFlowEventCuts.h"\r
  62 #include "AliFlowEvent.h"\r
  63 #include "AliFlowCommonConstants.h"\r
  64 \r
  65 #include "AliAnalysisTaskFlowCascade.h"\r
  66 \r
  67 ClassImp(AliAnalysisTaskFlowCascade)\r
  68 \r
  69 //_____________________________________________________________________________\r
  70   AliAnalysisTaskFlowCascade::AliAnalysisTaskFlowCascade() :\r
  71     AliAnalysisTaskSE(),\r
  72     //    fMinCent(0), fMaxCent(0),\r
  73     fSpecie(0),\r
  74     fMassBins(0),\r
  75     fMinMass(0.0),\r
  76     fMaxMass(0.0),\r
  77     fCutsEvent(NULL),\r
  78     fCutsRPTPC(NULL),\r
  79     fCutsRPVZE(NULL),\r
  80     fCutsPOI(NULL),\r
  81     fCutsDau(NULL),\r
  82     fPIDResponse(NULL),\r
  83     fFlowEventTPC(NULL),\r
  84     fFlowEventVZE(NULL),\r
  85     fCandidates(NULL),\r
  86     fQAList(NULL)\r
  87 {\r
  88   //ctor                                                                       \r
  89   for (Int_t i=0; i!=8; ++i)\r
  90     fCascadeCuts[i] = 0;\r
  91   \r
  92 }\r
  93 \r
  94 //_____________________________________________________________________________\r
  95 AliAnalysisTaskFlowCascade\r
  96 ::AliAnalysisTaskFlowCascade(const char *name,\r
  97                              AliFlowEventCuts *cutsEvent, \r
  98                              AliFlowTrackCuts *cutsRPTPC,\r
  99                              AliFlowTrackCuts *cutsRPVZE,\r
 100                              /* AliESDtrackCuts */ AliFlowTrackCuts *cutsDau ) :\r
 101   AliAnalysisTaskSE(name),\r
 102   //fMinCent(minCent), fMaxCent(maxCent),\r
 103   fSpecie(0),\r
 104   fMassBins(0),\r
 105   fMinMass(0.0),\r
 106   fMaxMass(0.0),\r
 107   fCutsEvent(cutsEvent),\r
 108   fCutsRPTPC(cutsRPTPC),\r
 109   fCutsRPVZE(cutsRPVZE),\r
 110   fCutsPOI(NULL),\r
 111   fCutsDau(cutsDau),\r
 112   fPIDResponse(NULL),\r
 113   fFlowEventTPC(NULL),\r
 114   fFlowEventVZE(NULL),\r
 115   fCandidates(NULL),\r
 116   fQAList(NULL)\r
 117 {\r
 118   //ctor                                                                       \r
 119   for (Int_t i=0; i!=8; ++i)\r
 120     fCascadeCuts[i] = 0;\r
 121 \r
 122   DefineInput( 0,TChain::Class());\r
 123   DefineOutput(1,AliFlowEventSimple::Class()); // TPC object\r
 124   DefineOutput(2,AliFlowEventSimple::Class()); // VZE object\r
 125   DefineOutput(3,TList::Class());\r
 126 }\r
 127 \r
 128 //_____________________________________________________________________________\r
 129 AliAnalysisTaskFlowCascade::~AliAnalysisTaskFlowCascade()\r
 130 {\r
 131   if(fQAList) delete fQAList;\r
 132   if (fFlowEventTPC) delete fFlowEventTPC;\r
 133   if (fFlowEventVZE) delete fFlowEventVZE;\r
 134   if (fCandidates)   delete fCandidates;\r
 135   if (fCutsDau)      delete fCutsDau;\r
 136   if (fCutsPOI)      delete fCutsPOI;\r
 137   if (fCutsRPTPC)   delete fCutsRPTPC;\r
 138   if (fCutsRPVZE)   delete fCutsRPVZE;\r
 139 }\r
 140 \r
 141 //_____________________________________________________________________________\r
 142 void AliAnalysisTaskFlowCascade::UserCreateOutputObjects()\r
 143 {\r
 144 \r
 145   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();\r
 146   AliInputEventHandler* inputHandler\r
 147     = (AliInputEventHandler*) (man->GetInputEventHandler());\r
 148   fPIDResponse = inputHandler->GetPIDResponse();\r
 149   \r
 150   fQAList = new TList();\r
 151   fQAList->SetOwner();\r
 152   AddQAEvents();\r
 153   AddQACandidates();\r
 154   PostData(3,fQAList);\r
 155 \r
 156   AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();\r
 157   cc->SetNbinsMult(1);\r
 158   cc->SetMultMin(0);\r
 159   cc->SetMultMax(1);\r
 160 \r
 161   cc->SetNbinsPt(20);\r
 162   cc->SetPtMin(0.0);\r
 163   cc->SetPtMax(10.0);\r
 164 \r
 165   cc->SetNbinsPhi(1);\r
 166   cc->SetPhiMin(0.0);\r
 167   cc->SetPhiMax(TMath::TwoPi());\r
 168 \r
 169   cc->SetNbinsEta(1);\r
 170   cc->SetEtaMin(-2.0);\r
 171   cc->SetEtaMax(+2.0);\r
 172 \r
 173   cc->SetNbinsQ(3);\r
 174   cc->SetQMin(0.0);\r
 175   cc->SetQMax(3.0);\r
 176 \r
 177   cc->SetNbinsMass(fMassBins);\r
 178   cc->SetMassMin(fMinMass);\r
 179   cc->SetMassMax(fMaxMass);\r
 180 \r
 181   fCutsPOI = new AliFlowTrackCuts("null_cuts");\r
 182   fCutsPOI->SetParamType(fCutsRPTPC->GetParamType());\r
 183   fCutsPOI->SetPtRange(+1,-1); // select nothing QUICK   \r
 184   fCutsPOI->SetEtaRange(+1,-1); // select nothing VZERO  \r
 185 \r
 186   fFlowEventTPC = new AliFlowEvent(3000);\r
 187   fFlowEventVZE = new AliFlowEvent(1000);\r
 188   fCandidates = new TObjArray(100);\r
 189   fCandidates->SetOwner();\r
 190 \r
 191   PostData(1,fFlowEventTPC);\r
 192   PostData(2,fFlowEventVZE);\r
 193   PostData(3,fQAList);\r
 194 }\r
 195 \r
 196 //_____________________________________________________________________________\r
 197 void AliAnalysisTaskFlowCascade::AddQAEvents()\r
 198 {\r
 199   TList *tQAEvents = new TList();\r
 200   tQAEvents->SetName("Events");\r
 201   tQAEvents->SetOwner();\r
 202   TH1I* tEvent = new TH1I("Event","Number of Events",   3,0,3);\r
 203   tQAEvents->Add(tEvent);\r
 204   \r
 205   TH1D *tTPCRFP = new TH1D("RFPTPC",\r
 206                            "TPC Reference Flow Particles;multiplicity",\r
 207                            100, 0, 3000); \r
 208   tQAEvents->Add(tTPCRFP);\r
 209   TH1D *tVZERFP = new TH1D("RFPVZE", \r
 210                            "VZERO Reference Flow Particles;multiplicity",\r
 211                            100, 0, 30000); \r
 212   tQAEvents->Add(tVZERFP);\r
 213 \r
 214   TProfile *tCuts = new TProfile("Cuts","Analysis Cuts",10,0,10);\r
 215   tCuts->Fill(0.5,fCascadeCuts[0],1); \r
 216   tCuts->GetXaxis()->SetBinLabel(1,"dcaXiDau");\r
 217   tCuts->Fill(1.5,fCascadeCuts[1],1); \r
 218   tCuts->GetXaxis()->SetBinLabel(2,"XiCPA");\r
 219   tCuts->Fill(2.5,fCascadeCuts[2],1); \r
 220   tCuts->GetXaxis()->SetBinLabel(3,"dcaV0Vtx");\r
 221   tCuts->Fill(3.5,fCascadeCuts[3],1); \r
 222   tCuts->GetXaxis()->SetBinLabel(4,"dcaBachVtx");\r
 223   tCuts->Fill(4.5,fCascadeCuts[4],1); \r
 224   tCuts->GetXaxis()->SetBinLabel(5,"dcaV0Dau");\r
 225   tCuts->Fill(5.5,fCascadeCuts[5],1); \r
 226   tCuts->GetXaxis()->SetBinLabel(6,"V0CPA");\r
 227   tCuts->Fill(6.5,fCascadeCuts[6],1); \r
 228   tCuts->GetXaxis()->SetBinLabel(7,"dcaV0DauVtx");\r
 229   tCuts->Fill(7.5,fCascadeCuts[7],1); \r
 230   tCuts->GetXaxis()->SetBinLabel(8,"V0Mass");\r
 231   tQAEvents->Add(tCuts);\r
 232 \r
 233   fQAList->Add(tQAEvents);\r
 234 }\r
 235 \r
 236 //_____________________________________________________________________________\r
 237 void AliAnalysisTaskFlowCascade::AddQACandidates()\r
 238 {\r
 239   TList *tQACandidates;\r
 240   \r
 241   tQACandidates = new TList();\r
 242   tQACandidates->SetOwner();\r
 243   tQACandidates->SetName("Candidates");\r
 244 \r
 245   TH1F* tChi2Xi = new TH1F("Chi2Xi", \r
 246                        "Cascade #chi^{2}; #chi^{2}; Number of Cascades", \r
 247                        160, 0, 160);\r
 248   tQACandidates->Add(tChi2Xi);\r
 249 \r
 250   TH1F* tDCAXiDaughters \r
 251     = new TH1F( "DcaXiDaughters",  \r
 252                 "DCA between Xi Daughters; DCA (cm); Number of Cascades", \r
 253                 100, 0., 0.5);\r
 254   tQACandidates->Add(tDCAXiDaughters);\r
 255 \r
 256   TH1F * tDCABachToPrimVertex\r
 257     = new TH1F("DcaBachToPrimVertex", \r
 258                "DCA of Bach. to Prim. Vertex; DCA (cm);Number of Cascades", \r
 259                250, 0., 2.5);\r
 260   tQACandidates->Add(tDCABachToPrimVertex);\r
 261   \r
 262   TH1F * tXiCosOfPointingAngle\r
 263     = new TH1F("XiCosOfPointingAngle",\r
 264                "Cos of Xi Pointing Angle; Cos (Xi Point.Angl);Number of Xis", \r
 265                200, 0.99, 1.0);\r
 266   tQACandidates->Add(tXiCosOfPointingAngle);\r
 267 \r
 268   TH1F * tXiRadius = new TH1F("XiRadius",  \r
 269                               "Casc. decay transv. radius; r (cm); Counts" , \r
 270                               1050, 0., 105.0 );\r
 271   tQACandidates->Add(tXiRadius);\r
 272   \r
 273   TH1F *tMassLambda \r
 274     = new TH1F("MassLambdaAsCascDghter",\r
 275                "#Lambda assoc. to Casc. candidates; Eff. Mass (GeV/c^{2}); Counts", \r
 276                300,1.00,1.3);\r
 277   tQACandidates->Add(tMassLambda);\r
 278 \r
 279   TH1F *tV0Chi2 = new TH1F("V0Chi2Xi", \r
 280                        "V0 #chi^{2}, in cascade; #chi^{2};Counts", \r
 281                        160, 0, 40);\r
 282   tQACandidates->Add(tV0Chi2);\r
 283   \r
 284   TH1F * tV0CosOfPointingAngle \r
 285     = new TH1F("V0CosOfPointingAngleXi", \r
 286                "Cos of V0 Pointing Angle, in cascade;Cos(V0 Point. Angl); Counts", \r
 287                200, 0.98, 1.0);\r
 288   tQACandidates->Add(tV0CosOfPointingAngle);\r
 289 \r
 290   TH1F *tV0Radius  = new TH1F("V0RadiusXi", \r
 291                           "V0 decay radius, in cascade; radius (cm); Counts", \r
 292                           1050, 0., 105.0);\r
 293   tQACandidates->Add(tV0Radius);\r
 294   \r
 295   TH1F * tDcaV0DaughtersXi \r
 296     = new TH1F("DcaV0DaughtersXi", \r
 297                "DCA between V0 daughters, in cascade;DCA (cm);Number of V0s", \r
 298                120, 0., 0.6);\r
 299   tQACandidates->Add(tDcaV0DaughtersXi);\r
 300 \r
 301   TH1F * tDcaV0ToPrimVertex \r
 302     = new TH1F("DcaV0ToPrimVertexXi", \r
 303                "DCA of V0 to Prim. Vertex, in cascade;DCA (cm);Number of Cascades", 200, 0., 1.);\r
 304   tQACandidates->Add(tDcaV0ToPrimVertex);\r
 305 \r
 306   TH1F * tDCAPosToPrimVertex =\r
 307     new TH1F("DcaPosToPrimVertexXi", \r
 308              "DCA of V0 pos daughter to Prim. Vertex;DCA (cm);Counts", \r
 309              300, 0, 3);\r
 310   tQACandidates->Add(tDCAPosToPrimVertex);\r
 311 \r
 312   TH1F * tDCANegToPrimVertex \r
 313     =  new TH1F("DcaNegToPrimVertexXi", \r
 314                 "DCA of V0 neg daughter to Prim. Vertex;DCA (cm);Counts", \r
 315                 300, 0, 3);\r
 316   tQACandidates->Add(tDCANegToPrimVertex);\r
 317 \r
 318   TH1F *tV0toXiCosOfPointingAngle \r
 319     = new TH1F("V0toXiCosOfPointingAngle", \r
 320                "Cos. of V0 Ptng Angl Xi vtx; Cos(V0 Point. Angl / Xi vtx); Counts", \r
 321                100, 0.99, 1.0);\r
 322   tQACandidates->Add(tV0toXiCosOfPointingAngle);\r
 323 \r
 324   TH2F *th2Armenteros \r
 325     = new TH2F("Armenteros", \r
 326                "#alpha_{Arm}(casc. cand.) Vs Pt_{Arm}(casc. cand.); #alpha_{Arm} ; Pt_{Arm} (GeV/c)", \r
 327                140, -1.2, 1.2, 300, 0., 0.3);\r
 328   tQACandidates->Add(th2Armenteros);\r
 329 \r
 330   TH2F *th2TPCdEdxOfCascDghters\r
 331     = new TH2F( "TPCdEdxOfCascDghters",\r
 332                 "TPC dE/dx of the cascade daughters; charge x || #vec{p}_{TPC inner wall}(Casc. daughter) || (GeV/c); TPC signal (ADC) ", \r
 333                 200, -10.0, 10.0, 450, 0., 900.);\r
 334   tQACandidates->Add(th2TPCdEdxOfCascDghters);\r
 335 \r
 336   TH2F *th2MassVsPtAll \r
 337     = new TH2F("MassVsPtAll", \r
 338                "M_{candidates} vs Pt; Pt (GeV/c); M (GeV/c^{2})", \r
 339                100, 0., 10., fMassBins, fMinMass, fMaxMass);\r
 340   tQACandidates->Add(th2MassVsPtAll);\r
 341 \r
 342   TH1F *tDistToVtxZAfter \r
 343     = new TH1F("DistToVtxZAfter", \r
 344                "Distance to vtx z after propagation to vtx; z [cm]", \r
 345                100, -5., 5.);\r
 346   tQACandidates->Add(tDistToVtxZAfter);\r
 347 \r
 348   TH1F * tDistToVtxXYAfter \r
 349     = new TH1F("DistToVtxXYAfter", \r
 350                "Distance to vtx xy after propagation to vtx",\r
 351                500, 0., 50.);\r
 352   tQACandidates->Add(tDistToVtxXYAfter);\r
 353 \r
 354   TH2F *th2DistToVtxZBeforeVsAfter \r
 355     = new TH2F("DistToVtxZBeforeVsAfter", \r
 356                "Distance to vtx z before vs after propagation; Distance before [cm]; Distance after [cm]", \r
 357                500, -50., 50., 100, -5., 5.);\r
 358   tQACandidates->Add(th2DistToVtxZBeforeVsAfter);\r
 359 \r
 360   TH2F *th2DistToVtxXYBeforeVsAfter\r
 361     = new TH2F("DistToVtxXYBeforeVsAfter",\r
 362                "Distance to vtx xy before vs after propagation; Distance before [cm]; Distance after [cm]", \r
 363                500, 0., 50, 500, 0., 50);\r
 364   tQACandidates->Add(th2DistToVtxXYBeforeVsAfter);\r
 365 \r
 366   TH2F * th2PxBeforeVsAfter \r
 367     = new TH2F("PxBeforeVsAfter", \r
 368                "Px before vs after propagation; Px [GeV/c]; Px' [GeV/c]", \r
 369                200, -10., 10, 200, -10., 10.);\r
 370   tQACandidates->Add(th2PxBeforeVsAfter);\r
 371 \r
 372   TH2F * th2PyBeforeVsAfter\r
 373     = new TH2F("PyBeforeVsAfter",\r
 374                "Py before vs after propagation; Py [GeV/c]; Py' [GeV/c]",\r
 375                200, -10., 10, 200, -10., 10.);\r
 376   tQACandidates->Add(th2PyBeforeVsAfter);\r
 377   \r
 378   TH2F * th2PhiPosBeforeVsAfter\r
 379     = new TH2F("PhiPosBeforeVsAfter", \r
 380                "Phi for positively charged candidates before vs after propagation; #phi; #phi'", \r
 381                360, 0., 2.0*TMath::Pi(), 360, 0., 2.0*TMath::Pi());\r
 382   tQACandidates->Add(th2PhiPosBeforeVsAfter);\r
 383 \r
 384   TH2F *th2PhiNegBeforeVsAfter\r
 385     = new TH2F("PhiNegBeforeVsAfter",\r
 386                "Phi for negatively charged candidates before vs after propagation; #phi; #phi'",\r
 387                360, 0., 2.0*TMath::Pi(), 360, 0., 2.0*TMath::Pi());\r
 388   tQACandidates->Add(th2PhiNegBeforeVsAfter);\r
 389 \r
 390   fQAList->Add(tQACandidates);\r
 391 }\r
 392 \r
 393 //_____________________________________________________________________________\r
 394 void AliAnalysisTaskFlowCascade::NotifyRun()\r
 395 {\r
 396 }\r
 397 //_____________________________________________________________________________\r
 398 void AliAnalysisTaskFlowCascade::UserExec(Option_t *)\r
 399 {\r
 400   AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());\r
 401   AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());\r
 402   Bool_t acceptEvent=kFALSE; \r
 403   fCandidates->SetLast(-1);\r
 404 \r
 405   if(fESD) {\r
 406     // recorrecting VZERO (for pass 1 only)\r
 407     /*\r
 408     Float_t *vChCorr = new Float_t[64];\r
 409     Float_t dummy;\r
 410     AliESDUtils::GetCorrV0(fESD,dummy,NULL,vChCorr);\r
 411     AliESDVZERO *vzero = (AliESDVZERO*) fESD->GetVZEROData();\r
 412     vzero->SetMultiplicity( vChCorr );\r
 413     delete [] vChCorr;\r
 414     */\r
 415     //\r
 416 \r
 417     ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(0);\r
 418     \r
 419     const AliVVertex *vtxGlb = fESD->GetPrimaryVertexTracks();\r
 420     const AliVVertex *vtxSPD = fESD->GetPrimaryVertexSPD();\r
 421     if(!vtxGlb || !vtxSPD) return;\r
 422     if( fCutsEvent->IsSelected(fESD, 0) && (TMath::Abs(vtxSPD->GetZ()-vtxGlb->GetZ()) <= 0.5) ) {\r
 423       acceptEvent = kTRUE;\r
 424     ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);\r
 425       ReadFromESDv0(fESD);\r
 426     }\r
 427   } else if(fAOD) {\r
 428     const AliVVertex *vtxGlb = fAOD->GetPrimaryVertex();\r
 429     const AliVVertex *vtxSPD = fAOD->GetPrimaryVertexSPD();\r
 430     if(!vtxGlb || !vtxSPD) return;\r
 431 \r
 432     ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(0);\r
 433         \r
 434     if(fCutsEvent->IsSelected(fAOD, 0) && (TMath::Abs(vtxSPD->GetZ()-vtxGlb->GetZ()) <= 0.5) ) {\r
 435       acceptEvent = kTRUE;\r
 436       ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);\r
 437       ReadFromAODv0(fAOD);\r
 438     }\r
 439 \r
 440     \r
 441     /*\r
 442 \r
 443     AliAODHeader *aodHeader = fAOD->GetHeader();\r
 444     if(!aodHeader) return;\r
 445     AliCentrality *centrality = aodHeader->GetCentralityP();\r
 446     if(!centrality) return;\r
 447     Double_t cent = centrality->GetCentralityPercentile("V0M" );\r
 448     Double_t cent1 = centrality->GetCentralityPercentile("TRK" );\r
 449     if(TMath::Abs(cent-cent1) >= 5.) return;\r
 450     \r
 451     if(cent<fMinCent||cent>=fMaxCent) return; //centrality cut\r
 452     \r
 453     Double_t zvtx = fAOD->GetPrimaryVertex()->GetZ();\r
 454     if(TMath::Abs(zvtx)>10.) return; //vertex cut\r
 455     \r
 456     ((TH1I*)((TList*)fQAList->FindObject("Events"))->FindObject("Event"))->Fill(2);\r
 457     ReadFromAODv0(fAOD);\r
 458     */  \r
 459   }\r
 460   \r
 461   if(!acceptEvent) return;\r
 462 \r
 463   ((TH1D*)((TList*)fQAList->FindObject("Events"))->FindObject("RFPTPC"))\r
 464     ->Fill(fFlowEventTPC->GetNumberOfRPs() );\r
 465   Double_t mult=0;\r
 466   for(Int_t i=0; i != fFlowEventVZE->GetNumberOfRPs(); ++i) {\r
 467     AliFlowTrackSimple *pTrack = fFlowEventVZE->GetTrack(i);\r
 468     mult += pTrack->Weight();\r
 469   }\r
 470   ((TH1D*)((TList*)fQAList->FindObject("Events"))->FindObject("RFPVZE"))\r
 471     ->Fill( mult );\r
 472 \r
 473   if(fDebug) printf("TPCevent %d | VZEevent %d\n",\r
 474                     fFlowEventTPC->NumberOfTracks(),\r
 475                     fFlowEventVZE->NumberOfTracks() );\r
 476   AddCandidates();\r
 477 \r
 478   PostData(1,fFlowEventTPC);\r
 479   PostData(2,fFlowEventVZE);\r
 480   PostData(3,fQAList);\r
 481 \r
 482   return;\r
 483 }\r
 484 \r
 485 //_____________________________________________________________________________\r
 486 void AliAnalysisTaskFlowCascade::AddCandidates(){\r
 487   \r
 488   if(fDebug) printf("I received %d candidates\n",\r
 489                     fCandidates->GetEntriesFast());\r
 490   for(int iCand=0; iCand!=fCandidates->GetEntriesFast(); ++iCand ) {\r
 491     AliFlowCandidateTrack *cand \r
 492       = dynamic_cast<AliFlowCandidateTrack*>(fCandidates->At(iCand));\r
 493     if(!cand) continue;\r
 494     if(fDebug) \r
 495       printf(" >Checking at candidate %d with %d daughters: mass %f\n",\r
 496              iCand, cand->GetNDaughters(), cand->Mass());\r
 497     // untagging ===>                      \r
 498     for(int iDau=0; iDau != cand->GetNDaughters(); ++iDau) {\r
 499       if(fDebug) \r
 500         printf(" >Daughter %d with fID %d", iDau, cand->GetIDDaughter(iDau));\r
 501       for(int iRPs=0; iRPs != fFlowEventTPC->NumberOfTracks(); ++iRPs ) {\r
 502         AliFlowTrack *iRP \r
 503           = dynamic_cast<AliFlowTrack*>(fFlowEventTPC->GetTrack( iRPs ));\r
 504         if (!iRP) continue;\r
 505         if( !iRP->InRPSelection() ) continue;\r
 506         if( cand->GetIDDaughter(iDau) == iRP->GetID() ) {\r
 507           if(fDebug) printf(" was in RP set");\r
 508           iRP->SetForRPSelection(kFALSE);\r
 509           fFlowEventTPC->SetNumberOfRPs( fFlowEventTPC->GetNumberOfRPs() -1 );\r
 510         }\r
 511       }\r
 512       if(fDebug) printf("\n");\r
 513     }\r
 514     // <=== untagging\r
 515     cand->SetForPOISelection(kTRUE);\r
 516     fFlowEventTPC->InsertTrack( ((AliFlowTrack*) cand) );\r
 517     fFlowEventVZE->InsertTrack( ((AliFlowTrack*) cand) );\r
 518   }\r
 519 \r
 520   if(fDebug) printf("TPCevent %d | VZEevent %d\n",\r
 521                     fFlowEventTPC->NumberOfTracks(),\r
 522                     fFlowEventVZE->NumberOfTracks() );\r
 523 \r
 524 }\r
 525 \r
 526 //______________________________________________________________________________\r
 527 void AliAnalysisTaskFlowCascade::ReadFromESDv0(AliESDEvent *fESD)\r
 528 {\r
 529   //AliFlowTrackCuts* cutsPOI = new AliFlowTrackCuts("null_cuts");\r
 530   //cutsPOI->SetParamType( fCutsRP->GetParamType() );\r
 531   //cutsPOI->SetParamType( AliFlowTrackCuts::kGlobal );\r
 532   // cutsPOI->SetPtRange(+1,-1); // select nothing\r
 533   //cutsPOI->SetEtaRange(+1,-1); // select nothing VZERO\r
 534 \r
 535   fCutsRPTPC->SetEvent(fESD,MCEvent());\r
 536   fCutsRPVZE->SetEvent(fESD,MCEvent());\r
 537 \r
 538   fCutsPOI->SetEvent(fESD,MCEvent());\r
 539 \r
 540   fFlowEventTPC->Fill(fCutsRPTPC,fCutsPOI);\r
 541   fFlowEventVZE->Fill(fCutsRPVZE,fCutsPOI);\r
 542   \r
 543   Double_t trkPrimaryVtxPos[3] = {-100., -100., -100.};\r
 544   Double_t bestPrimaryVtxPos[3] = {-100., -100., -100.};\r
 545   int nCascades=fESD->GetNumberOfCascades();\r
 546   \r
 547   const AliESDVertex *primaryTrackingESDVtx = fESD->GetPrimaryVertexTracks();\r
 548   primaryTrackingESDVtx->GetXYZ(trkPrimaryVtxPos);\r
 549   \r
 550   const AliESDVertex *primaryBestESDVtx = fESD->GetPrimaryVertex();\r
 551   primaryBestESDVtx->GetXYZ(bestPrimaryVtxPos);\r
 552 \r
 553   Double_t b = fESD->GetMagneticField();\r
 554 \r
 555   for(int i = 0; i != nCascades; ++i) {\r
 556     \r
 557     // Double_t trkPrimaryVtxRadius3D = -500.;\r
 558     // Double_t bestPrimaryVtxRadius3D = -500.;\r
 559     Double_t effMassXi = 0.;\r
 560     Double_t chi2Xi = -1.;\r
 561     Double_t dcaXiDaughters = -1.;\r
 562     Double_t XiCosOfPointingAngle = -1.;\r
 563     Double_t posXi[3] = {-1000., -1000., -1000.};\r
 564     Double_t XiRadius = -1000.;\r
 565 \r
 566     // Double_t innerWallMomCascDghters[3] = {-100., -100., -100.};\r
 567     //Double_t tpcSignalCascDghters[3] = {-100., -100., -100.};\r
 568 \r
 569     Double_t invMassLambdaAsCascDghter = 0.;\r
 570     Double_t V0Chi2Xi = -1.;\r
 571     Double_t dcaV0DaughtersXi = -1.;\r
 572     \r
 573     Double_t dcaBachToPrimaryVtxXi = -1.;\r
 574     Double_t dcaV0ToPrimaryVtxXi = -1.;\r
 575     Double_t dcaPosToPrimaryVtxXi = -1.;\r
 576     Double_t dcaNegToPrimaryVtxXi = -1.;\r
 577     Double_t V0CosOfPointingAngleXi = -1.;\r
 578     Double_t posV0Xi[3] = {-1000., -1000., -1000.};\r
 579     Double_t V0RadiusXi = -1000.;\r
 580     Double_t V0quality = 0.;\r
 581 \r
 582     Double_t invMassXiMinus = 0.;\r
 583     Double_t invMassXiPlus = 0.;\r
 584     Double_t invMassOmegaMinus = 0.;\r
 585     Double_t invMassOmegaPlus = 0.;\r
 586 \r
 587     /*\r
 588     Bool_t isPosInXiProton = kFALSE;\r
 589     Bool_t isPosInXiPion = kFALSE;\r
 590     Bool_t isPosInOmegaProton = kFALSE;\r
 591     Bool_t isPosInOmegaPion = kFALSE;\r
 592     \r
 593     Bool_t isNegInXiProton = kFALSE;\r
 594     Bool_t isNegInXiPion = kFALSE;\r
 595     Bool_t isNegInOmegaProton = kFALSE;\r
 596     Bool_t isNegInOmegaPion = kFALSE;\r
 597 \r
 598     Bool_t isBachelorKaon = kFALSE;\r
 599     Bool_t isBachelorPion = kFALSE;\r
 600     */\r
 601 \r
 602     Bool_t isBachelorKaonForTPC = kFALSE;\r
 603     Bool_t isBachelorPionForTPC = kFALSE;\r
 604     Bool_t isNegPionForTPC = kFALSE;\r
 605     Bool_t isPosPionForTPC = kFALSE;\r
 606     Bool_t isNegProtonForTPC = kFALSE;\r
 607     Bool_t isPosProtonForTPC = kFALSE;\r
 608 \r
 609     Double_t XiPx = 0., XiPy = 0., XiPz = 0.;\r
 610     Double_t XiPt = 0.;\r
 611     Double_t XiPtot = 0.;\r
 612     \r
 613     Double_t bachPx = 0., bachPy = 0., bachPz = 0.;\r
 614     Double_t bachPt = 0.;\r
 615     Double_t bachPtot = 0.;\r
 616     \r
 617     //Short_t chargeXi = -2;\r
 618     Double_t V0toXiCosOfPointingAngle = 0.;\r
 619     \r
 620     Double_t rapXi = -20.;\r
 621     Double_t rapOmega = -20.;\r
 622     Double_t phi = 6.3;\r
 623     Double_t alphaXi = -200.;\r
 624     Double_t ptArmXi = -200.;\r
 625     //    TLorentzVector lv1, lv2, lv3, lv12, lvXi;\r
 626 \r
 627     Double_t distToVtxZBefore = -999.;\r
 628     Double_t distToVtxZAfter = -999.;\r
 629     Double_t distToVtxXYBefore = -999.;\r
 630     Double_t distToVtxXYAfter = -999.;\r
 631     Double_t XiPAfter[3] = {-999., -999., -999.};\r
 632     Double_t phiAfter = -999.;\r
 633     \r
 634     AliESDcascade *xi = fESD->GetCascade(i);\r
 635     if(!xi) continue;\r
 636     \r
 637     if(xi->Charge()<0)\r
 638       xi->ChangeMassHypothesis(V0quality, 3312); // Xi- hypothesis\r
 639     else if(xi->Charge() > 0)\r
 640       xi->ChangeMassHypothesis(V0quality, -3312);\r
 641     else continue;\r
 642 \r
 643     effMassXi = xi->GetEffMassXi();\r
 644     chi2Xi = xi->GetChi2Xi();\r
 645     dcaXiDaughters = xi->GetDcaXiDaughters();\r
 646     XiCosOfPointingAngle \r
 647       = xi->GetCascadeCosineOfPointingAngle(bestPrimaryVtxPos[0],\r
 648                                             bestPrimaryVtxPos[1],\r
 649                                             bestPrimaryVtxPos[2]);\r
 650     xi->GetXYZcascade(posXi[0], posXi[1], posXi[2]);\r
 651     XiRadius = TMath::Sqrt(posXi[0]*posXi[0]\r
 652                            +posXi[1]*posXi[1]\r
 653                            +posXi[2]*posXi[2]);\r
 654     \r
 655     UInt_t idxPosXi = (UInt_t)TMath::Abs(xi->GetPindex());\r
 656     UInt_t idxNegXi = (UInt_t)TMath::Abs(xi->GetNindex());\r
 657     UInt_t idxBach = (UInt_t)TMath::Abs(xi->GetBindex());\r
 658 \r
 659     if(idxBach == idxPosXi || idxBach == idxNegXi) continue;\r
 660 \r
 661     AliESDtrack *pTrkXi = fESD->GetTrack(idxPosXi);\r
 662     AliESDtrack *nTrkXi = fESD->GetTrack(idxNegXi);\r
 663     AliESDtrack *bTrkXi = fESD->GetTrack(idxBach);\r
 664 \r
 665     if( !pTrkXi || !nTrkXi || !bTrkXi ) continue;\r
 666 \r
 667     if( !fCutsDau->IsSelected(pTrkXi) \r
 668         || !fCutsDau->IsSelected(nTrkXi)\r
 669         || !fCutsDau->IsSelected(bTrkXi) ) continue;\r
 670 \r
 671     const AliExternalTrackParam *pExtTrk = pTrkXi->GetInnerParam();\r
 672     const AliExternalTrackParam *nExtTrk = nTrkXi->GetInnerParam();\r
 673     const AliExternalTrackParam *bExtTrk = bTrkXi->GetInnerParam();\r
 674     \r
 675     if(pExtTrk && pTrkXi->IsOn(AliESDtrack::kTPCin)){\r
 676       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(pExtTrk->GetP()*pExtTrk->Charge(), pTrkXi->GetTPCsignal());\r
 677     }\r
 678     if(nExtTrk && nTrkXi->IsOn(AliESDtrack::kTPCin)){\r
 679       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(nExtTrk->GetP()*nExtTrk->Charge(), nTrkXi->GetTPCsignal());\r
 680     }\r
 681     if(bExtTrk && bTrkXi->IsOn(AliESDtrack::kTPCin)){\r
 682       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(bExtTrk->GetP()*bExtTrk->Charge(), bTrkXi->GetTPCsignal());\r
 683     }\r
 684     \r
 685     invMassLambdaAsCascDghter = xi->GetEffMass(); // from V0\r
 686     dcaV0DaughtersXi = xi->GetDcaV0Daughters();\r
 687     V0Chi2Xi = xi->GetChi2V0();\r
 688     V0CosOfPointingAngleXi \r
 689       = xi->GetV0CosineOfPointingAngle(bestPrimaryVtxPos[0],\r
 690                                        bestPrimaryVtxPos[1],\r
 691                                        bestPrimaryVtxPos[2]);\r
 692     dcaV0ToPrimaryVtxXi = xi->GetD(bestPrimaryVtxPos[0], \r
 693                                    bestPrimaryVtxPos[1],\r
 694                                    bestPrimaryVtxPos[2]);\r
 695     dcaBachToPrimaryVtxXi = TMath::Abs(bTrkXi->GetD(bestPrimaryVtxPos[0],\r
 696                                                     bestPrimaryVtxPos[1],\r
 697                                                     bestPrimaryVtxPos[2]));\r
 698     \r
 699     //V0\r
 700     xi->GetXYZ(posV0Xi[0], posV0Xi[1], posV0Xi[2]);\r
 701     V0RadiusXi = TMath::Sqrt(posV0Xi[0]*posV0Xi[0]\r
 702                              +posV0Xi[1]*posV0Xi[1]\r
 703                              +posV0Xi[2]*posV0Xi[2]);\r
 704     dcaPosToPrimaryVtxXi = TMath::Abs(pTrkXi->GetD(bestPrimaryVtxPos[0],\r
 705                                                    bestPrimaryVtxPos[1],\r
 706                                                    bestPrimaryVtxPos[2]));\r
 707     dcaNegToPrimaryVtxXi = TMath::Abs(nTrkXi->GetD(bestPrimaryVtxPos[0],\r
 708                                                    bestPrimaryVtxPos[1],\r
 709                                                    bestPrimaryVtxPos[2]));\r
 710 \r
 711     //apply cuts\r
 712     //if(XiRadius < 0.9 || XiRadius > 100.) continue;\r
 713     //if(dcaXiDaughters > 0.2) continue;\r
 714     //if(XiCosOfPointingAngle < 0.99) continue;\r
 715     //if(dcaV0ToPrimaryVtxXi < 0.03) continue;\r
 716     //if(dcaBachToPrimaryVtxXi < 0.01) continue;\r
 717     \r
 718     if(dcaXiDaughters > fCascadeCuts[0]) continue;\r
 719     if(XiCosOfPointingAngle < fCascadeCuts[1]) continue;\r
 720     if(dcaV0ToPrimaryVtxXi < fCascadeCuts[2]) continue;\r
 721     if(dcaBachToPrimaryVtxXi < fCascadeCuts[3]) continue;\r
 722 \r
 723     //V0 mass cut?\r
 724     //  if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > 0.01) continue;\r
 725     if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > fCascadeCuts[7]) \r
 726       continue;\r
 727     \r
 728     //if(dcaV0DaughtersXi > 1.) continue;\r
 729     //if(V0CosOfPointingAngleXi > 0.9999) continue;\r
 730     //if(dcaPosToPrimaryVtxXi < 0.1) continue;\r
 731     //if(dcaNegToPrimaryVtxXi < 0.1) continue;\r
 732 \r
 733     if(dcaV0DaughtersXi > fCascadeCuts[4]) continue;\r
 734     if(V0CosOfPointingAngleXi > fCascadeCuts[5]) continue;\r
 735     if(dcaPosToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
 736     if(dcaNegToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
 737 \r
 738     //if(V0RadiusXi < 1.0 || V0RadiusXi > 100) continue;\r
 739     \r
 740     //other cuts?\r
 741     // change mass hypothesis to cover all the possibilities\r
 742     if(bTrkXi->Charge()<0){\r
 743       V0quality = 0.;\r
 744       xi->ChangeMassHypothesis(V0quality, 3312); //Xi- hyp.\r
 745       invMassXiMinus = xi->GetEffMassXi();\r
 746 \r
 747       V0quality = 0.;\r
 748       xi->ChangeMassHypothesis(V0quality, 3334); //Omega- hyp.\r
 749       invMassOmegaMinus = xi->GetEffMassXi();\r
 750 \r
 751       V0quality = 0.;\r
 752       xi->ChangeMassHypothesis(V0quality, 3312); //back to default hyp.\r
 753     }\r
 754 \r
 755     if(bTrkXi->Charge() > 0){\r
 756       V0quality = 0.;\r
 757       xi->ChangeMassHypothesis(V0quality, -3312); //anti-Xi- hyp.\r
 758       invMassXiPlus = xi->GetEffMassXi();\r
 759 \r
 760       V0quality = 0.;\r
 761       xi->ChangeMassHypothesis(V0quality, -3334); //anti-Omega- hyp.\r
 762       invMassOmegaPlus = xi->GetEffMassXi();\r
 763 \r
 764       V0quality = 0.;\r
 765       xi->ChangeMassHypothesis(V0quality, -3312); //back to default hyp.\r
 766     }\r
 767 \r
 768     //PID on the daughter tracks\r
 769     /*\r
 770     //A - Combined PID\r
 771     //Resonable guess the priors for the cascade track sample\r
 772     //(e, mu, pi, K, p)\r
 773     Double_t priorsGuessXi[5] = {0, 0, 2, 0, 1};\r
 774     Double_t priorsGuessOmega[5] = {0, 0, 1, 1, 1};\r
 775 \r
 776     //Combined bachelor-daughter PID\r
 777     AliPID pidXi;\r
 778     pidXi.SetPriors(priorsGuessXi);\r
 779     AliPID pidOmega;\r
 780     pidOmega.SetPriors(priorsGuessOmega);\r
 781 \r
 782     if(pTrkXi->IsOn(AliESDtrack::kESDpid)){// combined PID exists\r
 783       Double_t r[10] = {0.};\r
 784       pTrkXi->GetESDpid(r);\r
 785       pidXi.SetProbabilities(r);\r
 786       pidOmega.SetProbabilities(r);\r
 787 \r
 788       //Check if the V0 postive track is proton (case for Xi-)\r
 789       Double_t pProton = pidXi.GetProbability(AliPID::kProton);\r
 790       if(pProton > pidXi.GetProbability(AliPID::kElectron)\r
 791          && pProton > pidXi.GetProbability(AliPID::kMuon)\r
 792          && pProton > pidXi.GetProbability(AliPID::kPion)\r
 793          && pProton > pidXi.GetProbability(AliPID::kKaon))\r
 794         isPosInXiProton = kTRUE;\r
 795         \r
 796       //Check if the V0 postive track is a pi+ (case for Xi+)\r
 797       Double_t pPion = pidXi.GetProbability(AliPID::kPion);\r
 798       if(pPion > pidXi.GetProbability(AliPID::kElectron)\r
 799          && pPion > pidXi.GetProbability(AliPID::kMuon)\r
 800          && pPion > pidXi.GetProbability(AliPID::kKaon)\r
 801          && pPion > pidXi.GetProbability(AliPID::kProton))\r
 802         isPosInXiPion = kTRUE;\r
 803       // Check if the V0 positive track is a proton (case for Omega-)\r
 804       pProton = pidOmega.GetProbability(AliPID::kProton);\r
 805       if(pProton > pidOmega.GetProbability(AliPID::kElectron)\r
 806          && pProton > pidOmega.GetProbability(AliPID::kMuon)\r
 807          && pProton > pidOmega.GetProbability(AliPID::kPion)\r
 808          && pProton > pidOmega.GetProbability(AliPID::kKaon))\r
 809         isPosInOmegaProton = kTRUE;\r
 810     \r
 811       // Check if the V0 positive track is a pi+ (case for Omega+)\r
 812       pPion =  pidOmega.GetProbability(AliPID::kPion);\r
 813       if(pPion > pidOmega.GetProbability(AliPID::kElectron)\r
 814          && pPion > pidOmega.GetProbability(AliPID::kMuon)\r
 815          && pPion > pidOmega.GetProbability(AliPID::kKaon)\r
 816          && pPion > pidOmega.GetProbability(AliPID::kProton))\r
 817         isPosInOmegaPion = kTRUE;\r
 818     }\r
 819 \r
 820     //Combined V0-negative-daughter PID\r
 821     pidXi.SetPriors(priorsGuessXi);\r
 822     pidOmega.SetPriors(priorsGuessOmega);\r
 823     if(nTrkXi->IsOn(AliESDtrack::kESDpid)){\r
 824       Double_t r[10] = {0.};\r
 825       nTrkXi->GetESDpid(r);\r
 826       pidXi.SetProbabilities(r);\r
 827       pidOmega.SetProbabilities(r);\r
 828       \r
 829       // Check if the V0 negative track is a pi- (case for Xi-)\r
 830       Double_t pPion = pidXi.GetProbability(AliPID::kPion);\r
 831       if(pPion > pidXi.GetProbability(AliPID::kElectron)\r
 832          && pPion > pidXi.GetProbability(AliPID::kMuon)\r
 833          && pPion >  pidXi.GetProbability(AliPID::kKaon)\r
 834          && pPion > pidXi.GetProbability(AliPID::kProton))\r
 835         isNegInXiPion = kTRUE;\r
 836 \r
 837       // Check if the V0 negative track is an anti-proton (case for Xi+)\r
 838       Double_t pProton = pidXi.GetProbability(AliPID::kProton);\r
 839       if(pProton > pidXi.GetProbability(AliPID::kElectron)\r
 840          && pProton > pidXi.GetProbability(AliPID::kMuon)\r
 841          && pProton > pidXi.GetProbability(AliPID::kPion) \r
 842          && pProton > pidXi.GetProbability(AliPID::kKaon))\r
 843         isNegInXiProton = kTRUE;\r
 844       \r
 845       // Check if the V0 negative track is a pi- (case for Omega-)\r
 846       pPion = pidOmega.GetProbability(AliPID::kPion);\r
 847       if(pPion > pidOmega.GetProbability(AliPID::kElectron)\r
 848          && pPion > pidOmega.GetProbability(AliPID::kMuon)\r
 849          && pPion > pidOmega.GetProbability(AliPID::kKaon)\r
 850          && pPion > pidOmega.GetProbability(AliPID::kProton))\r
 851         isNegInOmegaPion = kTRUE;\r
 852       \r
 853       // Check if the V0 negative track is an anti-proton (case for Omega+)   \r
 854       pProton =  pidOmega.GetProbability(AliPID::kProton);\r
 855       if(pProton > pidOmega.GetProbability(AliPID::kElectron)\r
 856          && pProton > pidOmega.GetProbability(AliPID::kMuon)\r
 857          && pProton > pidOmega.GetProbability(AliPID::kPion)\r
 858          && pProton > pidOmega.GetProbability(AliPID::kKaon))\r
 859         isNegInOmegaProton = kTRUE;\r
 860       \r
 861     }\r
 862 \r
 863     // Combined bachelor PID\r
 864     pidXi.SetPriors(priorsGuessXi);\r
 865     pidOmega.SetPriors(priorsGuessOmega);\r
 866     if(bTrkXi->IsOn(AliESDtrack::kESDpid)){//Combined PID exists\r
 867       Double_t r[10] = {0.};\r
 868       bTrkXi->GetESDpid(r);\r
 869       pidXi.SetProbabilities(r);\r
 870       pidOmega.SetProbabilities(r);\r
 871       \r
 872       //Check if the bachelor track is a pion\r
 873       Double_t pPion = pidXi.GetProbability(AliPID::kPion);\r
 874       if(pPion >  pidXi.GetProbability(AliPID::kElectron)\r
 875          && pPion >  pidXi.GetProbability(AliPID::kMuon)\r
 876          && pPion > pidXi.GetProbability(AliPID::kKaon)\r
 877          && pPion > pidXi.GetProbability(AliPID::kProton))\r
 878         isBachelorPion = kTRUE;\r
 879       \r
 880       // Check if the bachelor track is a kaon\r
 881       Double_t pKaon = pidOmega.GetProbability(AliPID::kKaon);\r
 882       if(pKaon > pidOmega.GetProbability(AliPID::kElectron)\r
 883          && pKaon > pidOmega.GetProbability(AliPID::kMuon)\r
 884          && pKaon > pidOmega.GetProbability(AliPID::kPion)\r
 885          && pKaon > pidOmega.GetProbability(AliPID::kProton))\r
 886         isBachelorKaon = kTRUE;\r
 887     }\r
 888     */\r
 889 \r
 890 \r
 891     //B - TPC PID: 3-sigma bands on Bethe-Bloch curve\r
 892     //Bachelor\r
 893     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kKaon))<3.)\r
 894       isBachelorKaonForTPC = kTRUE;\r
 895     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kPion))<3.)\r
 896       isBachelorPionForTPC = kTRUE;\r
 897 \r
 898     //Negative V0 daughter\r
 899     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kPion))<3.)\r
 900       isNegPionForTPC = kTRUE;\r
 901     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kProton))<3.)\r
 902       isNegProtonForTPC = kTRUE;\r
 903     \r
 904     //Positive V0 daughter\r
 905     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kPion))<3.)\r
 906       isPosPionForTPC = kTRUE;\r
 907     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kProton))<3.)\r
 908       isPosProtonForTPC = kTRUE;\r
 909  \r
 910    \r
 911     //Extra QA information\r
 912     xi->GetPxPyPz(XiPx, XiPy, XiPz);\r
 913     XiPt = TMath::Sqrt(XiPx*XiPx + XiPy*XiPy);\r
 914     XiPtot= TMath::Sqrt(XiPx*XiPx + XiPy*XiPy + XiPz*XiPz);\r
 915     \r
 916     XiPAfter[0] = XiPx;\r
 917     XiPAfter[1] = XiPy;\r
 918     XiPAfter[2] = XiPz;\r
 919 \r
 920     xi->GetBPxPyPz(bachPx, bachPy, bachPz);\r
 921     bachPt = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy);\r
 922     bachPtot = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy + bachPz*bachPz);\r
 923 \r
 924     //chargeXi = xi->Charge();\r
 925     \r
 926     V0toXiCosOfPointingAngle \r
 927       = xi->GetV0CosineOfPointingAngle(posXi[0], posXi[1], posXi[2]);\r
 928     rapXi = xi->RapXi();\r
 929     rapOmega = xi->RapOmega();\r
 930     phi = xi->Phi();\r
 931     alphaXi = xi->AlphaXi();\r
 932     ptArmXi = xi->PtArmXi();\r
 933 \r
 934     distToVtxZBefore = posXi[2]-bestPrimaryVtxPos[2];\r
 935     distToVtxXYBefore \r
 936       = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
 937                     *(posXi[0] - bestPrimaryVtxPos[0])\r
 938                     +(posXi[1] - bestPrimaryVtxPos[1])\r
 939                     *(posXi[1] - bestPrimaryVtxPos[1]));\r
 940     \r
 941     //propagation to the best primary vertex to determine the momentum\r
 942     Propagate(bestPrimaryVtxPos, posXi, XiPAfter, b, xi->Charge());\r
 943     distToVtxZAfter = posXi[2] - bestPrimaryVtxPos[2];\r
 944     distToVtxXYAfter = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
 945                                    *(posXi[0] - bestPrimaryVtxPos[0])\r
 946                                    +(posXi[1] - bestPrimaryVtxPos[1])\r
 947                                    *(posXi[1] - bestPrimaryVtxPos[1]));\r
 948     phiAfter = TMath::Pi() + TMath::ATan2(-XiPAfter[1],-XiPAfter[0]);\r
 949     \r
 950     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZAfter"))->Fill(distToVtxZAfter);\r
 951     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYAfter"))->Fill(distToVtxXYAfter);\r
 952     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZBeforeVsAfter"))->Fill(distToVtxZBefore, distToVtxZAfter);\r
 953     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYBeforeVsAfter"))->Fill(distToVtxXYBefore, distToVtxXYAfter);\r
 954     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PxBeforeVsAfter"))->Fill(XiPx, XiPAfter[0]);\r
 955     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PyBeforeVsAfter"))->Fill(XiPy, XiPAfter[1]);\r
 956     if(xi->Charge()>0)\r
 957       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiPosBeforeVsAfter"))->Fill(phi, phiAfter);\r
 958     else if(xi->Charge()<0)\r
 959       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiNegBeforeVsAfter"))->Fill(phi, phiAfter);\r
 960 \r
 961     \r
 962     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Chi2Xi"))->Fill(chi2Xi);\r
 963     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaXiDaughters"))->Fill(dcaXiDaughters);\r
 964     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaBachToPrimVertex"))->Fill(dcaBachToPrimaryVtxXi);\r
 965     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiCosOfPointingAngle"))->Fill(XiCosOfPointingAngle);\r
 966     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiRadius"))->Fill(XiRadius);\r
 967     \r
 968     //V0\r
 969     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassLambdaAsCascDghter"))->Fill(invMassLambdaAsCascDghter);\r
 970     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0Chi2Xi"))->Fill(V0Chi2Xi);\r
 971     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0DaughtersXi"))->Fill(dcaV0DaughtersXi);\r
 972     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0CosOfPointingAngleXi"))->Fill(V0CosOfPointingAngleXi);\r
 973     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0RadiusXi"))->Fill(V0RadiusXi);\r
 974     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0ToPrimVertexXi"))->Fill(dcaV0ToPrimaryVtxXi);\r
 975     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaPosToPrimVertexXi"))->Fill(dcaPosToPrimaryVtxXi);\r
 976     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaNegToPrimVertexXi"))->Fill(dcaNegToPrimaryVtxXi);\r
 977     \r
 978     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0toXiCosOfPointingAngle"))->Fill(V0toXiCosOfPointingAngle);\r
 979     \r
 980     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Armenteros"))->Fill(alphaXi, ptArmXi);\r
 981   \r
 982     //PID cuts with TPC cuts\r
 983     if(xi->Charge() < 0){\r
 984       if(isPosProtonForTPC\r
 985          && isNegPionForTPC){\r
 986         \r
 987         switch(fSpecie) {\r
 988         case 0:\r
 989           if(isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
 990             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiMinus);\r
 991           \r
 992           //candidate inserting\r
 993             MakeTrack(invMassXiMinus, XiPt, /*xi->Phi()*/\r
 994                       phiAfter, xi->Eta(),  pTrkXi->GetID(),\r
 995                       nTrkXi->GetID(), bTrkXi->GetID());\r
 996           }\r
 997           break;\r
 998           \r
 999         case 1:\r
1000           if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8){\r
1001             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaMinus);\r
1002             MakeTrack(invMassOmegaMinus, XiPt, /*xi->Phi()*/\r
1003                       phiAfter, xi->Eta(),\r
1004                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1005           }\r
1006           break;\r
1007         }\r
1008       }\r
1009     }\r
1010 \r
1011     if(xi->Charge() > 0){\r
1012       if(isNegProtonForTPC\r
1013          && isPosPionForTPC){\r
1014         \r
1015         switch (fSpecie){\r
1016         case 0:\r
1017           if(isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
1018             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiPlus);\r
1019           \r
1020           //candidate inserting                                              \r
1021             MakeTrack(invMassXiPlus, XiPt, /*xi->Phi()*/\r
1022                       phiAfter, xi->Eta(),\r
1023                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1024           }\r
1025           break;\r
1026 \r
1027         case 1:\r
1028           if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8){\r
1029             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaPlus);\r
1030             MakeTrack(invMassOmegaPlus, XiPt, /*xi->Phi()*/\r
1031                       phiAfter, xi->Eta(),\r
1032                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1033           }\r
1034           break;\r
1035         }\r
1036       }\r
1037     }\r
1038   \r
1039   }\r
1040 \r
1041   return;\r
1042 }\r
1043 \r
1044 //______________________________________________________________________________\r
1045 void AliAnalysisTaskFlowCascade::ReadFromAODv0(AliAODEvent *fAOD)\r
1046 {\r
1047   \r
1048   fCutsRPTPC->SetEvent(fAOD, MCEvent());\r
1049   fCutsRPVZE->SetEvent(fAOD, MCEvent());\r
1050   fCutsPOI->SetEvent(fAOD, MCEvent());\r
1051   fFlowEventTPC->Fill(fCutsRPTPC, fCutsPOI);\r
1052   fFlowEventVZE->Fill(fCutsRPVZE, fCutsPOI);\r
1053 \r
1054   //  Double_t trkPrimaryVtxPos[3] = {-100., -100., -100.};\r
1055   Double_t bestPrimaryVtxPos[3] = {-100., -100., -100.};\r
1056 \r
1057   Double_t b = fAOD->GetMagneticField();\r
1058 \r
1059   int nCascades=fAOD->GetNumberOfCascades();\r
1060   const AliAODVertex *primaryBestAODVtx = fAOD->GetPrimaryVertex();\r
1061   primaryBestAODVtx->GetXYZ(bestPrimaryVtxPos);\r
1062   \r
1063   // calculation part dedicated to Xi vertices\r
1064   for(Int_t iXi = 0; iXi < nCascades; iXi++){\r
1065     Double_t effMassXi = 0.;\r
1066     Double_t chi2Xi = -1.;\r
1067     Double_t dcaXiDaughters = -1.;\r
1068     Double_t XiCosOfPointingAngle = -1.;\r
1069     Double_t posXi[3] = {-1000., -1000., -1000.};\r
1070     Double_t XiRadius = -1000.;\r
1071     \r
1072     Double_t invMassLambdaAsCascDghter = 0.;\r
1073     Double_t V0Chi2Xi = -1.;\r
1074     Double_t dcaV0DaughtersXi = -1.;\r
1075     \r
1076     Double_t dcaBachToPrimaryVtxXi = -1.;\r
1077     Double_t dcaV0ToPrimaryVtxXi = -1.;\r
1078     Double_t dcaPosToPrimaryVtxXi = -1.;\r
1079     Double_t dcaNegToPrimaryVtxXi = -1.;\r
1080     Double_t V0CosOfPointingAngleXi = -1.;\r
1081     Double_t posV0Xi[3] = {-1000., -1000., -1000.};\r
1082     Double_t V0RadiusXi = -1000.;\r
1083     //    Double_t V0quality = 0.;\r
1084 \r
1085     Double_t invMassXiMinus = 0.;\r
1086     Double_t invMassXiPlus = 0.;\r
1087     Double_t invMassOmegaMinus = 0.;\r
1088     Double_t invMassOmegaPlus = 0.;\r
1089     \r
1090     /*\r
1091     Bool_t isPosInXiProton = kFALSE;\r
1092     Bool_t isPosInXiPion = kFALSE;\r
1093     Bool_t isPosInOmegaProton = kFALSE;\r
1094     Bool_t isPosInOmegaPion = kFALSE;\r
1095     \r
1096     Bool_t isNegInXiProton = kFALSE;\r
1097     Bool_t isNegInXiPion = kFALSE;\r
1098     Bool_t isNegInOmegaProton = kFALSE;\r
1099     Bool_t isNegInOmegaPion = kFALSE;\r
1100 \r
1101     Bool_t isBachelorKaon = kFALSE;\r
1102     Bool_t isBachelorPion = kFALSE;\r
1103     */\r
1104 \r
1105 \r
1106     Bool_t isBachelorKaonForTPC = kFALSE;\r
1107     Bool_t isBachelorPionForTPC = kFALSE;\r
1108     Bool_t isNegPionForTPC = kFALSE;\r
1109     Bool_t isPosPionForTPC = kFALSE;\r
1110     Bool_t isNegProtonForTPC = kFALSE;\r
1111     Bool_t isPosProtonForTPC = kFALSE;\r
1112     \r
1113     Double_t XiPx = 0., XiPy = 0., XiPz = 0.;\r
1114     Double_t XiPt = 0.;\r
1115     Double_t XiPtot = 0.;\r
1116     \r
1117     Double_t bachPx = 0., bachPy = 0., bachPz = 0.;\r
1118     Double_t bachPt = 0.;\r
1119     Double_t bachPtot = 0.;\r
1120     \r
1121     //Short_t chargeXi = -2;\r
1122     Double_t V0toXiCosOfPointingAngle = 0.;\r
1123     \r
1124     Double_t rapXi = -20.;\r
1125     Double_t rapOmega = -20.;\r
1126     Double_t phi = 6.3;\r
1127     Double_t alphaXi = -200.;\r
1128     Double_t ptArmXi = -200.;\r
1129 \r
1130     Double_t distToVtxZBefore = -999.;\r
1131     Double_t distToVtxZAfter = -999.;\r
1132     Double_t distToVtxXYBefore = -999.;\r
1133     Double_t distToVtxXYAfter = -999.;\r
1134     Double_t XiPAfter[3] = {-999., -999., -999.};\r
1135     Double_t phiAfter = -999.;\r
1136 \r
1137     const AliAODcascade *xi = fAOD->GetCascade(iXi);\r
1138     if (!xi) continue;\r
1139 \r
1140     effMassXi = xi->MassXi(); //default working hypothesis: Xi- decay\r
1141     chi2Xi = xi->Chi2Xi();\r
1142     dcaXiDaughters = xi->DcaXiDaughters();\r
1143     XiCosOfPointingAngle = xi->CosPointingAngleXi(bestPrimaryVtxPos[0],\r
1144                                                   bestPrimaryVtxPos[1],\r
1145                                                   bestPrimaryVtxPos[2]);\r
1146     posXi[0] = xi->DecayVertexXiX();\r
1147     posXi[1] = xi->DecayVertexXiY();\r
1148     posXi[2] = xi->DecayVertexXiZ();\r
1149     XiRadius = TMath::Sqrt(posXi[0]*posXi[0]\r
1150                            +posXi[1]*posXi[1]\r
1151                            +posXi[2]*posXi[2]);\r
1152 \r
1153     AliAODTrack *pTrkXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(0) );\r
1154     AliAODTrack *nTrkXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(1) );\r
1155     AliAODTrack *bTrkXi \r
1156       = dynamic_cast<AliAODTrack*>( xi->GetDecayVertexXi()->GetDaughter(0) );\r
1157 \r
1158     if(!pTrkXi || !nTrkXi || !bTrkXi) continue;\r
1159 \r
1160     UInt_t idxPosXi  = (UInt_t) TMath::Abs( pTrkXi->GetID() );\r
1161     UInt_t idxNegXi  = (UInt_t) TMath::Abs( nTrkXi->GetID() );\r
1162     UInt_t idxBach   = (UInt_t) TMath::Abs( bTrkXi->GetID() );\r
1163 \r
1164     if(idxBach == idxNegXi || idxBach == idxPosXi) continue;\r
1165 \r
1166     if( !fCutsDau->IsSelected(pTrkXi) \r
1167         || !fCutsDau->IsSelected(nTrkXi)\r
1168         || !fCutsDau->IsSelected(bTrkXi) ) continue;\r
1169 \r
1170     \r
1171     if(pTrkXi->IsOn(AliESDtrack::kTPCin)){\r
1172       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(pTrkXi->P()*pTrkXi->Charge(), pTrkXi->GetTPCsignal());\r
1173     }\r
1174     if( nTrkXi->IsOn(AliESDtrack::kTPCin) ){\r
1175       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(nTrkXi->P()*nTrkXi->Charge(), nTrkXi->GetTPCsignal());\r
1176     }\r
1177     if(bTrkXi->IsOn(AliESDtrack::kTPCin)){\r
1178       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("TPCdEdxOfCascDghters"))->Fill(bTrkXi->P()*bTrkXi->Charge(), bTrkXi->GetTPCsignal());\r
1179     }\r
1180     \r
1181     if(xi->ChargeXi() < 0)\r
1182       invMassLambdaAsCascDghter = xi->MassLambda();\r
1183     else\r
1184       invMassLambdaAsCascDghter = xi->MassAntiLambda();\r
1185     \r
1186     dcaV0DaughtersXi = xi->DcaV0Daughters();\r
1187     V0Chi2Xi = xi->Chi2V0();\r
1188     V0CosOfPointingAngleXi \r
1189       = xi->CosPointingAngle(bestPrimaryVtxPos);\r
1190     dcaV0ToPrimaryVtxXi = xi->DcaV0ToPrimVertex();\r
1191     dcaBachToPrimaryVtxXi = xi->DcaBachToPrimVertex();\r
1192     \r
1193     //V0\r
1194     posV0Xi[0] = xi->DecayVertexV0X();\r
1195     posV0Xi[1] = xi->DecayVertexV0Y();\r
1196     posV0Xi[2] = xi->DecayVertexV0Z();\r
1197     V0RadiusXi = TMath::Sqrt(posV0Xi[0]*posV0Xi[0]\r
1198                              +posV0Xi[1]*posV0Xi[1]\r
1199                              +posV0Xi[2]*posV0Xi[2]);\r
1200     dcaPosToPrimaryVtxXi = xi->DcaPosToPrimVertex();\r
1201     dcaNegToPrimaryVtxXi = xi->DcaNegToPrimVertex();\r
1202 \r
1203     //apply cuts ?\r
1204     // if(XiRadius < 1. || XiRadius > 100.) continue;\r
1205     //if(dcaXiDaughters > 0.1) continue;\r
1206     //if(XiCosOfPointingAngle < 0.999) continue;\r
1207     //if(dcaV0ToPrimaryVtxXi < 0.05) continue;\r
1208     //if(dcaBachToPrimaryVtxXi < 0.03) continue;\r
1209 \r
1210     if(dcaXiDaughters > fCascadeCuts[0]) continue;\r
1211     if(XiCosOfPointingAngle < fCascadeCuts[1]) continue;\r
1212     if(dcaV0ToPrimaryVtxXi < fCascadeCuts[2]) continue;\r
1213     if(dcaBachToPrimaryVtxXi < fCascadeCuts[3]) continue;\r
1214     \r
1215     //V0 mass cut?\r
1216     //if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > 0.006) continue;\r
1217     if(TMath::Abs(invMassLambdaAsCascDghter-1.11568) > fCascadeCuts[7]) \r
1218       continue;\r
1219 \r
1220     //if(dcaV0DaughtersXi > 1.) continue;\r
1221     //if(V0CosOfPointingAngleXi > 0.9999) continue;\r
1222     //if(dcaPosToPrimaryVtxXi < 0.1) continue;\r
1223     //if(dcaNegToPrimaryVtxXi < 0.1) continue;\r
1224     if(dcaV0DaughtersXi > fCascadeCuts[4]) continue;\r
1225     if(V0CosOfPointingAngleXi > fCascadeCuts[5]) continue;\r
1226     if(dcaPosToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
1227     if(dcaNegToPrimaryVtxXi < fCascadeCuts[6]) continue;\r
1228     \r
1229     // if(V0RadiusXi < 1.0 || V0RadiusXi > 100) continue;\r
1230     \r
1231     //other cuts?\r
1232 \r
1233 \r
1234     //???\r
1235     if(xi->ChargeXi()<0){\r
1236       invMassXiMinus = xi->MassXi();\r
1237       invMassOmegaMinus = xi->MassOmega();\r
1238     }else{\r
1239       invMassXiPlus = xi->MassXi();\r
1240       invMassOmegaPlus = xi->MassOmega();\r
1241     }\r
1242 \r
1243     /*\r
1244     if(pTrkXi->GetMostProbablePID() == AliAODTrack::kProton) {\r
1245       isPosInXiProton = kTRUE;\r
1246       isPosInOmegaProton = kTRUE;\r
1247     }\r
1248     if(pTrkXi->GetMostProbablePID() == AliAODTrack::kPion){\r
1249       isPosInXiPion = kTRUE;\r
1250       isPosInOmegaPion = kTRUE;\r
1251     }\r
1252     \r
1253     if(nTrkXi->GetMostProbablePID() == AliAODTrack::kPion){\r
1254       isNegInXiPion = kTRUE;\r
1255       isNegInOmegaPion = kTRUE;\r
1256     }\r
1257     if(nTrkXi->GetMostProbablePID() == AliAODTrack::kProton){\r
1258       isNegInXiProton = kTRUE;\r
1259       isNegInOmegaProton = kTRUE;\r
1260     }\r
1261 \r
1262     if(bTrkXi->GetMostProbablePID() == AliAODTrack::kPion)\r
1263       isBachelorPion = kTRUE;\r
1264     if(bTrkXi->GetMostProbablePID() == AliAODTrack::kKaon)\r
1265       isBachelorKaon = kTRUE;\r
1266     */\r
1267 \r
1268     //PID with TPC only: \r
1269     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kKaon))<3.)\r
1270       isBachelorKaonForTPC = kTRUE;\r
1271     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(bTrkXi, AliPID::kPion))<3.)\r
1272       isBachelorPionForTPC = kTRUE;\r
1273 \r
1274     //Negative V0 daughter\r
1275     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kPion))<3.)\r
1276       isNegPionForTPC = kTRUE;\r
1277     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrkXi, AliPID::kProton))<3.)\r
1278       isNegProtonForTPC = kTRUE;\r
1279     \r
1280     //Positive V0 daughter\r
1281     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kPion))<3.)\r
1282       isPosPionForTPC = kTRUE;\r
1283     if(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrkXi, AliPID::kProton))<3.)\r
1284       isPosProtonForTPC = kTRUE;\r
1285 \r
1286     //Extra QA information\r
1287     XiPx = xi->MomXiX();\r
1288     XiPy = xi->MomXiY();\r
1289     XiPz = xi->MomXiZ();\r
1290     XiPt = TMath::Sqrt(XiPx*XiPx + XiPy*XiPy);\r
1291     XiPtot= TMath::Sqrt(XiPx*XiPx + XiPy*XiPy + XiPz*XiPz);\r
1292     \r
1293     bachPx = xi->MomBachX();\r
1294     bachPy = xi->MomBachY();\r
1295     bachPz = xi->MomBachZ();\r
1296     \r
1297     bachPt = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy);\r
1298     bachPtot = TMath::Sqrt(bachPx*bachPx + bachPy*bachPy + bachPz*bachPz);\r
1299   \r
1300     V0toXiCosOfPointingAngle = xi->CosPointingAngle( xi->GetDecayVertexXi() );\r
1301     \r
1302     rapXi = xi->RapXi();\r
1303     rapOmega = xi->RapOmega();\r
1304     phi = xi->Phi();\r
1305     alphaXi = xi->AlphaXi();\r
1306     ptArmXi = xi->PtArmXi();\r
1307 \r
1308     distToVtxZBefore = posXi[2]-bestPrimaryVtxPos[2];\r
1309     distToVtxXYBefore\r
1310       = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
1311                     *(posXi[0] - bestPrimaryVtxPos[0])\r
1312                     +(posXi[1] - bestPrimaryVtxPos[1])\r
1313                     *(posXi[1] - bestPrimaryVtxPos[1]));\r
1314 \r
1315 \r
1316     XiPAfter[0] = XiPx;\r
1317     XiPAfter[1] = XiPy;\r
1318     XiPAfter[2] = XiPz;\r
1319     //propagation to the best primary vertex to determine the momentum\r
1320     Propagate(bestPrimaryVtxPos, posXi, XiPAfter, b, xi->ChargeXi());\r
1321     distToVtxZAfter = posXi[2] - bestPrimaryVtxPos[2];\r
1322     distToVtxXYAfter = TMath::Sqrt((posXi[0] - bestPrimaryVtxPos[0])\r
1323                                    *(posXi[0] - bestPrimaryVtxPos[0])\r
1324                                    +(posXi[1] - bestPrimaryVtxPos[1])\r
1325                                    *(posXi[1] - bestPrimaryVtxPos[1]));\r
1326     phiAfter = TMath::Pi() + TMath::ATan2(-XiPAfter[1],-XiPAfter[0]);\r
1327 \r
1328     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZAfter"))->Fill(distToVtxZAfter);\r
1329     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYAfter"))->Fill(distToVtxXYAfter);\r
1330     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxZBeforeVsAfter"))->Fill(distToVtxZBefore, distToVtxZAfter);\r
1331     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DistToVtxXYBeforeVsAfter"))->Fill(distToVtxXYBefore, distToVtxXYAfter);\r
1332     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PxBeforeVsAfter"))->Fill(XiPx, XiPAfter[0]);\r
1333     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PyBeforeVsAfter"))->Fill(XiPy, XiPAfter[1]);\r
1334     if(xi->ChargeXi()>0)\r
1335       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiPosBeforeVsAfter"))->Fill(phi, phiAfter);\r
1336     else if(xi->ChargeXi()<0)\r
1337       ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("PhiNegBeforeVsAfter"))->Fill(phi, phiAfter);\r
1338     \r
1339     //for default hypothesis\r
1340     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Chi2Xi"))->Fill(chi2Xi);\r
1341     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaXiDaughters"))->Fill(dcaXiDaughters);\r
1342     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaBachToPrimVertex"))->Fill(dcaBachToPrimaryVtxXi);\r
1343     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiCosOfPointingAngle"))->Fill(XiCosOfPointingAngle);\r
1344     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("XiRadius"))->Fill(XiRadius);\r
1345     \r
1346     //V0\r
1347     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassLambdaAsCascDghter"))->Fill(invMassLambdaAsCascDghter);\r
1348     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0Chi2Xi"))->Fill(V0Chi2Xi);\r
1349     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0DaughtersXi"))->Fill(dcaV0DaughtersXi);\r
1350     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0CosOfPointingAngleXi"))->Fill(V0CosOfPointingAngleXi);\r
1351     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0RadiusXi"))->Fill(V0RadiusXi);\r
1352     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaV0ToPrimVertexXi"))->Fill(dcaV0ToPrimaryVtxXi);\r
1353     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaPosToPrimVertexXi"))->Fill(dcaPosToPrimaryVtxXi);\r
1354     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("DcaNegToPrimVertexXi"))->Fill(dcaNegToPrimaryVtxXi);\r
1355     \r
1356     ((TH1F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("V0toXiCosOfPointingAngle"))->Fill(V0toXiCosOfPointingAngle);\r
1357     \r
1358     ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("Armenteros"))->Fill(alphaXi, ptArmXi);\r
1359   \r
1360     //with PID cuts\r
1361     if(xi->ChargeXi()<0){\r
1362       if(isPosProtonForTPC && isNegPionForTPC){\r
1363         switch (fSpecie){\r
1364         case 0:\r
1365           if( isBachelorPionForTPC && TMath::Abs(rapXi) < 0.8){\r
1366             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiMinus);\r
1367             MakeTrack(invMassXiMinus, XiPt, /*xi->Phi(),*/\r
1368                       phiAfter, xi->Eta(), \r
1369                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1370           }// endif\r
1371       \r
1372           break;\r
1373 \r
1374         case 1:\r
1375           if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8\r
1376              && (invMassXiMinus > 1.32486 || invMassXiMinus < 1.30486)){\r
1377             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaMinus);      \r
1378             \r
1379             MakeTrack(invMassOmegaMinus, XiPt, /* xi->Phi(),*/\r
1380                       phiAfter,   xi->Eta(),\r
1381                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1382           }//endif\r
1383           break;\r
1384         }\r
1385       }\r
1386     }//end if ChargeXi()<0\r
1387     \r
1388     if(xi->ChargeXi() > 0){\r
1389       if(isNegProtonForTPC && isPosPionForTPC){ \r
1390         switch(fSpecie){\r
1391         case 0:\r
1392           if (isBachelorPionForTPC  && TMath::Abs(rapXi) < 0.8){\r
1393             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassXiPlus);\r
1394           \r
1395             //candidate inserting                                              \r
1396             MakeTrack(invMassXiPlus, XiPt, /* xi->Phi(),*/\r
1397                       phiAfter, xi->Eta(),\r
1398                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1399           }//endif particle id\r
1400           break;\r
1401           \r
1402         case 1:\r
1403           if(isBachelorKaonForTPC && TMath::Abs(rapOmega) < 0.8\r
1404              && (invMassXiPlus > 1.32486 || invMassXiPlus < 1.30486)){\r
1405             ((TH2F*)((TList*)fQAList->FindObject("Candidates"))->FindObject("MassVsPtAll"))->Fill(XiPt, invMassOmegaPlus);\r
1406             MakeTrack(invMassOmegaPlus, XiPt, /* xi->Phi(),*/\r
1407                       phiAfter, xi->Eta(),\r
1408                       pTrkXi->GetID(), nTrkXi->GetID(), bTrkXi->GetID());\r
1409           }//endif particle id\r
1410         }\r
1411       }\r
1412     }//endif ChargeXi()>0  \r
1413   }//for Xi candidate loop\r
1414 \r
1415   return;\r
1416 \r
1417 }\r
1418 \r
1419 \r
1420 void AliAnalysisTaskFlowCascade::MakeTrack( double mass, \r
1421                                             double pt, \r
1422                                             double phi, \r
1423                                             double eta, \r
1424                                             int iid, \r
1425                                             int jid,\r
1426                                             int kid) {\r
1427   // create track for flow tasks        \r
1428   if(fCandidates->GetLast()+1>=fCandidates->GetSize()) {\r
1429     fCandidates->Expand( 2*fCandidates->GetSize() );\r
1430   }\r
1431   Bool_t overwrite = kTRUE;\r
1432 \r
1433   AliFlowCandidateTrack *sTrack \r
1434     = (static_cast<AliFlowCandidateTrack*> (fCandidates->At( fCandidates->GetLast()+1 )));\r
1435   if( !sTrack ) { // creates new\r
1436     sTrack = new AliFlowCandidateTrack();\r
1437     overwrite = kFALSE;\r
1438   } else { // overwrites\r
1439     sTrack->ClearMe();\r
1440   }\r
1441 \r
1442 \r
1443   sTrack->SetMass(mass);\r
1444   sTrack->SetPt(pt);\r
1445   sTrack->SetPhi(phi);\r
1446   sTrack->SetEta(eta);\r
1447   sTrack->AddDaughter(iid);\r
1448   sTrack->AddDaughter(jid);\r
1449   sTrack->AddDaughter(kid);\r
1450   sTrack->SetForPOISelection(kTRUE);\r
1451   sTrack->SetForRPSelection(kFALSE);\r
1452  \r
1453   if(overwrite) {\r
1454     fCandidates->SetLast( fCandidates->GetLast()+1 );\r
1455   } else {\r
1456     fCandidates->AddLast(sTrack);\r
1457   }\r
1458   \r
1459   return;\r
1460 }\r
1461 //_____________________________________________________________________________\r
1462 void AliAnalysisTaskFlowCascade::Terminate(Option_t *)\r
1463 {\r
1464 \r
1465 }\r
1466 \r
1467 void AliAnalysisTaskFlowCascade::Propagate(Double_t vv[3], \r
1468                                            Double_t x[3], \r
1469                                            Double_t p[3], \r
1470                                            Double_t bz, \r
1471                                            Short_t sign){\r
1472   //Propagation to the primary vertex to determine the px and py\r
1473   //x, p are the position and momentum as input and output\r
1474   //bz is the magnetic field along z direction\r
1475   //sign is the charge of particle for propagation\r
1476 \r
1477   Double_t pp = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);\r
1478   Double_t len = (vv[2]-x[2])*pp/p[2];\r
1479   Double_t a = -kB2C*bz*sign;  \r
1480 \r
1481   Double_t rho = a/pp;\r
1482   x[0] += p[0]*TMath::Sin(rho*len)/a - p[1]*(1-TMath::Cos(rho*len))/a;\r
1483   x[1] += p[1]*TMath::Sin(rho*len)/a + p[0]*(1-TMath::Cos(rho*len))/a;\r
1484   x[2] += p[2]*len/pp;\r
1485 \r
1486   Double_t p0=p[0];\r
1487   p[0] = p0  *TMath::Cos(rho*len) - p[1]*TMath::Sin(rho*len);\r
1488   p[1] = p[1]*TMath::Cos(rho*len) + p0  *TMath::Sin(rho*len);\r
1489 }\r
1490 \r
1491 \r
1492 //=====================================================================       \r
1493 void AliAnalysisTaskFlowCascade::SetCommonConstants(Int_t massBins, \r
1494                                                     Double_t minMass, \r
1495                                                     Double_t maxMass)\r
1496 {\r
1497   // setter for mass bins                          \r
1498                           \r
1499   fMassBins = massBins;\r
1500   fMinMass = minMass;\r
1501   fMaxMass = maxMass;\r
1502 }\r
1503 \r
1504 \r
1505 //====================================================================         \r
1506 void AliAnalysisTaskFlowCascade::SetCuts2010(int set) {\r
1507 \r
1508   // fCascadeCuts[0]: DcaXiDaughter; fCascadeCuts[1]: XiCosOfPointingAngle\r
1509   // fCascadeCuts[2]: DcaV0ToPrimaryVtxXi; fCascadeCuts[3]: DcaBachToPrimaryVtxXi\r
1510   // fCascadeCuts[4]: DcaV0DaughersXi; fCascadeCuts[5]: V0CosOfPointingAngleXi\r
1511   // fCascadeCuts[6]: DcaV0DaughterToPrimaryVtxXi; fCascadeCuts[7]: V0MassWidth\r
1512   \r
1513   switch(set){\r
1514   \r
1515   case 0: //tighter\r
1516     fCascadeCuts[0] = 0.2; fCascadeCuts[1] = 0.999;\r
1517     fCascadeCuts[2] = 0.03; fCascadeCuts[3] = 0.05;\r
1518     fCascadeCuts[4] = .5; fCascadeCuts[5] = 0.9998;\r
1519     fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.006;\r
1520     break;\r
1521     \r
1522   case 1: //middle\r
1523     fCascadeCuts[0] = 0.3; fCascadeCuts[1] = 0.99;\r
1524     fCascadeCuts[2] = 0.01; fCascadeCuts[3] = 0.03;\r
1525     fCascadeCuts[4] = .6; fCascadeCuts[5] = 0.9999;\r
1526     fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.008;\r
1527     break;\r
1528 \r
1529   case 2: //looser\r
1530     fCascadeCuts[0] = 0.3; fCascadeCuts[1] = 0.99;\r
1531     fCascadeCuts[2] = 0.01; fCascadeCuts[3] = 0.03;\r
1532     fCascadeCuts[4] = 1.; fCascadeCuts[5] = 1.;\r
1533     fCascadeCuts[6] = 0.1; fCascadeCuts[7] = 0.01;\r
1534     break;\r
1535   }\r
1536   \r
1537 }\r