[u/mrichter/AliRoot.git] / PWG0 / base / AliAnalysisTaskESDfilter.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */
 
#include <TChain.h>
#include <TFile.h>

#include "AliAnalysisTaskESDfilter.h"
#include "AliAnalysisManager.h"
#include "AliESDEvent.h"
#include "AliAODEvent.h"
#include "AliAODHandler.h"
#include "AliAnalysisFilter.h"
#include "AliESDtrack.h"
#include "AliESDMuonTrack.h"
#include "AliESDVertex.h"
#include "AliESDv0.h"
#include "AliESDkink.h"
#include "AliESDcascade.h"
#include "AliLog.h"

ClassImp(AliAnalysisTaskESDfilter)

////////////////////////////////////////////////////////////////////////

AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
    fDebug(0),
    fChain(0x0),
    fESD(0x0),
    fAOD(0x0),
    fTreeA(0x0),
    fTrackFilter(0x0),
    fKinkFilter(0x0),
    fV0Filter(0x0)
{
  // Default constructor
}

AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
    AliAnalysisTask(name, "AnalysisTaskESDfilter"),
    fDebug(0),
    fChain(0x0),
    fESD(0x0),
    fAOD(0x0),
    fTreeA(0x0),
    fTrackFilter(0x0),
    fKinkFilter(0x0),
    fV0Filter(0x0)
{
  // Default constructor
    DefineInput (0, TChain::Class());
    DefineOutput(0, TTree::Class());
}

void AliAnalysisTaskESDfilter::CreateOutputObjects()
{
// Create the output container
    if (fDebug > 1) AliInfo("CreateOutPutData() \n");
    AliAODHandler* handler = (AliAODHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
    
    fAOD   = handler->GetAOD();
    fTreeA = handler->GetTree();
    fTreeA->GetUserInfo()->Add(fTrackFilter);
}

void AliAnalysisTaskESDfilter::Init()
{
    // Initialization
    if (fDebug > 1) AliInfo("Init() \n");
    // Call configuration file
}

void AliAnalysisTaskESDfilter::ConnectInputData(Option_t */*option*/)
{
// Connect the input data
//
    if (fDebug > 1) AliInfo("ConnectInputData() \n");
    fChain = (TChain*)GetInputData(0);
    fESD = new AliESDEvent();
    fESD->ReadFromTree(fChain);
}

void AliAnalysisTaskESDfilter::Exec(Option_t */*option*/)
{
// Execute analysis for current event
//
    AliESD* old = fESD->GetAliESDOld();
    if (old) fESD->CopyFromOldESD();
    
    Long64_t ientry = fChain->GetReadEntry();
    AliInfo(Form("ESD Filter: Analysing event # %5d\n", (Int_t) ientry));
    
    
    // ESD Filter analysis task executed for each event
    
    // set arrays and pointers
    Double_t pos[3];
    Double_t p[3];
    Double_t covVtx[6];
    Double_t covTr[21];
    Double_t pid[10];

    for (Int_t i = 0; i < 6; i++)  covVtx[i] = 0.;
    for (Int_t i = 0; i < 21; i++) covTr [i] = 0.;

    
  // loop over events and fill them
  
  // Multiplicity information needed by the header (to be revised!)
    Int_t nTracks    = fESD->GetNumberOfTracks();
    Int_t nPosTracks = 0;
    for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack) 
	if (fESD->GetTrack(iTrack)->GetSign()> 0) nPosTracks++;
    
    // Update the header

    AliAODHeader* header = fAOD->GetHeader();
    header->SetRunNumber(fESD->GetRunNumber());
    if (old) {
	header->SetBunchCrossNumber(0);
	header->SetOrbitNumber(0);
	header->SetPeriodNumber(0);
	header->SetEventType(0);
	header->SetMuonMagFieldScale(-999.);
	header->SetCentrality(-999.);
    } else {
	header->SetBunchCrossNumber(fESD->GetBunchCrossNumber());
	header->SetOrbitNumber(fESD->GetOrbitNumber());
	header->SetPeriodNumber(fESD->GetPeriodNumber());
	header->SetEventType(fESD->GetEventType());
	header->SetMuonMagFieldScale(-999.);
	header->SetCentrality(-999.);
    }
    
    header->SetTriggerMask(fESD->GetTriggerMask()); 
    header->SetTriggerCluster(fESD->GetTriggerCluster());
    header->SetMagneticField(fESD->GetMagneticField());
    header->SetZDCN1Energy(fESD->GetZDCN1Energy());
    header->SetZDCP1Energy(fESD->GetZDCP1Energy());
    header->SetZDCN2Energy(fESD->GetZDCN2Energy());
    header->SetZDCP2Energy(fESD->GetZDCP2Energy());
    header->SetZDCEMEnergy(fESD->GetZDCEMEnergy(0),fESD->GetZDCEMEnergy(1));
    header->SetRefMultiplicity(nTracks);
    header->SetRefMultiplicityPos(nPosTracks);
    header->SetRefMultiplicityNeg(nTracks - nPosTracks);
//
//    
    Int_t nV0s      = fESD->GetNumberOfV0s();
    Int_t nCascades = fESD->GetNumberOfCascades();
    Int_t nKinks    = fESD->GetNumberOfKinks();
    
    printf("   NV0=%d  NCASCADES=%d  NKINKS=%d\n", nV0s, nCascades, nKinks);

    Int_t nVertices = nV0s + nCascades + nKinks;
    
    fAOD->ResetStd(nTracks, nVertices);
    AliAODTrack *aodTrack;
    
    
    // Array to take into account the tracks already added to the AOD
    Bool_t * usedTrack = NULL;
    if (nTracks>0) {
	usedTrack = new Bool_t[nTracks];
	for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
    }
    // Array to take into account the V0s already added to the AOD
    Bool_t * usedV0 = NULL;
    if (nV0s>0) {
	usedV0 = new Bool_t[nV0s];
	for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
    }
    // Array to take into account the kinks already added to the AOD
    Bool_t * usedKink = NULL;
    if (nKinks>0) {
	usedKink = new Bool_t[nKinks];
	for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
    }
    
    // Access to the AOD container of vertices
    TClonesArray &vertices = *(fAOD->GetVertices());
    Int_t jVertices=0;
    
    // Access to the AOD container of tracks
    TClonesArray &tracks = *(fAOD->GetTracks());
    Int_t jTracks=0; 
    
    // Add primary vertex. The primary tracks will be defined
    // after the loops on the composite objects (V0, cascades, kinks)
    const AliESDVertex *vtx = fESD->GetPrimaryVertex();
    
    vtx->GetXYZ(pos); // position
    vtx->GetCovMatrix(covVtx); //covariance matrix
    
    AliAODVertex * primary = new(vertices[jVertices++])
	AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, AliAODVertex::kPrimary);
    primary->Print();


    // Create vertices starting from the most complex objects
    Double_t chi2 = 0.;
    
    // Cascades
    for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
	AliESDcascade *cascade = fESD->GetCascade(nCascade);
	
	cascade->GetXYZ(pos[0], pos[1], pos[2]);

	if (!old) {
	    chi2 = vtx->GetChi2toNDF();
	    cascade->GetPosCovXi(covVtx);
	} else {
	    chi2 = -999.;
	}
	// Add the cascade vertex
	AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
									  covVtx,
									  chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
									  primary,
									  AliAODVertex::kCascade);
	
	primary->AddDaughter(vcascade);
	
	// Add the V0 from the cascade. The ESD class have to be optimized...
	// Now we have to search for the corresponding Vo in the list of V0s
	// using the indeces of the positive and negative tracks
	
	Int_t posFromV0 = cascade->GetPindex();
	Int_t negFromV0 = cascade->GetNindex();
	
      
	AliESDv0 * v0 = 0x0;
	Int_t indV0 = -1;
	
	for (Int_t iV0=0; iV0<nV0s; ++iV0) {
	    
	    v0 = fESD->GetV0(iV0);
	    Int_t posV0 = v0->GetPindex();
	    Int_t negV0 = v0->GetNindex();
	    
	    if (posV0==posFromV0 && negV0==negFromV0) {
		indV0 = iV0;
		break;
	    }
	}
	
	AliAODVertex * vV0FromCascade = 0x0;
	
	if (indV0>-1 && !usedV0[indV0] ) {
	    
	    // the V0 exists in the array of V0s and is not used
	    
	    usedV0[indV0] = kTRUE;
	    
	    v0->GetXYZ(pos[0], pos[1], pos[2]);
	    if (!old) {
		chi2 = v0->GetChi2V0();
		v0->GetPosCov(covVtx);
	    } else {
		chi2 = -999.;
	    }

	    vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
								     covVtx,
								     chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
								     vcascade,
								     AliAODVertex::kV0);
	} else {
	    
	    // the V0 doesn't exist in the array of V0s or was used
//	    cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
//		 << " The V0 " << indV0 
//		 << " doesn't exist in the array of V0s or was used!" << endl;
	    
	    cascade->GetXYZ(pos[0], pos[1], pos[2]);
	    
	    if (!old) {
		chi2 = v0->GetChi2V0();
		cascade->GetPosCov(covVtx);
	    } else {
		chi2 = -999.;
	    }

	    vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
								     covVtx,
								     chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
								     vcascade,
								     AliAODVertex::kV0);
	    vcascade->AddDaughter(vV0FromCascade);
	}
	
	// Add the positive tracks from the V0
	
	if (posFromV0>-1 && !usedTrack[posFromV0]) {
	    
	    usedTrack[posFromV0] = kTRUE;
	    
	    AliESDtrack *esdTrack = fESD->GetTrack(posFromV0);
	    esdTrack->GetPxPyPz(p);
	    esdTrack->GetXYZ(pos);
	    esdTrack->GetCovarianceXYZPxPyPz(covTr);
	    esdTrack->GetESDpid(pid);
	    
	    vV0FromCascade->AddDaughter(aodTrack =
					new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
									   esdTrack->GetLabel(), 
									   p, 
									   kTRUE,
									   pos,
									   kFALSE,
									   covTr, 
									   (Short_t)esdTrack->GetSign(),
									   esdTrack->GetITSClusterMap(), 
									   pid,
									   vV0FromCascade,
									   kTRUE,  // check if this is right
									   kFALSE, // check if this is right
									   AliAODTrack::kSecondary)
		);
	    aodTrack->ConvertAliPIDtoAODPID();
	}
	else {
//	    cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
//		 << " track " << posFromV0 << " has already been used!" << endl;
	}
	
	// Add the negative tracks from the V0
	
	if (negFromV0>-1 && !usedTrack[negFromV0]) {
	    
	    usedTrack[negFromV0] = kTRUE;
	    
	    AliESDtrack *esdTrack = fESD->GetTrack(negFromV0);
	    esdTrack->GetPxPyPz(p);
	    esdTrack->GetXYZ(pos);
	    esdTrack->GetCovarianceXYZPxPyPz(covTr);
	    esdTrack->GetESDpid(pid);
	    
	    vV0FromCascade->AddDaughter(aodTrack =
					new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
									   esdTrack->GetLabel(),
									   p,
									   kTRUE,
									   pos,
									   kFALSE,
									   covTr, 
									   (Short_t)esdTrack->GetSign(),
									   esdTrack->GetITSClusterMap(), 
									   pid,
									   vV0FromCascade,
									   kTRUE,  // check if this is right
									   kFALSE, // check if this is right
									   AliAODTrack::kSecondary)
		);
	    aodTrack->ConvertAliPIDtoAODPID();
	}
	else {
//	    cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
//		 << " track " << negFromV0 << " has already been used!" << endl;
	}
	
	// Add the bachelor track from the cascade
	
	Int_t bachelor = cascade->GetBindex();
	
	if(bachelor>-1 && !usedTrack[bachelor]) {
	    
	    usedTrack[bachelor] = kTRUE;
	    
	    AliESDtrack *esdTrack = fESD->GetTrack(bachelor);
	    esdTrack->GetPxPyPz(p);
	    esdTrack->GetXYZ(pos);
	    esdTrack->GetCovarianceXYZPxPyPz(covTr);
	    esdTrack->GetESDpid(pid);
	    
	    vcascade->AddDaughter(aodTrack =
				  new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
								     esdTrack->GetLabel(),
								     p,
								     kTRUE,
								     pos,
								     kFALSE,
								     covTr, 
								     (Short_t)esdTrack->GetSign(),
								     esdTrack->GetITSClusterMap(), 
								     pid,
								     vcascade,
								     kTRUE,  // check if this is right
								     kFALSE, // check if this is right
								     AliAODTrack::kSecondary)
		);
	    aodTrack->ConvertAliPIDtoAODPID();
	}
	else {
//	    cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
//		 << " track " << bachelor << " has already been used!" << endl;
	}
	
	// Add the primary track of the cascade (if any)
	
    } // end of the loop on cascades
   
    // V0s
    
    for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
	
	if (usedV0[nV0]) continue; // skip if aready added to the AOD
	
	AliESDv0 *v0 = fESD->GetV0(nV0);
	
	v0->GetXYZ(pos[0], pos[1], pos[2]);

	if (!old) {
	    chi2 = v0->GetChi2V0();
	    v0->GetPosCov(covVtx);
	} else {
	    chi2 = -999.;
	}


	AliAODVertex * vV0 = 
	    new(vertices[jVertices++]) AliAODVertex(pos,
						    covVtx,
						    chi2, // = chi2/NDF since NDF = 2*2-3
						    primary,
						    AliAODVertex::kV0);
	primary->AddDaughter(vV0);
	
	Int_t posFromV0 = v0->GetPindex();
	Int_t negFromV0 = v0->GetNindex();
	
	// Add the positive tracks from the V0
	
	if (posFromV0>-1 && !usedTrack[posFromV0]) {
	    
	    usedTrack[posFromV0] = kTRUE;
	    
	    AliESDtrack *esdTrack = fESD->GetTrack(posFromV0);
	    esdTrack->GetPxPyPz(p);
	    esdTrack->GetXYZ(pos);
	    esdTrack->GetCovarianceXYZPxPyPz(covTr);
	    esdTrack->GetESDpid(pid);
	    
	    vV0->AddDaughter(aodTrack =
			     new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
								esdTrack->GetLabel(), 
								p, 
								kTRUE,
								pos,
								kFALSE,
								covTr, 
								(Short_t)esdTrack->GetSign(),
								esdTrack->GetITSClusterMap(), 
								pid,
								vV0,
								kTRUE,  // check if this is right
								kFALSE, // check if this is right
								AliAODTrack::kSecondary)
		);
	    aodTrack->ConvertAliPIDtoAODPID();
	}
	else {
//	    cerr << "Error: event " << fESD->GetEventNumberInFile() << " V0 " << nV0
//		 << " track " << posFromV0 << " has already been used!" << endl;
	}
	
	// Add the negative tracks from the V0
	
	if (negFromV0>-1 && !usedTrack[negFromV0]) {
	    
	    usedTrack[negFromV0] = kTRUE;
	    
	    AliESDtrack *esdTrack = fESD->GetTrack(negFromV0);
	    esdTrack->GetPxPyPz(p);
	    esdTrack->GetXYZ(pos);
	    esdTrack->GetCovarianceXYZPxPyPz(covTr);
	    esdTrack->GetESDpid(pid);
	    
	    vV0->AddDaughter(aodTrack =
			     new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
								esdTrack->GetLabel(),
								p,
								kTRUE,
								pos,
								kFALSE,
								covTr, 
								(Short_t)esdTrack->GetSign(),
								esdTrack->GetITSClusterMap(), 
								pid,
								vV0,
								kTRUE,  // check if this is right
								kFALSE, // check if this is right
								AliAODTrack::kSecondary)
		);
	    aodTrack->ConvertAliPIDtoAODPID();
	}
	else {
//	    cerr << "Error: event " << fESD->GetEventNumberInFile() << " V0 " << nV0
//		 << " track " << negFromV0 << " has already been used!" << endl;
	}
	
    } // end of the loop on V0s
    
    // Kinks: it is a big mess the access to the information in the kinks
    // The loop is on the tracks in order to find the mother and daugther of each kink
    
    
    for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
	
	
	AliESDtrack * esdTrack = fESD->GetTrack(iTrack);
	
	Int_t ikink = esdTrack->GetKinkIndex(0);
	
	if (ikink && nKinks) {
	    // Negative kink index: mother, positive: daughter
	    
	    // Search for the second track of the kink
	    
	    for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
		
		AliESDtrack * esdTrack1 = fESD->GetTrack(jTrack);
		
		Int_t jkink = esdTrack1->GetKinkIndex(0);
		
		if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
		    
		    // The two tracks are from the same kink
		    
		    if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
		    
		    Int_t imother = -1;
		    Int_t idaughter = -1;
		    
		    if (ikink<0 && jkink>0) {
			
			imother = iTrack;
			idaughter = jTrack;
		    }
		    else if (ikink>0 && jkink<0) {
			
			imother = jTrack;
			idaughter = iTrack;
		    }
		    else {
//			cerr << "Error: Wrong combination of kink indexes: "
//			     << ikink << " " << jkink << endl;
			continue;
		    }
		    
		    // Add the mother track
		    
		    AliAODTrack * mother = NULL;
		    
		    if (!usedTrack[imother]) {
			
			usedTrack[imother] = kTRUE;
			
			AliESDtrack *esdTrack = fESD->GetTrack(imother);
			esdTrack->GetPxPyPz(p);
			esdTrack->GetXYZ(pos);
			esdTrack->GetCovarianceXYZPxPyPz(covTr);
			esdTrack->GetESDpid(pid);
			
			mother = 
			    new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
							       esdTrack->GetLabel(),
							       p,
							       kTRUE,
							       pos,
							       kFALSE,
							       covTr, 
							       (Short_t)esdTrack->GetSign(),
							       esdTrack->GetITSClusterMap(), 
							       pid,
							       primary,
							       kTRUE, // check if this is right
							       kTRUE, // check if this is right
							       AliAODTrack::kPrimary);
			primary->AddDaughter(mother);
			mother->ConvertAliPIDtoAODPID();
		    }
		    else {
//			cerr << "Error: event " << fESD->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
//			     << " track " << imother << " has already been used!" << endl;
		  }
		    
		    // Add the kink vertex
		    AliESDkink * kink = fESD->GetKink(TMath::Abs(ikink)-1);
		    
		    AliAODVertex * vkink = 
			new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
								NULL,
								0.,
								mother,
								AliAODVertex::kKink);
		    // Add the daughter track
		  
		    AliAODTrack * daughter = NULL;
		    
		    if (!usedTrack[idaughter]) {
			
			usedTrack[idaughter] = kTRUE;
			
			AliESDtrack *esdTrack = fESD->GetTrack(idaughter);
			esdTrack->GetPxPyPz(p);
			esdTrack->GetXYZ(pos);
			esdTrack->GetCovarianceXYZPxPyPz(covTr);
			esdTrack->GetESDpid(pid);
			
			daughter = 
			    new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
							       esdTrack->GetLabel(),
							       p,
							       kTRUE,
							       pos,
							       kFALSE,
							       covTr, 
							       (Short_t)esdTrack->GetSign(),
							       esdTrack->GetITSClusterMap(), 
							       pid,
							       vkink,
							       kTRUE, // check if this is right
							       kTRUE, // check if this is right
							       AliAODTrack::kPrimary);
			vkink->AddDaughter(daughter);
			daughter->ConvertAliPIDtoAODPID();
		    }
		    else {
//			cerr << "Error: event " << fESD->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
//			     << " track " << idaughter << " has already been used!" << endl;
		    }
		}
	    }
	}      
    }
    
  
    // Tracks (primary and orphan)

    printf("NUMBER OF TRACKS %5d\n", nTracks);
    
    for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
	
	
	if (usedTrack[nTrack]) continue;
	AliESDtrack *esdTrack = fESD->GetTrack(nTrack);
	UInt_t selectInfo = 0;
	//
	// Track selection
	if (fTrackFilter) {
	    selectInfo = fTrackFilter->IsSelected(esdTrack);
	    if (!selectInfo) continue;
	}
	
	//
	esdTrack->GetPxPyPz(p);
	esdTrack->GetXYZ(pos);
	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	esdTrack->GetESDpid(pid);
	
	Float_t impactXY, impactZ;
	
	esdTrack->GetImpactParameters(impactXY,impactZ);
	
	if (impactXY<3) {
	    // track inside the beam pipe
	    
	    primary->AddDaughter(aodTrack =
				 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
								    esdTrack->GetLabel(),
								    p,
								    kTRUE,
								    pos,
								    kFALSE,
								    covTr, 
								    (Short_t)esdTrack->GetSign(),
								    esdTrack->GetITSClusterMap(), 
								    pid,
								    primary,
								    kTRUE, // check if this is right
								    kTRUE, // check if this is right
								    AliAODTrack::kPrimary, 
								    selectInfo)
		);
	    aodTrack->ConvertAliPIDtoAODPID();
	}
	else {
	    // outside the beam pipe: orphan track
	    aodTrack =
		new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
						   esdTrack->GetLabel(),
						   p,
						   kTRUE,
						   pos,
						   kFALSE,
						   covTr, 
						   (Short_t)esdTrack->GetSign(),
						   esdTrack->GetITSClusterMap(), 
						   pid,
						   NULL,
						   kFALSE, // check if this is right
						   kFALSE, // check if this is right
						   AliAODTrack::kOrphan,
						   selectInfo);
	    aodTrack->ConvertAliPIDtoAODPID();
	}	
    } // end of loop on tracks
    
    // muon tracks
    Int_t nMuTracks = fESD->GetNumberOfMuonTracks();
    for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
	
	AliESDMuonTrack *esdMuTrack = fESD->GetMuonTrack(nMuTrack);     
	p[0] = esdMuTrack->Px(); 
	p[1] = esdMuTrack->Py(); 
	p[2] = esdMuTrack->Pz();
	pos[0] = primary->GetX(); 
	pos[1] = primary->GetY(); 
	pos[2] = primary->GetZ();
	
	// has to be changed once the muon pid is provided by the ESD
	for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
	
	primary->AddDaughter(
	    new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
					       0, // no label provided
					       p,
					       kTRUE,
					       pos,
					       kFALSE,
					       NULL, // no covariance matrix provided
					       (Short_t)-99, // no charge provided
					       0, // no ITSClusterMap
					       pid,
					       primary,
					       kTRUE,  // check if this is right
					       kTRUE,  // not used for vertex fit
					       AliAODTrack::kPrimary)
	    );
    }
  
    delete [] usedTrack;
    delete [] usedV0;
    delete [] usedKink;


    //
    PostData(0, fTreeA);
    return;
}

void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
{
// Terminate analysis
//
    if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
}
Commit	Line	Data
	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, 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	/* $Id$ */
	17
	18	#include <TChain.h>
	19	#include <TFile.h>
	20
	21	#include "AliAnalysisTaskESDfilter.h"
	22	#include "AliAnalysisManager.h"
	23	#include "AliESDEvent.h"
	24	#include "AliAODEvent.h"
	25	#include "AliAODHandler.h"
	26	#include "AliAnalysisFilter.h"
	27	#include "AliESDtrack.h"
	28	#include "AliESDMuonTrack.h"
	29	#include "AliESDVertex.h"
	30	#include "AliESDv0.h"
	31	#include "AliESDkink.h"
	32	#include "AliESDcascade.h"
	33	#include "AliLog.h"
	34
	35	ClassImp(AliAnalysisTaskESDfilter)
	36
	37	////////////////////////////////////////////////////////////////////////
	38
	39	AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
	40	fDebug(0),
	41	fChain(0x0),
	42	fESD(0x0),
	43	fAOD(0x0),
	44	fTreeA(0x0),
	45	fTrackFilter(0x0),
	46	fKinkFilter(0x0),
	47	fV0Filter(0x0)
	48	{
	49	// Default constructor
	50	}
	51
	52	AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
	53	AliAnalysisTask(name, "AnalysisTaskESDfilter"),
	54	fDebug(0),
	55	fChain(0x0),
	56	fESD(0x0),
	57	fAOD(0x0),
	58	fTreeA(0x0),
	59	fTrackFilter(0x0),
	60	fKinkFilter(0x0),
	61	fV0Filter(0x0)
	62	{
	63	// Default constructor
	64	DefineInput (0, TChain::Class());
	65	DefineOutput(0, TTree::Class());
	66	}
	67
	68	void AliAnalysisTaskESDfilter::CreateOutputObjects()
	69	{
	70	// Create the output container
	71	if (fDebug > 1) AliInfo("CreateOutPutData() \n");
	72	AliAODHandler* handler = (AliAODHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
	73
	74	fAOD = handler->GetAOD();
	75	fTreeA = handler->GetTree();
	76	fTreeA->GetUserInfo()->Add(fTrackFilter);
	77	}
	78
	79	void AliAnalysisTaskESDfilter::Init()
	80	{
	81	// Initialization
	82	if (fDebug > 1) AliInfo("Init() \n");
	83	// Call configuration file
	84	}
	85
	86	void AliAnalysisTaskESDfilter::ConnectInputData(Option_t /option*/)
	87	{
	88	// Connect the input data
	89	//
	90	if (fDebug > 1) AliInfo("ConnectInputData() \n");
	91	fChain = (TChain*)GetInputData(0);
	92	fESD = new AliESDEvent();
	93	fESD->ReadFromTree(fChain);
	94	}
	95
	96	void AliAnalysisTaskESDfilter::Exec(Option_t /option*/)
	97	{
	98	// Execute analysis for current event
	99	//
	100	AliESD* old = fESD->GetAliESDOld();
	101	if (old) fESD->CopyFromOldESD();
	102
	103	Long64_t ientry = fChain->GetReadEntry();
	104	AliInfo(Form("ESD Filter: Analysing event # %5d\n", (Int_t) ientry));
	105
	106
	107	// ESD Filter analysis task executed for each event
	108
	109	// set arrays and pointers
	110	Double_t pos[3];
	111	Double_t p[3];
	112	Double_t covVtx[6];
	113	Double_t covTr[21];
	114	Double_t pid[10];
	115
	116	for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
	117	for (Int_t i = 0; i < 21; i++) covTr [i] = 0.;
	118
	119
	120	// loop over events and fill them
	121
	122	// Multiplicity information needed by the header (to be revised!)
	123	Int_t nTracks = fESD->GetNumberOfTracks();
	124	Int_t nPosTracks = 0;
	125	for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack)
	126	if (fESD->GetTrack(iTrack)->GetSign()> 0) nPosTracks++;
	127
	128	// Update the header
	129
	130	AliAODHeader* header = fAOD->GetHeader();
	131	header->SetRunNumber(fESD->GetRunNumber());
	132	if (old) {
	133	header->SetBunchCrossNumber(0);
	134	header->SetOrbitNumber(0);
	135	header->SetPeriodNumber(0);
	136	header->SetEventType(0);
	137	header->SetMuonMagFieldScale(-999.);
	138	header->SetCentrality(-999.);
	139	} else {
	140	header->SetBunchCrossNumber(fESD->GetBunchCrossNumber());
	141	header->SetOrbitNumber(fESD->GetOrbitNumber());
	142	header->SetPeriodNumber(fESD->GetPeriodNumber());
	143	header->SetEventType(fESD->GetEventType());
	144	header->SetMuonMagFieldScale(-999.);
	145	header->SetCentrality(-999.);
	146	}
	147
	148	header->SetTriggerMask(fESD->GetTriggerMask());
	149	header->SetTriggerCluster(fESD->GetTriggerCluster());
	150	header->SetMagneticField(fESD->GetMagneticField());
	151	header->SetZDCN1Energy(fESD->GetZDCN1Energy());
	152	header->SetZDCP1Energy(fESD->GetZDCP1Energy());
	153	header->SetZDCN2Energy(fESD->GetZDCN2Energy());
	154	header->SetZDCP2Energy(fESD->GetZDCP2Energy());
	155	header->SetZDCEMEnergy(fESD->GetZDCEMEnergy(0),fESD->GetZDCEMEnergy(1));
	156	header->SetRefMultiplicity(nTracks);
	157	header->SetRefMultiplicityPos(nPosTracks);
	158	header->SetRefMultiplicityNeg(nTracks - nPosTracks);
	159	//
	160	//
	161	Int_t nV0s = fESD->GetNumberOfV0s();
	162	Int_t nCascades = fESD->GetNumberOfCascades();
	163	Int_t nKinks = fESD->GetNumberOfKinks();
	164
	165	printf(" NV0=%d NCASCADES=%d NKINKS=%d\n", nV0s, nCascades, nKinks);
	166
	167	Int_t nVertices = nV0s + nCascades + nKinks;
	168
	169	fAOD->ResetStd(nTracks, nVertices);
	170	AliAODTrack *aodTrack;
	171
	172
	173	// Array to take into account the tracks already added to the AOD
	174	Bool_t * usedTrack = NULL;
	175	if (nTracks>0) {
	176	usedTrack = new Bool_t[nTracks];
	177	for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
	178	}
	179	// Array to take into account the V0s already added to the AOD
	180	Bool_t * usedV0 = NULL;
	181	if (nV0s>0) {
	182	usedV0 = new Bool_t[nV0s];
	183	for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
	184	}
	185	// Array to take into account the kinks already added to the AOD
	186	Bool_t * usedKink = NULL;
	187	if (nKinks>0) {
	188	usedKink = new Bool_t[nKinks];
	189	for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
	190	}
	191
	192	// Access to the AOD container of vertices
	193	TClonesArray &vertices = *(fAOD->GetVertices());
	194	Int_t jVertices=0;
	195
	196	// Access to the AOD container of tracks
	197	TClonesArray &tracks = *(fAOD->GetTracks());
	198	Int_t jTracks=0;
	199
	200	// Add primary vertex. The primary tracks will be defined
	201	// after the loops on the composite objects (V0, cascades, kinks)
	202	const AliESDVertex *vtx = fESD->GetPrimaryVertex();
	203
	204	vtx->GetXYZ(pos); // position
	205	vtx->GetCovMatrix(covVtx); //covariance matrix
	206
	207	AliAODVertex * primary = new(vertices[jVertices++])
	208	AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, AliAODVertex::kPrimary);
	209	primary->Print();
	210
	211
	212	// Create vertices starting from the most complex objects
	213	Double_t chi2 = 0.;
	214
	215	// Cascades
	216	for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
	217	AliESDcascade *cascade = fESD->GetCascade(nCascade);
	218
	219	cascade->GetXYZ(pos[0], pos[1], pos[2]);
	220
	221	if (!old) {
	222	chi2 = vtx->GetChi2toNDF();
	223	cascade->GetPosCovXi(covVtx);
	224	} else {
	225	chi2 = -999.;
	226	}
	227	// Add the cascade vertex
	228	AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
	229	covVtx,
	230	chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
	231	primary,
	232	AliAODVertex::kCascade);
	233
	234	primary->AddDaughter(vcascade);
	235
	236	// Add the V0 from the cascade. The ESD class have to be optimized...
	237	// Now we have to search for the corresponding Vo in the list of V0s
	238	// using the indeces of the positive and negative tracks
	239
	240	Int_t posFromV0 = cascade->GetPindex();
	241	Int_t negFromV0 = cascade->GetNindex();
	242
	243
	244	AliESDv0 * v0 = 0x0;
	245	Int_t indV0 = -1;
	246
	247	for (Int_t iV0=0; iV0<nV0s; ++iV0) {
	248
	249	v0 = fESD->GetV0(iV0);
	250	Int_t posV0 = v0->GetPindex();
	251	Int_t negV0 = v0->GetNindex();
	252
	253	if (posV0==posFromV0 && negV0==negFromV0) {
	254	indV0 = iV0;
	255	break;
	256	}
	257	}
	258
	259	AliAODVertex * vV0FromCascade = 0x0;
	260
	261	if (indV0>-1 && !usedV0[indV0] ) {
	262
	263	// the V0 exists in the array of V0s and is not used
	264
	265	usedV0[indV0] = kTRUE;
	266
	267	v0->GetXYZ(pos[0], pos[1], pos[2]);
	268	if (!old) {
	269	chi2 = v0->GetChi2V0();
	270	v0->GetPosCov(covVtx);
	271	} else {
	272	chi2 = -999.;
	273	}
	274
	275	vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
	276	covVtx,
	277	chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
	278	vcascade,
	279	AliAODVertex::kV0);
	280	} else {
	281
	282	// the V0 doesn't exist in the array of V0s or was used
	283	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
	284	// << " The V0 " << indV0
	285	// << " doesn't exist in the array of V0s or was used!" << endl;
	286
	287	cascade->GetXYZ(pos[0], pos[1], pos[2]);
	288
	289	if (!old) {
	290	chi2 = v0->GetChi2V0();
	291	cascade->GetPosCov(covVtx);
	292	} else {
	293	chi2 = -999.;
	294	}
	295
	296	vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
	297	covVtx,
	298	chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
	299	vcascade,
	300	AliAODVertex::kV0);
	301	vcascade->AddDaughter(vV0FromCascade);
	302	}
	303
	304	// Add the positive tracks from the V0
	305
	306	if (posFromV0>-1 && !usedTrack[posFromV0]) {
	307
	308	usedTrack[posFromV0] = kTRUE;
	309
	310	AliESDtrack *esdTrack = fESD->GetTrack(posFromV0);
	311	esdTrack->GetPxPyPz(p);
	312	esdTrack->GetXYZ(pos);
	313	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	314	esdTrack->GetESDpid(pid);
	315
	316	vV0FromCascade->AddDaughter(aodTrack =
	317	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	318	esdTrack->GetLabel(),
	319	p,
	320	kTRUE,
	321	pos,
	322	kFALSE,
	323	covTr,
	324	(Short_t)esdTrack->GetSign(),
	325	esdTrack->GetITSClusterMap(),
	326	pid,
	327	vV0FromCascade,
	328	kTRUE, // check if this is right
	329	kFALSE, // check if this is right
	330	AliAODTrack::kSecondary)
	331	);
	332	aodTrack->ConvertAliPIDtoAODPID();
	333	}
	334	else {
	335	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
	336	// << " track " << posFromV0 << " has already been used!" << endl;
	337	}
	338
	339	// Add the negative tracks from the V0
	340
	341	if (negFromV0>-1 && !usedTrack[negFromV0]) {
	342
	343	usedTrack[negFromV0] = kTRUE;
	344
	345	AliESDtrack *esdTrack = fESD->GetTrack(negFromV0);
	346	esdTrack->GetPxPyPz(p);
	347	esdTrack->GetXYZ(pos);
	348	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	349	esdTrack->GetESDpid(pid);
	350
	351	vV0FromCascade->AddDaughter(aodTrack =
	352	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	353	esdTrack->GetLabel(),
	354	p,
	355	kTRUE,
	356	pos,
	357	kFALSE,
	358	covTr,
	359	(Short_t)esdTrack->GetSign(),
	360	esdTrack->GetITSClusterMap(),
	361	pid,
	362	vV0FromCascade,
	363	kTRUE, // check if this is right
	364	kFALSE, // check if this is right
	365	AliAODTrack::kSecondary)
	366	);
	367	aodTrack->ConvertAliPIDtoAODPID();
	368	}
	369	else {
	370	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
	371	// << " track " << negFromV0 << " has already been used!" << endl;
	372	}
	373
	374	// Add the bachelor track from the cascade
	375
	376	Int_t bachelor = cascade->GetBindex();
	377
	378	if(bachelor>-1 && !usedTrack[bachelor]) {
	379
	380	usedTrack[bachelor] = kTRUE;
	381
	382	AliESDtrack *esdTrack = fESD->GetTrack(bachelor);
	383	esdTrack->GetPxPyPz(p);
	384	esdTrack->GetXYZ(pos);
	385	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	386	esdTrack->GetESDpid(pid);
	387
	388	vcascade->AddDaughter(aodTrack =
	389	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	390	esdTrack->GetLabel(),
	391	p,
	392	kTRUE,
	393	pos,
	394	kFALSE,
	395	covTr,
	396	(Short_t)esdTrack->GetSign(),
	397	esdTrack->GetITSClusterMap(),
	398	pid,
	399	vcascade,
	400	kTRUE, // check if this is right
	401	kFALSE, // check if this is right
	402	AliAODTrack::kSecondary)
	403	);
	404	aodTrack->ConvertAliPIDtoAODPID();
	405	}
	406	else {
	407	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " cascade " << nCascade
	408	// << " track " << bachelor << " has already been used!" << endl;
	409	}
	410
	411	// Add the primary track of the cascade (if any)
	412
	413	} // end of the loop on cascades
	414
	415	// V0s
	416
	417	for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
	418
	419	if (usedV0[nV0]) continue; // skip if aready added to the AOD
	420
	421	AliESDv0 *v0 = fESD->GetV0(nV0);
	422
	423	v0->GetXYZ(pos[0], pos[1], pos[2]);
	424
	425	if (!old) {
	426	chi2 = v0->GetChi2V0();
	427	v0->GetPosCov(covVtx);
	428	} else {
	429	chi2 = -999.;
	430	}
	431
	432
	433	AliAODVertex * vV0 =
	434	new(vertices[jVertices++]) AliAODVertex(pos,
	435	covVtx,
	436	chi2, // = chi2/NDF since NDF = 2*2-3
	437	primary,
	438	AliAODVertex::kV0);
	439	primary->AddDaughter(vV0);
	440
	441	Int_t posFromV0 = v0->GetPindex();
	442	Int_t negFromV0 = v0->GetNindex();
	443
	444	// Add the positive tracks from the V0
	445
	446	if (posFromV0>-1 && !usedTrack[posFromV0]) {
	447
	448	usedTrack[posFromV0] = kTRUE;
	449
	450	AliESDtrack *esdTrack = fESD->GetTrack(posFromV0);
	451	esdTrack->GetPxPyPz(p);
	452	esdTrack->GetXYZ(pos);
	453	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	454	esdTrack->GetESDpid(pid);
	455
	456	vV0->AddDaughter(aodTrack =
	457	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	458	esdTrack->GetLabel(),
	459	p,
	460	kTRUE,
	461	pos,
	462	kFALSE,
	463	covTr,
	464	(Short_t)esdTrack->GetSign(),
	465	esdTrack->GetITSClusterMap(),
	466	pid,
	467	vV0,
	468	kTRUE, // check if this is right
	469	kFALSE, // check if this is right
	470	AliAODTrack::kSecondary)
	471	);
	472	aodTrack->ConvertAliPIDtoAODPID();
	473	}
	474	else {
	475	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " V0 " << nV0
	476	// << " track " << posFromV0 << " has already been used!" << endl;
	477	}
	478
	479	// Add the negative tracks from the V0
	480
	481	if (negFromV0>-1 && !usedTrack[negFromV0]) {
	482
	483	usedTrack[negFromV0] = kTRUE;
	484
	485	AliESDtrack *esdTrack = fESD->GetTrack(negFromV0);
	486	esdTrack->GetPxPyPz(p);
	487	esdTrack->GetXYZ(pos);
	488	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	489	esdTrack->GetESDpid(pid);
	490
	491	vV0->AddDaughter(aodTrack =
	492	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	493	esdTrack->GetLabel(),
	494	p,
	495	kTRUE,
	496	pos,
	497	kFALSE,
	498	covTr,
	499	(Short_t)esdTrack->GetSign(),
	500	esdTrack->GetITSClusterMap(),
	501	pid,
	502	vV0,
	503	kTRUE, // check if this is right
	504	kFALSE, // check if this is right
	505	AliAODTrack::kSecondary)
	506	);
	507	aodTrack->ConvertAliPIDtoAODPID();
	508	}
	509	else {
	510	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " V0 " << nV0
	511	// << " track " << negFromV0 << " has already been used!" << endl;
	512	}
	513
	514	} // end of the loop on V0s
	515
	516	// Kinks: it is a big mess the access to the information in the kinks
	517	// The loop is on the tracks in order to find the mother and daugther of each kink
	518
	519
	520	for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
	521
	522
	523	AliESDtrack * esdTrack = fESD->GetTrack(iTrack);
	524
	525	Int_t ikink = esdTrack->GetKinkIndex(0);
	526
	527	if (ikink && nKinks) {
	528	// Negative kink index: mother, positive: daughter
	529
	530	// Search for the second track of the kink
	531
	532	for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
	533
	534	AliESDtrack * esdTrack1 = fESD->GetTrack(jTrack);
	535
	536	Int_t jkink = esdTrack1->GetKinkIndex(0);
	537
	538	if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
	539
	540	// The two tracks are from the same kink
	541
	542	if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
	543
	544	Int_t imother = -1;
	545	Int_t idaughter = -1;
	546
	547	if (ikink<0 && jkink>0) {
	548
	549	imother = iTrack;
	550	idaughter = jTrack;
	551	}
	552	else if (ikink>0 && jkink<0) {
	553
	554	imother = jTrack;
	555	idaughter = iTrack;
	556	}
	557	else {
	558	// cerr << "Error: Wrong combination of kink indexes: "
	559	// << ikink << " " << jkink << endl;
	560	continue;
	561	}
	562
	563	// Add the mother track
	564
	565	AliAODTrack * mother = NULL;
	566
	567	if (!usedTrack[imother]) {
	568
	569	usedTrack[imother] = kTRUE;
	570
	571	AliESDtrack *esdTrack = fESD->GetTrack(imother);
	572	esdTrack->GetPxPyPz(p);
	573	esdTrack->GetXYZ(pos);
	574	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	575	esdTrack->GetESDpid(pid);
	576
	577	mother =
	578	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	579	esdTrack->GetLabel(),
	580	p,
	581	kTRUE,
	582	pos,
	583	kFALSE,
	584	covTr,
	585	(Short_t)esdTrack->GetSign(),
	586	esdTrack->GetITSClusterMap(),
	587	pid,
	588	primary,
	589	kTRUE, // check if this is right
	590	kTRUE, // check if this is right
	591	AliAODTrack::kPrimary);
	592	primary->AddDaughter(mother);
	593	mother->ConvertAliPIDtoAODPID();
	594	}
	595	else {
	596	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
	597	// << " track " << imother << " has already been used!" << endl;
	598	}
	599
	600	// Add the kink vertex
	601	AliESDkink * kink = fESD->GetKink(TMath::Abs(ikink)-1);
	602
	603	AliAODVertex * vkink =
	604	new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
	605	NULL,
	606	0.,
	607	mother,
	608	AliAODVertex::kKink);
	609	// Add the daughter track
	610
	611	AliAODTrack * daughter = NULL;
	612
	613	if (!usedTrack[idaughter]) {
	614
	615	usedTrack[idaughter] = kTRUE;
	616
	617	AliESDtrack *esdTrack = fESD->GetTrack(idaughter);
	618	esdTrack->GetPxPyPz(p);
	619	esdTrack->GetXYZ(pos);
	620	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	621	esdTrack->GetESDpid(pid);
	622
	623	daughter =
	624	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	625	esdTrack->GetLabel(),
	626	p,
	627	kTRUE,
	628	pos,
	629	kFALSE,
	630	covTr,
	631	(Short_t)esdTrack->GetSign(),
	632	esdTrack->GetITSClusterMap(),
	633	pid,
	634	vkink,
	635	kTRUE, // check if this is right
	636	kTRUE, // check if this is right
	637	AliAODTrack::kPrimary);
	638	vkink->AddDaughter(daughter);
	639	daughter->ConvertAliPIDtoAODPID();
	640	}
	641	else {
	642	// cerr << "Error: event " << fESD->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
	643	// << " track " << idaughter << " has already been used!" << endl;
	644	}
	645	}
	646	}
	647	}
	648	}
	649
	650
	651	// Tracks (primary and orphan)
	652
	653	printf("NUMBER OF TRACKS %5d\n", nTracks);
	654
	655	for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
	656
	657
	658	if (usedTrack[nTrack]) continue;
	659	AliESDtrack *esdTrack = fESD->GetTrack(nTrack);
	660	UInt_t selectInfo = 0;
	661	//
	662	// Track selection
	663	if (fTrackFilter) {
	664	selectInfo = fTrackFilter->IsSelected(esdTrack);
	665	if (!selectInfo) continue;
	666	}
	667
	668	//
	669	esdTrack->GetPxPyPz(p);
	670	esdTrack->GetXYZ(pos);
	671	esdTrack->GetCovarianceXYZPxPyPz(covTr);
	672	esdTrack->GetESDpid(pid);
	673
	674	Float_t impactXY, impactZ;
	675
	676	esdTrack->GetImpactParameters(impactXY,impactZ);
	677
	678	if (impactXY<3) {
	679	// track inside the beam pipe
	680
	681	primary->AddDaughter(aodTrack =
	682	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	683	esdTrack->GetLabel(),
	684	p,
	685	kTRUE,
	686	pos,
	687	kFALSE,
	688	covTr,
	689	(Short_t)esdTrack->GetSign(),
	690	esdTrack->GetITSClusterMap(),
	691	pid,
	692	primary,
	693	kTRUE, // check if this is right
	694	kTRUE, // check if this is right
	695	AliAODTrack::kPrimary,
	696	selectInfo)
	697	);
	698	aodTrack->ConvertAliPIDtoAODPID();
	699	}
	700	else {
	701	// outside the beam pipe: orphan track
	702	aodTrack =
	703	new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
	704	esdTrack->GetLabel(),
	705	p,
	706	kTRUE,
	707	pos,
	708	kFALSE,
	709	covTr,
	710	(Short_t)esdTrack->GetSign(),
	711	esdTrack->GetITSClusterMap(),
	712	pid,
	713	NULL,
	714	kFALSE, // check if this is right
	715	kFALSE, // check if this is right
	716	AliAODTrack::kOrphan,
	717	selectInfo);
	718	aodTrack->ConvertAliPIDtoAODPID();
	719	}
	720	} // end of loop on tracks
	721
	722	// muon tracks
	723	Int_t nMuTracks = fESD->GetNumberOfMuonTracks();
	724	for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
	725
	726	AliESDMuonTrack *esdMuTrack = fESD->GetMuonTrack(nMuTrack);
	727	p[0] = esdMuTrack->Px();
	728	p[1] = esdMuTrack->Py();
	729	p[2] = esdMuTrack->Pz();
	730	pos[0] = primary->GetX();
	731	pos[1] = primary->GetY();
	732	pos[2] = primary->GetZ();
	733
	734	// has to be changed once the muon pid is provided by the ESD
	735	for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
	736
	737	primary->AddDaughter(
	738	new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
	739	0, // no label provided
	740	p,
	741	kTRUE,
	742	pos,
	743	kFALSE,
	744	NULL, // no covariance matrix provided
	745	(Short_t)-99, // no charge provided
	746	0, // no ITSClusterMap
	747	pid,
	748	primary,
	749	kTRUE, // check if this is right
	750	kTRUE, // not used for vertex fit
	751	AliAODTrack::kPrimary)
	752	);
	753	}
	754
	755	delete [] usedTrack;
	756	delete [] usedV0;
	757	delete [] usedKink;
	758
	759
	760	//
	761	PostData(0, fTreeA);
	762	return;
	763	}
	764
	765	void AliAnalysisTaskESDfilter::Terminate(Option_t /option*/)
	766	{
	767	// Terminate analysis
	768	//
	769	if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
	770	}
	771