[u/mrichter/AliRoot.git] / PWG / FLOW / Tasks / AliAnalysisTaskFlowK0Candidates.cxx

/*************************************************************************
* 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.                  *
**************************************************************************/

/////////////////////////////////////////////////////
// AliAnalysisTaskFlowK0Candidates:
// Analysis task to select K0 candidates for flow analysis.
// Uses one AliESDtrackCuts object for both daughters and
// QA histograms to monitor the reconstruction.
// Author: Carlos Perez (cperez@cern.ch)
//////////////////////////////////////////////////////

#include "TChain.h"
#include "TList.h"
#include "TH1D.h"
#include "TVector3.h"

#include "AliAnalysisTaskSE.h"
#include "AliAnalysisManager.h"

#include "AliESDEvent.h"
#include "AliESDInputHandler.h"
#include "AliESDtrack.h"
#include "AliESDtrackCuts.h"

#include "AliAODEvent.h"
#include "AliAODInputHandler.h"
#include "AliAODTrack.h"

#include "AliCFManager.h"

#include "TMath.h"
#include "TObjArray.h"
#include "AliFlowCandidateTrack.h"
#include "AliFlowEventCuts.h"

#include "AliAnalysisTaskFlowK0Candidates.h"


ClassImp(AliAnalysisTaskFlowK0Candidates)

//_____________________________________________________________________________
AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates() :
  AliAnalysisTaskSE(),
  fCutsEvent(NULL),
  fCuts(NULL),
  fQAList(NULL),
  fMassMin(0),
  fMassMax(0),
  fDLcut(0),
  fEvent(NULL),
  fMulti(NULL)
{
  //
  for(int i=0; i!=4; ++i) {
    fMass[i] = NULL;
    fDCA[i]  = NULL;
    fDL[i]   = NULL;
    fCTP[i]  = NULL;
    fd0d0[i] = NULL;
    fPhi[i]  = NULL;
    fEta[i]  = NULL;
    fPt[i]   = NULL;
    fAPhi[i] = NULL;
    fAEta[i] = NULL;
    fAPt[i]  = NULL;
    fBPhi[i] = NULL;
    fBEta[i] = NULL;
    fBPt[i]  = NULL;
  }
}

//_____________________________________________________________________________
AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates(const char *name, AliFlowEventCuts *cutsEvent, AliESDtrackCuts *cuts, Double_t MassMin, Double_t MassMax) :
  AliAnalysisTaskSE(name),
  fCutsEvent(cutsEvent),
  fCuts(cuts),
  fQAList(NULL),
  fMassMin(MassMin),
  fMassMax(MassMax),
  fDLcut(0),
  fEvent(NULL),
  fMulti(NULL)
{
  //
  for(int i=0; i!=4; ++i) {
    fMass[i] = NULL;
    fDCA[i]  = NULL;
    fDL[i]   = NULL;
    fCTP[i]  = NULL;
    fd0d0[i] = NULL;
    fPhi[i]  = NULL;
    fEta[i]  = NULL;
    fPt[i]   = NULL;
    fAPhi[i] = NULL;
    fAEta[i] = NULL;
    fAPt[i]  = NULL;
    fBPhi[i] = NULL;
    fBEta[i] = NULL;
    fBPt[i]  = NULL;
  }
  DefineInput(  0, TChain::Class() );
  DefineOutput( 1, TObjArray::Class() );
  DefineOutput( 2, TList::Class() );
}

//_____________________________________________________________________________
AliAnalysisTaskFlowK0Candidates::~AliAnalysisTaskFlowK0Candidates()
{
  //
  if(fQAList) delete fQAList;
}

//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::UserCreateOutputObjects()
{
  //
  fQAList = new TList();
  fQAList->SetOwner();
  AddQAEvents();
  AddQACandidates();
  PostData( 2, fQAList);
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::AddQAEvents()
{
  //
  TList *tQAEvents = new TList();
  tQAEvents->SetName("Events");
  tQAEvents->SetOwner();
  fEvent = new TH1D("Event", "Number of Events",    2, 0, 2);
  tQAEvents->Add(fEvent);
  fMulti = new TH1D("Multiplicity", "Multiplicity", 180, 0, 10000);
  tQAEvents->Add(fMulti);
  fQAList->Add(tQAEvents);
}  
//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::AddQACandidates()
{
  //
  TList *tQACandidates[4];
  TList *tQADaughters[4];
  for(int i=0; i!=4; ++i) {
    tQACandidates[i] = new TList();
    tQACandidates[i]->SetOwner();
    tQACandidates[i]->SetName(Form("Candidates%d",i));
    fMass[i] = new TH1D( Form("Mass%i",i), "Mass;M_{#pi#pi} [GeV];Counts per MeV", 180, 0.41, 0.59); tQACandidates[i]->Add( fMass[i] );
    fDCA[i]  = new TH1D( Form("DCA%i" ,i), "DCA;[cm];Counts per 10 um",            180, 0.00, 0.18); tQACandidates[i]->Add( fDCA[i]  );
    fDL[i]   = new TH1D( Form("DL%i"  ,i), "Decay Length;[cm];Counts per 0.1 mm",  180, 0.00, 1.80); tQACandidates[i]->Add( fDL[i]   );
    fCTP[i]  = new TH1D( Form("CTP%i" ,i), "Cos#theta_{p}",                        180,-1.10, 1.10); tQACandidates[i]->Add( fCTP[i]  );
    fd0d0[i] = new TH1D( Form("d0d0%i",i), "d_{0}xd_{0};[cm^{2}];Cnts 0.01 mm^{2}",180,-0.009,0.009);tQACandidates[i]->Add( fd0d0[i] );
    fPhi[i]  = new TH1D( Form("Phi%i" ,i), "Phi;[rad];Counts per degree",     180,0,TMath::TwoPi()); tQACandidates[i]->Add( fPhi[i]  );
    fEta[i]  = new TH1D( Form("Eta%i" ,i), "Eta;;Counts per 0.04",                 180,-3.60, 3.60); tQACandidates[i]->Add( fEta[i]  );
    fPt[i]   = new TH1D( Form("Pt%i"  ,i), "Pt;[GeV];Counts per 0.1 GeV",          180, 0.00,18.00); tQACandidates[i]->Add( fPt[i]   );
    tQADaughters[i] = new TList();
    tQADaughters[i]->SetOwner();
    tQADaughters[i]->SetName(Form("Daughters%d",i));
    fAPhi[i] = new TH1D( Form("PhiBef%i",i), "Phi prePropagation;[rad];Cnts per degree",180,0,TMath::TwoPi()); tQADaughters[i]->Add( fAPhi[i] );
    fAEta[i] = new TH1D( Form("EtaBef%i",i), "Eta prePropagation;;Counts per 0.04"     ,180,-3.6,3.6); tQADaughters[i]->Add( fAEta[i] );
    fAPt[i]  = new TH1D( Form("PtBef%i" ,i), "Pt prePropagation;[GeV];Counts per 0.1 GeV",180, 0, 18); tQADaughters[i]->Add( fAPt[i]  );
    fBPhi[i] = new TH1D( Form("PhiAft%i",i), "Phi posPropagation;[rad];Cnts per degree",180,0,TMath::TwoPi()); tQADaughters[i]->Add( fBPhi[i] );
    fBPhi[i] = new TH1D( Form("EtaAft%i",i), "Eta posPropagation;;Counts per 0.04"     ,180,-3.6,3.6); tQADaughters[i]->Add( fBPhi[i] );
    fBPt[i]  = new TH1D( Form("PtAft%i" ,i), "Pt posPropagation;[GeV];Counts per 0.1 GeV",180, 0, 18); tQADaughters[i]->Add( fBPt[i]  );
    tQACandidates[i]->Add(tQADaughters[i]);
    fQAList->Add(tQACandidates[i]);
  }
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::NotifyRun()
{
  //
//  AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());
//  if(!fESD) return;
}

//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::UserExec(Option_t *)
{
  //
  AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());
  AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
  if( (!fESD)&&(!fAOD) )
    return;
//  printf("\nEvent found");
  fEvent->Fill( 0 );
  if(fCutsEvent)
    if( !(fCutsEvent->IsSelected(InputEvent())) )
      return;
  fEvent->Fill( 1 );
  if(fESD)
    ReadFromESD(fESD);
  else if(fAOD)
    ReadFromAOD(fAOD);
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::ReadFromESD(const AliESDEvent *fESD)
{
  //
  TObjArray *gPOIselection =  new TObjArray();
  gPOIselection->SetOwner();
  int nTracks = fESD->GetNumberOfTracks();
  fMulti->Fill( nTracks );
  for(int i=0; i!=nTracks; ++i) {
    AliESDtrack *ioT = fESD->GetTrack(i);
    if(fCuts)
      if( !(fCuts->IsSelected(ioT)) )
        continue;
    // printf("\n first particle OK...");
    for(int j=i+1; j!=nTracks; ++j) {
      AliESDtrack *joT = fESD->GetTrack(j);
      if( (ioT->Charge()*joT->Charge()) > 0 )
        continue;
      if(fCuts)
        if( !(fCuts->IsSelected(joT)) )
          continue;
      // printf("\n second also...");
      AliESDtrack *iT = new AliESDtrack(*ioT);
      AliESDtrack *jT = new AliESDtrack(*joT);
      // getting distance of closest approach
      double DCA = iT->PropagateToDCA(jT,fESD->GetMagneticField());
      // printf(" propagated...");
      // getting decay length
      TVector3 vp, vv, vl;
      vp = TVector3( (iT->Xv()+jT->Xv())/2, (iT->Yv()+jT->Yv())/2, (iT->Zv()+jT->Zv())/2 ); // candidate position
      vv = TVector3( fESD->GetPrimaryVertex()->GetX(), fESD->GetPrimaryVertex()->GetY(), fESD->GetPrimaryVertex()->GetZ() ); // vertex position
      vl = vp - vv; // flight line
      double gDL;
      gDL  = vl.Mag();
      // getting cos pointing angle
      double gCTP;
      TVector3 vi, vj, vs;
      vi = TVector3( iT->Px(), iT->Py(), iT->Pz() );
      vj = TVector3( jT->Px(), jT->Py(), jT->Pz() );
      vs = vi + vj;
      gCTP = TMath::Cos( vl.Angle(vs) ); // Marta says: "Carlos, it is faster if you calculate this by hand!"
      double gD0D0;
      gD0D0 = iT->GetD(vv.X(),vv.Y(),fESD->GetMagneticField())*jT->GetD(vv.X(),vv.Y(),fESD->GetMagneticField());
      // getting invariant mass
      double sum12 = iT->P()*iT->P()+jT->P()*jT->P();
      double pro12 = iT->P()*iT->P()*jT->P()*jT->P();
      double gInvMass = 0;
      gInvMass += TMath::Power(0.13957018,2);
      gInvMass += sqrt( TMath::Power(0.13957018,4) + TMath::Power(0.13957018,2)*(sum12) + pro12 );
      gInvMass -= ( iT->Px()*jT->Px()+iT->Py()*jT->Py()+iT->Pz()*jT->Pz() );
      gInvMass *= 2;
      gInvMass = sqrt(gInvMass);
      // filtering
      int iLevel = 3;
      if( gCTP<0.8 )
        iLevel = 2;
      if( gDL<fDLcut ) // 0.5
        iLevel = 1;
      if( DCA>0.05 )
        iLevel = 0;
      // printf(" candidate at level %d...",iLevel);
      for(int h=0; h!=iLevel+1; ++h) {
        // printf(" %d",h);
        // candidate
        fDCA[h]->Fill( DCA );
        fDL[h]->Fill( gDL );
        fCTP[h]->Fill( gCTP );
        fd0d0[h]->Fill( gD0D0 );
        fMass[h]->Fill( gInvMass );
        fPhi[h]->Fill( vs.Phi()+TMath::Pi() );
        fEta[h]->Fill( vs.Eta() );
        fPt[h]->Fill( vs.Pt() );
        // daughters
        fAPhi[h]->Fill( iT->Phi() );
        fAEta[h]->Fill( iT->Eta() );
        fAPt[h]->Fill(  iT->Pt()  );
        fAPhi[h]->Fill( jT->Phi() );
        fAEta[h]->Fill( jT->Eta() );
        fAPt[h]->Fill(  jT->Pt()  );
        fBPhi[h]->Fill( ioT->Phi());
        fBPhi[h]->Fill( ioT->Eta());
        fBPt[h]->Fill(  ioT->Pt() );
        fBPhi[h]->Fill( joT->Phi());
        fBPhi[h]->Fill( joT->Eta());
        fBPt[h]->Fill(  joT->Pt() );
      }
      // building candidate
      if( (gInvMass>fMassMin) && (gInvMass<fMassMax) && (iLevel==3) ) {
        AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
        sTrack->SetMass( gInvMass );
        sTrack->SetPt( vs.Pt() );
        sTrack->SetPhi( vs.Phi()+TMath::Pi() );
        sTrack->SetEta( vs.Eta() );
        sTrack->SetCharge( 0 );
        sTrack->AddDaughter( iT->GetID() );
        sTrack->AddDaughter( jT->GetID() );
        gPOIselection->AddLast( sTrack );
      }
      delete iT;
      delete jT;
    } // endfor j
  } // endfor i
  PostData( 1, gPOIselection );
  PostData( 2, fQAList);
}

//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::ReadFromAOD(const AliAODEvent *fAOD)
{
  //
  TObjArray *gPOIselection =  new TObjArray();
  gPOIselection->SetOwner();
  int nTracks = fAOD->GetNumberOfTracks();
  fMulti->Fill( nTracks );
  for(int i=0; i!=fAOD->GetNumberOfV0s(); ++i) {
    AliAODv0 *iV0 = (AliAODv0*) fAOD->GetV0(i);
    if ( iV0->GetOnFlyStatus() ) continue;
    if ( (iV0->ChargeProng(0)*iV0->ChargeProng(1))>0 ) continue;
    // getting distance of closest approach
    double DCA = iV0->DcaV0Daughters();
    // getting decay length
    double gDL;
    double vv[3];
    vv[0] = fAOD->GetPrimaryVertex()->GetX();
    vv[1] = fAOD->GetPrimaryVertex()->GetY();
    vv[2] = fAOD->GetPrimaryVertex()->GetZ();
    gDL = iV0->DecayLengthV0(vv);
    // cos pointing angle
    double gCTP;
    gCTP = iV0->CosPointingAngle(vv);
    double gD0D0;
    gD0D0 = iV0->Prodd0d0();
    // getting invariant mass
    double gInvMass = iV0->MassK0Short();
    // filtering
    int iLevel = 3;
    if( gCTP<0.8 )
      iLevel = 2;
    if( gDL<fDLcut ) // 0.5
      iLevel = 1;
    if( DCA>0.05 )
      iLevel = 0;
    for(int h=0; h!=iLevel+1; ++h) {
      // candidate
      fDCA[h]->Fill( DCA );
      fDL[h]->Fill( gDL );
      fCTP[h]->Fill( gCTP );
      fd0d0[h]->Fill( gD0D0 );
      fMass[h]->Fill( gInvMass );
      fPhi[h]->Fill( iV0->Phi() );
      fEta[h]->Fill( iV0->Eta() );
      fPt[h]->Fill( iV0->Pt() );
      // daughters
      fAPhi[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Phi() );
      fAEta[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Eta() );
      fAPt[h]->Fill(  ( (AliAODTrack*) iV0->GetDaughter(0) )->Pt()  );
      fAPhi[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Phi() );
      fAEta[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Eta() );
      fAPt[h]->Fill(  ( (AliAODTrack*) iV0->GetDaughter(1) )->Pt()  );
      fBPhi[h]->Fill( iV0->PhiProng(0) );
      fBPhi[h]->Fill( iV0->EtaProng(0) );
      fBPt[h]->Fill(  iV0->PtProng(0)  );
      fBPhi[h]->Fill( iV0->PhiProng(1) );
      fBPhi[h]->Fill( iV0->EtaProng(1) );
      fBPt[h]->Fill(  iV0->PtProng(1)  );
    }
    AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
    sTrack->SetMass( gInvMass );
    sTrack->SetPt( iV0->Pt() );
    sTrack->SetPhi( iV0->Phi() );
    sTrack->SetEta( iV0->Eta() );
    sTrack->SetCharge( 0 );
    sTrack->AddDaughter( iV0->GetPosID() );
    sTrack->AddDaughter( iV0->GetNegID() );
    gPOIselection->AddLast( sTrack );
  }
  PostData( 1, gPOIselection );
  PostData( 2, fQAList);
}

//_____________________________________________________________________________
void AliAnalysisTaskFlowK0Candidates::Terminate(Option_t *)
{
  //
}
Commit	Line	Data
6e214c87	1	/*************************************************************************
	2	* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	/////////////////////////////////////////////////////
	17	// AliAnalysisTaskFlowK0Candidates:
	18	// Analysis task to select K0 candidates for flow analysis.
	19	// Uses one AliESDtrackCuts object for both daughters and
	20	// QA histograms to monitor the reconstruction.
	21	// Author: Carlos Perez (cperez@cern.ch)
	22	//////////////////////////////////////////////////////
	23
	24	#include "TChain.h"
	25	#include "TList.h"
	26	#include "TH1D.h"
	27	#include "TVector3.h"
	28
	29	#include "AliAnalysisTaskSE.h"
	30	#include "AliAnalysisManager.h"
	31
	32	#include "AliESDEvent.h"
	33	#include "AliESDInputHandler.h"
	34	#include "AliESDtrack.h"
	35	#include "AliESDtrackCuts.h"
	36
	37	#include "AliAODEvent.h"
	38	#include "AliAODInputHandler.h"
	39	#include "AliAODTrack.h"
	40
	41	#include "AliCFManager.h"
	42
41dc4195	43	#include "TMath.h"
6e214c87	44	#include "TObjArray.h"
	45	#include "AliFlowCandidateTrack.h"
	46	#include "AliFlowEventCuts.h"
	47
	48	#include "AliAnalysisTaskFlowK0Candidates.h"
	49
	50
	51	ClassImp(AliAnalysisTaskFlowK0Candidates)
	52
	53	//_____________________________________________________________________________
	54	AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates() :
	55	AliAnalysisTaskSE(),
	56	fCutsEvent(NULL),
	57	fCuts(NULL),
41dc4195	58	fQAList(NULL),
6e214c87	59	fMassMin(0),
41dc4195	60	fMassMax(0),
41dc4195	61	fDLcut(0),
26baa5fb	62	fEvent(NULL),
26baa5fb	63	fMulti(NULL)
6e214c87	64	{
26baa5fb	65	//
41dc4195	66	for(int i=0; i!=4; ++i) {
26baa5fb	67	fMass[i] = NULL;
	68	fDCA[i] = NULL;
	69	fDL[i] = NULL;
	70	fCTP[i] = NULL;
	71	fd0d0[i] = NULL;
	72	fPhi[i] = NULL;
	73	fEta[i] = NULL;
	74	fPt[i] = NULL;
	75	fAPhi[i] = NULL;
	76	fAEta[i] = NULL;
	77	fAPt[i] = NULL;
	78	fBPhi[i] = NULL;
	79	fBEta[i] = NULL;
	80	fBPt[i] = NULL;
41dc4195	81	}
6e214c87	82	}
	83
	84	//_____________________________________________________________________________
	85	AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates(const char name, AliFlowEventCuts cutsEvent, AliESDtrackCuts *cuts, Double_t MassMin, Double_t MassMax) :
	86	AliAnalysisTaskSE(name),
	87	fCutsEvent(cutsEvent),
	88	fCuts(cuts),
41dc4195	89	fQAList(NULL),
6e214c87	90	fMassMin(MassMin),
41dc4195	91	fMassMax(MassMax),
41dc4195	92	fDLcut(0),
26baa5fb	93	fEvent(NULL),
26baa5fb	94	fMulti(NULL)
6e214c87	95	{
26baa5fb	96	//
41dc4195	97	for(int i=0; i!=4; ++i) {
26baa5fb	98	fMass[i] = NULL;
	99	fDCA[i] = NULL;
	100	fDL[i] = NULL;
	101	fCTP[i] = NULL;
	102	fd0d0[i] = NULL;
	103	fPhi[i] = NULL;
	104	fEta[i] = NULL;
	105	fPt[i] = NULL;
	106	fAPhi[i] = NULL;
	107	fAEta[i] = NULL;
	108	fAPt[i] = NULL;
	109	fBPhi[i] = NULL;
	110	fBEta[i] = NULL;
	111	fBPt[i] = NULL;
41dc4195	112	}
6e214c87	113	DefineInput( 0, TChain::Class() );
	114	DefineOutput( 1, TObjArray::Class() );
	115	DefineOutput( 2, TList::Class() );
	116	}
	117
	118	//_____________________________________________________________________________
	119	AliAnalysisTaskFlowK0Candidates::~AliAnalysisTaskFlowK0Candidates()
	120	{
26baa5fb	121	//
41dc4195	122	if(fQAList) delete fQAList;
6e214c87	123	}
	124
	125	//_____________________________________________________________________________
	126	void AliAnalysisTaskFlowK0Candidates::UserCreateOutputObjects()
	127	{
26baa5fb	128	//
41dc4195	129	fQAList = new TList();
	130	fQAList->SetOwner();
	131	AddQAEvents();
	132	AddQACandidates();
	133	PostData( 2, fQAList);
	134	}
	135	//_____________________________________________________________________________
	136	void AliAnalysisTaskFlowK0Candidates::AddQAEvents()
	137	{
26baa5fb	138	//
41dc4195	139	TList *tQAEvents = new TList();
	140	tQAEvents->SetName("Events");
	141	tQAEvents->SetOwner();
26baa5fb	142	fEvent = new TH1D("Event", "Number of Events", 2, 0, 2);
	143	tQAEvents->Add(fEvent);
	144	fMulti = new TH1D("Multiplicity", "Multiplicity", 180, 0, 10000);
	145	tQAEvents->Add(fMulti);
41dc4195	146	fQAList->Add(tQAEvents);
	147	}
	148	//_____________________________________________________________________________
	149	void AliAnalysisTaskFlowK0Candidates::AddQACandidates()
	150	{
26baa5fb	151	//
41dc4195	152	TList *tQACandidates[4];
	153	TList *tQADaughters[4];
	154	for(int i=0; i!=4; ++i) {
	155	tQACandidates[i] = new TList();
	156	tQACandidates[i]->SetOwner();
	157	tQACandidates[i]->SetName(Form("Candidates%d",i));
26baa5fb	158	fMass[i] = new TH1D( Form("Mass%i",i), "Mass;M_{#pi#pi} [GeV];Counts per MeV", 180, 0.41, 0.59); tQACandidates[i]->Add( fMass[i] );
	159	fDCA[i] = new TH1D( Form("DCA%i" ,i), "DCA;[cm];Counts per 10 um", 180, 0.00, 0.18); tQACandidates[i]->Add( fDCA[i] );
	160	fDL[i] = new TH1D( Form("DL%i" ,i), "Decay Length;[cm];Counts per 0.1 mm", 180, 0.00, 1.80); tQACandidates[i]->Add( fDL[i] );
	161	fCTP[i] = new TH1D( Form("CTP%i" ,i), "Cos#theta_{p}", 180,-1.10, 1.10); tQACandidates[i]->Add( fCTP[i] );
	162	fd0d0[i] = new TH1D( Form("d0d0%i",i), "d_{0}xd_{0};[cm^{2}];Cnts 0.01 mm^{2}",180,-0.009,0.009);tQACandidates[i]->Add( fd0d0[i] );
	163	fPhi[i] = new TH1D( Form("Phi%i" ,i), "Phi;[rad];Counts per degree", 180,0,TMath::TwoPi()); tQACandidates[i]->Add( fPhi[i] );
	164	fEta[i] = new TH1D( Form("Eta%i" ,i), "Eta;;Counts per 0.04", 180,-3.60, 3.60); tQACandidates[i]->Add( fEta[i] );
	165	fPt[i] = new TH1D( Form("Pt%i" ,i), "Pt;[GeV];Counts per 0.1 GeV", 180, 0.00,18.00); tQACandidates[i]->Add( fPt[i] );
41dc4195	166	tQADaughters[i] = new TList();
	167	tQADaughters[i]->SetOwner();
	168	tQADaughters[i]->SetName(Form("Daughters%d",i));
26baa5fb	169	fAPhi[i] = new TH1D( Form("PhiBef%i",i), "Phi prePropagation;[rad];Cnts per degree",180,0,TMath::TwoPi()); tQADaughters[i]->Add( fAPhi[i] );
	170	fAEta[i] = new TH1D( Form("EtaBef%i",i), "Eta prePropagation;;Counts per 0.04" ,180,-3.6,3.6); tQADaughters[i]->Add( fAEta[i] );
	171	fAPt[i] = new TH1D( Form("PtBef%i" ,i), "Pt prePropagation;[GeV];Counts per 0.1 GeV",180, 0, 18); tQADaughters[i]->Add( fAPt[i] );
	172	fBPhi[i] = new TH1D( Form("PhiAft%i",i), "Phi posPropagation;[rad];Cnts per degree",180,0,TMath::TwoPi()); tQADaughters[i]->Add( fBPhi[i] );
	173	fBPhi[i] = new TH1D( Form("EtaAft%i",i), "Eta posPropagation;;Counts per 0.04" ,180,-3.6,3.6); tQADaughters[i]->Add( fBPhi[i] );
	174	fBPt[i] = new TH1D( Form("PtAft%i" ,i), "Pt posPropagation;[GeV];Counts per 0.1 GeV",180, 0, 18); tQADaughters[i]->Add( fBPt[i] );
41dc4195	175	tQACandidates[i]->Add(tQADaughters[i]);
	176	fQAList->Add(tQACandidates[i]);
	177	}
6e214c87	178	}
6e214c87	179	//_____________________________________________________________________________
	180	void AliAnalysisTaskFlowK0Candidates::NotifyRun()
	181	{
26baa5fb	182	//
41dc4195	183	// AliESDEvent fESD = dynamic_cast<AliESDEvent>(InputEvent());
41dc4195	184	// if(!fESD) return;
6e214c87	185	}
	186
	187	//_____________________________________________________________________________
	188	void AliAnalysisTaskFlowK0Candidates::UserExec(Option_t *)
	189	{
26baa5fb	190	//
41dc4195	191	AliESDEvent fESD = dynamic_cast<AliESDEvent>(InputEvent());
	192	AliAODEvent fAOD = dynamic_cast<AliAODEvent>(InputEvent());
	193	if( (!fESD)&&(!fAOD) )
	194	return;
	195	// printf("\nEvent found");
26baa5fb	196	fEvent->Fill( 0 );
6e214c87	197	if(fCutsEvent)
41dc4195	198	if( !(fCutsEvent->IsSelected(InputEvent())) )
41dc4195	199	return;
26baa5fb	200	fEvent->Fill( 1 );
41dc4195	201	if(fESD)
	202	ReadFromESD(fESD);
	203	else if(fAOD)
	204	ReadFromAOD(fAOD);
	205	}
	206	//_____________________________________________________________________________
26baa5fb	207	void AliAnalysisTaskFlowK0Candidates::ReadFromESD(const AliESDEvent *fESD)
41dc4195	208	{
26baa5fb	209	//
	210	TObjArray *gPOIselection = new TObjArray();
	211	gPOIselection->SetOwner();
6e214c87	212	int nTracks = fESD->GetNumberOfTracks();
26baa5fb	213	fMulti->Fill( nTracks );
41dc4195	214	for(int i=0; i!=nTracks; ++i) {
6e214c87	215	AliESDtrack *ioT = fESD->GetTrack(i);
	216	if(fCuts)
	217	if( !(fCuts->IsSelected(ioT)) )
	218	continue;
41dc4195	219	// printf("\n first particle OK...");
41dc4195	220	for(int j=i+1; j!=nTracks; ++j) {
6e214c87	221	AliESDtrack *joT = fESD->GetTrack(j);
41dc4195	222	if( (ioT->Charge()*joT->Charge()) > 0 )
6e214c87	223	continue;
	224	if(fCuts)
	225	if( !(fCuts->IsSelected(joT)) )
	226	continue;
41dc4195	227	// printf("\n second also...");
6e214c87	228	AliESDtrack iT = new AliESDtrack(ioT);
	229	AliESDtrack jT = new AliESDtrack(joT);
	230	// getting distance of closest approach
	231	double DCA = iT->PropagateToDCA(jT,fESD->GetMagneticField());
41dc4195	232	// printf(" propagated...");
6e214c87	233	// getting decay length
41dc4195	234	TVector3 vp, vv, vl;
	235	vp = TVector3( (iT->Xv()+jT->Xv())/2, (iT->Yv()+jT->Yv())/2, (iT->Zv()+jT->Zv())/2 ); // candidate position
	236	vv = TVector3( fESD->GetPrimaryVertex()->GetX(), fESD->GetPrimaryVertex()->GetY(), fESD->GetPrimaryVertex()->GetZ() ); // vertex position
	237	vl = vp - vv; // flight line
26baa5fb	238	double gDL;
26baa5fb	239	gDL = vl.Mag();
41dc4195	240	// getting cos pointing angle
26baa5fb	241	double gCTP;
41dc4195	242	TVector3 vi, vj, vs;
6e214c87	243	vi = TVector3( iT->Px(), iT->Py(), iT->Pz() );
	244	vj = TVector3( jT->Px(), jT->Py(), jT->Pz() );
	245	vs = vi + vj;
26baa5fb	246	gCTP = TMath::Cos( vl.Angle(vs) ); // Marta says: "Carlos, it is faster if you calculate this by hand!"
	247	double gD0D0;
	248	gD0D0 = iT->GetD(vv.X(),vv.Y(),fESD->GetMagneticField())*jT->GetD(vv.X(),vv.Y(),fESD->GetMagneticField());
6e214c87	249	// getting invariant mass
	250	double sum12 = iT->P()iT->P()+jT->P()jT->P();
	251	double pro12 = iT->P()iT->P()jT->P()*jT->P();
26baa5fb	252	double gInvMass = 0;
	253	gInvMass += TMath::Power(0.13957018,2);
	254	gInvMass += sqrt( TMath::Power(0.13957018,4) + TMath::Power(0.13957018,2)*(sum12) + pro12 );
	255	gInvMass -= ( iT->Px()jT->Px()+iT->Py()jT->Py()+iT->Pz()*jT->Pz() );
	256	gInvMass *= 2;
	257	gInvMass = sqrt(gInvMass);
41dc4195	258	// filtering
41dc4195	259	int iLevel = 3;
26baa5fb	260	if( gCTP<0.8 )
41dc4195	261	iLevel = 2;
26baa5fb	262	if( gDL<fDLcut ) // 0.5
41dc4195	263	iLevel = 1;
	264	if( DCA>0.05 )
	265	iLevel = 0;
	266	// printf(" candidate at level %d...",iLevel);
	267	for(int h=0; h!=iLevel+1; ++h) {
	268	// printf(" %d",h);
	269	// candidate
26baa5fb	270	fDCA[h]->Fill( DCA );
	271	fDL[h]->Fill( gDL );
	272	fCTP[h]->Fill( gCTP );
	273	fd0d0[h]->Fill( gD0D0 );
	274	fMass[h]->Fill( gInvMass );
	275	fPhi[h]->Fill( vs.Phi()+TMath::Pi() );
	276	fEta[h]->Fill( vs.Eta() );
	277	fPt[h]->Fill( vs.Pt() );
41dc4195	278	// daughters
26baa5fb	279	fAPhi[h]->Fill( iT->Phi() );
	280	fAEta[h]->Fill( iT->Eta() );
	281	fAPt[h]->Fill( iT->Pt() );
	282	fAPhi[h]->Fill( jT->Phi() );
	283	fAEta[h]->Fill( jT->Eta() );
	284	fAPt[h]->Fill( jT->Pt() );
	285	fBPhi[h]->Fill( ioT->Phi());
	286	fBPhi[h]->Fill( ioT->Eta());
	287	fBPt[h]->Fill( ioT->Pt() );
	288	fBPhi[h]->Fill( joT->Phi());
	289	fBPhi[h]->Fill( joT->Eta());
	290	fBPt[h]->Fill( joT->Pt() );
41dc4195	291	}
41dc4195	292	// building candidate
26baa5fb	293	if( (gInvMass>fMassMin) && (gInvMass<fMassMax) && (iLevel==3) ) {
41dc4195	294	AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
26baa5fb	295	sTrack->SetMass( gInvMass );
41dc4195	296	sTrack->SetPt( vs.Pt() );
	297	sTrack->SetPhi( vs.Phi()+TMath::Pi() );
	298	sTrack->SetEta( vs.Eta() );
	299	sTrack->SetCharge( 0 );
	300	sTrack->AddDaughter( iT->GetID() );
	301	sTrack->AddDaughter( jT->GetID() );
26baa5fb	302	gPOIselection->AddLast( sTrack );
6e214c87	303	}
6e214c87	304	delete iT;
6e214c87	305	delete jT;
41dc4195	306	} // endfor j
41dc4195	307	} // endfor i
26baa5fb	308	PostData( 1, gPOIselection );
41dc4195	309	PostData( 2, fQAList);
	310	}
	311
	312	//_____________________________________________________________________________
26baa5fb	313	void AliAnalysisTaskFlowK0Candidates::ReadFromAOD(const AliAODEvent *fAOD)
41dc4195	314	{
26baa5fb	315	//
	316	TObjArray *gPOIselection = new TObjArray();
	317	gPOIselection->SetOwner();
41dc4195	318	int nTracks = fAOD->GetNumberOfTracks();
26baa5fb	319	fMulti->Fill( nTracks );
41dc4195	320	for(int i=0; i!=fAOD->GetNumberOfV0s(); ++i) {
	321	AliAODv0 iV0 = (AliAODv0) fAOD->GetV0(i);
	322	if ( iV0->GetOnFlyStatus() ) continue;
	323	if ( (iV0->ChargeProng(0)*iV0->ChargeProng(1))>0 ) continue;
	324	// getting distance of closest approach
	325	double DCA = iV0->DcaV0Daughters();
	326	// getting decay length
26baa5fb	327	double gDL;
41dc4195	328	double vv[3];
	329	vv[0] = fAOD->GetPrimaryVertex()->GetX();
	330	vv[1] = fAOD->GetPrimaryVertex()->GetY();
	331	vv[2] = fAOD->GetPrimaryVertex()->GetZ();
26baa5fb	332	gDL = iV0->DecayLengthV0(vv);
41dc4195	333	// cos pointing angle
26baa5fb	334	double gCTP;
	335	gCTP = iV0->CosPointingAngle(vv);
	336	double gD0D0;
	337	gD0D0 = iV0->Prodd0d0();
41dc4195	338	// getting invariant mass
26baa5fb	339	double gInvMass = iV0->MassK0Short();
41dc4195	340	// filtering
41dc4195	341	int iLevel = 3;
26baa5fb	342	if( gCTP<0.8 )
41dc4195	343	iLevel = 2;
26baa5fb	344	if( gDL<fDLcut ) // 0.5
41dc4195	345	iLevel = 1;
	346	if( DCA>0.05 )
	347	iLevel = 0;
	348	for(int h=0; h!=iLevel+1; ++h) {
	349	// candidate
26baa5fb	350	fDCA[h]->Fill( DCA );
	351	fDL[h]->Fill( gDL );
	352	fCTP[h]->Fill( gCTP );
	353	fd0d0[h]->Fill( gD0D0 );
	354	fMass[h]->Fill( gInvMass );
	355	fPhi[h]->Fill( iV0->Phi() );
	356	fEta[h]->Fill( iV0->Eta() );
	357	fPt[h]->Fill( iV0->Pt() );
41dc4195	358	// daughters
26baa5fb	359	fAPhi[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Phi() );
	360	fAEta[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Eta() );
	361	fAPt[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Pt() );
	362	fAPhi[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Phi() );
	363	fAEta[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Eta() );
	364	fAPt[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Pt() );
	365	fBPhi[h]->Fill( iV0->PhiProng(0) );
	366	fBPhi[h]->Fill( iV0->EtaProng(0) );
	367	fBPt[h]->Fill( iV0->PtProng(0) );
	368	fBPhi[h]->Fill( iV0->PhiProng(1) );
	369	fBPhi[h]->Fill( iV0->EtaProng(1) );
	370	fBPt[h]->Fill( iV0->PtProng(1) );
6e214c87	371	}
41dc4195	372	AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
26baa5fb	373	sTrack->SetMass( gInvMass );
41dc4195	374	sTrack->SetPt( iV0->Pt() );
	375	sTrack->SetPhi( iV0->Phi() );
	376	sTrack->SetEta( iV0->Eta() );
	377	sTrack->SetCharge( 0 );
	378	sTrack->AddDaughter( iV0->GetPosID() );
	379	sTrack->AddDaughter( iV0->GetNegID() );
26baa5fb	380	gPOIselection->AddLast( sTrack );
6e214c87	381	}
26baa5fb	382	PostData( 1, gPOIselection );
41dc4195	383	PostData( 2, fQAList);
6e214c87	384	}
	385
	386	//_____________________________________________________________________________
	387	void AliAnalysisTaskFlowK0Candidates::Terminate(Option_t *)
	388	{
26baa5fb	389	//
6e214c87	390	}
6e214c87	391