[u/mrichter/AliRoot.git] / PWGCF / Correlations / Base / AliAnalysisTaskCFTree.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.                  *
 **************************************************************************/

// Analysis task to produce trees of lightweight events
// evgeny.kryshen@cern.ch

// aliroot
#include "AliAnalysisTaskCFTree.h"
#include "AliCFParticle.h"
#include "AliAnalysisManager.h"
#include "AliInputEventHandler.h"
#include "AliESDEvent.h"
#include "AliAODEvent.h"
#include "AliMCEvent.h"
#include "AliHeader.h"
#include "AliESDHeader.h"
#include "AliAODHeader.h"
#include "AliAODMCHeader.h"
#include "AliGenEventHeader.h"
#include "AliGenCocktailEventHeader.h"
#include "AliAODTrack.h"
#include "AliESDtrack.h"
#include "AliMCParticle.h"
#include "AliAODMCParticle.h"
#include "AliExternalTrackParam.h"
#include "AliCentrality.h"
#include "AliAnalysisFilter.h"
#include "AliPIDResponse.h"
#include "AliTPCPIDResponse.h"
// root
#include "TMath.h"
#include "TFile.h"
#include "TList.h"
#include "TH1I.h"
#include "TChain.h"
#include "TTree.h"
#include "TClonesArray.h"
ClassImp(AliAnalysisTaskCFTree)

//-----------------------------------------------------------------------------
AliAnalysisTaskCFTree::AliAnalysisTaskCFTree(const char* name) :
  AliAnalysisTask(name,""),
  fDebug(0),
  fMode(0),
  fInputHandler(0x0),
  fMcHandler(0x0),
  fTrackFilter(0x0),
  fHybridConstrainedMask(0),
  fTPConlyConstrainedMask(0),
  fPIDResponse(0x0),
  fListOfHistos(0x0),
  fEventStatistics(0x0),
  fTree(0x0),
  fParticles(0x0),
  fField(0),
  fCentrality(0),
  fZvtx(0),
  fRunNumber(0),
  fPeriod(0),
  fOrbit(),
  fBc(),
  fSelectMask(0),
  fSelectBit(AliVEvent::kMB),
  fZVertexCut(10.),
  fnTracksVertex(1),
  fCentralityMethod("V0M"),
  fTrackFilterBit(0xffffffff),
  fTrackEtaCut(1.0),
  fPtMin(0.15),
  fSharedClusterCut(0.4),
  fCrossedRowsCut(100),
  fFoundFractionCut(0.8)
{
  Info("AliAnalysisTaskCFTree","Calling Constructor");
  DefineInput(0,TChain::Class());
  DefineOutput(0,TList::Class());
  DefineOutput(1,TTree::Class());
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
void AliAnalysisTaskCFTree::ConnectInputData(Option_t* /*option*/){
  fInputHandler = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
  fMcHandler    = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
void AliAnalysisTaskCFTree::CreateOutputObjects(){
  fListOfHistos = new TList();
  fListOfHistos->SetOwner();
  fEventStatistics = new TH1I("fEventStatistics","",10,0,10);
  fEventStatistics->SetBit(TH1::kCanRebin);

  fListOfHistos->Add(fEventStatistics);

  fParticles = new TClonesArray("AliCFParticle",2000);
  // create file-resident tree
  TDirectory *owd = gDirectory;
  OpenFile(1);
  fTree = new TTree("events","events");
  owd->cd();
  fTree->Branch("cent",&fCentrality);
  fTree->Branch("zvtx",&fZvtx);
  fTree->Branch("field",&fField);
  fTree->Branch("run",&fRunNumber);
  fTree->Branch("period",&fPeriod);
  fTree->Branch("orbit",&fOrbit);
  fTree->Branch("bc",&fBc);
  fTree->Branch("mask",&fSelectMask);
  fTree->Branch("particles",&fParticles);
  PostData(0,fListOfHistos);
  PostData(1,fTree);
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
void AliAnalysisTaskCFTree::Exec(Option_t *){
  fParticles->Clear();
  
  fEventStatistics->Fill("before cuts",1);
  
  AliVEvent* event = fInputHandler->GetEvent();
  if (!event) return;
  fEventStatistics->Fill("after event check",1);
  
  fPIDResponse = fInputHandler->GetPIDResponse();
  
  fSelectMask = fInputHandler->IsEventSelected();
  if (!(fSelectMask & fSelectBit)) return;
  fEventStatistics->Fill("after trigger selection",1);
  fRunNumber  = event->GetRunNumber();
  fPeriod     = event->GetPeriodNumber();
  fOrbit      = event->GetOrbitNumber();
  fBc         = event->GetBunchCrossNumber();
  fField      = event->GetMagneticField();
  fCentrality = event->GetCentrality()->GetCentralityPercentile(fCentralityMethod);
  
  if (fMode==0){ // Data analysis
    const AliVVertex* vertex  = event->GetPrimaryVertex();
    if (!vertex) return;
    fZvtx  = vertex->GetZ();
    TString name(vertex->GetName());
    if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex")) return; // Reject TPC only vertex
    fEventStatistics->Fill("after rejection of TPC only vertex",1);
    if (TMath::Abs(fZvtx) >= fZVertexCut || vertex->GetNContributors()<fnTracksVertex)  return;
    fEventStatistics->Fill("after vertex selection",1);
    
    for (Int_t ipart=0;ipart<event->GetNumberOfTracks();ipart++){
      AliVTrack* track = (AliVTrack*) event->GetTrack(ipart);
      if (!track) continue;
      UInt_t mask = GetFilterMap(track);
      if (!(mask & fTrackFilterBit)) continue;

      if (track->InheritsFrom("AliAODTrack")) AddTrack(track,mask,0);
      else if (track->InheritsFrom("AliESDtrack")) {
        if (mask)                           AddTrack(track,mask,1);
        if (mask & fHybridConstrainedMask)  AddTrack(track,mask,2);
        if (mask & fTPConlyConstrainedMask) AddTrack(track,mask,3);
      }
    }
  }
  else { // MC analysis
    AliMCEvent* mcEvent = 0;
    TClonesArray* mcTracks = 0;
    Int_t nPrimGen = 0;
    Int_t nProduced = 0;

    mcEvent = fMcHandler ? fMcHandler->MCEvent() : 0;
    mcTracks = (TClonesArray*) event->FindListObject(AliAODMCParticle::StdBranchName());
    
    if (!mcEvent && !mcTracks) { printf("No mc object found\n"); return; }
    fEventStatistics->Fill("after mc objects check",1);

    if (mcEvent) {
      AliHeader* header = (AliHeader*) mcEvent->Header();
      AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
      AliGenEventHeader* mcHeader = dynamic_cast<AliGenEventHeader*> (cocktailHeader ? cocktailHeader->GetHeaders()->First() : header->GenEventHeader());
      nProduced = mcEvent->GetNumberOfTracks();
      nPrimGen = mcHeader->NProduced();
      fZvtx = mcEvent->GetPrimaryVertex()->GetZ();
    } else if (mcTracks) {
      AliAODMCHeader* mcHeader = (AliAODMCHeader*) event->FindListObject(AliAODMCHeader::StdBranchName());
      if (!mcHeader) { printf("AliAODMCHeader not found\n"); return; }
      if (mcHeader->GetCocktailHeaders()) {
        AliGenEventHeader* header0 =  mcHeader->GetCocktailHeader(0);
        if (!header0) { printf("first header expected but not found\n"); return; }
        nPrimGen = header0->NProduced();
      } else  nPrimGen = nProduced;
      fZvtx = mcHeader->GetVtxZ();
    } else return;
    fEventStatistics->Fill("after mc header check",1);

    if (TMath::Abs(fZvtx)>fZVertexCut) return;
    fEventStatistics->Fill("after MC vertex cut",1);

    const AliVVertex* vertex = event->GetPrimaryVertex();
    if (!vertex) return;
    fEventStatistics->Fill("after check for vertex object",1);
    TString name(vertex->GetName());
    if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex")) fZvtx=1000; // Reject TPC only vertex
    fEventStatistics->Fill("after rejection of TPC only vertex",1);
    if (vertex->GetNContributors()<fnTracksVertex)  fZvtx=-1000;
    fEventStatistics->Fill("after check for vertex contributors",1);
    
    UInt_t* masks = new UInt_t[nProduced];
    for (Int_t i=0;i<nProduced;i++) masks[i]=0;
    
    // Loop over reconstructed tracks to set masks
    for (Int_t ipart=0;ipart<event->GetNumberOfTracks();ipart++){
      AliVTrack* part = (AliVTrack*) event->GetTrack(ipart);
      if (!part) continue;
      Int_t label = TMath::Abs(part->GetLabel());
      if (label>=nProduced) continue;
      masks[label] |= GetFilterMap(part);
    }
    
    // Loop over mc tracks to be stored
    for (Int_t ipart=0;ipart<nProduced;ipart++){
      AliVParticle* part = 0;
      Bool_t isPrimary = 0;
      if (mcTracks) { // AOD analysis
        AliAODMCParticle* particle = (AliAODMCParticle*) mcTracks->At(ipart);
        if (!particle) continue;
        // skip injected signals
        AliAODMCParticle* mother = particle;
        while (!mother->IsPrimary()) mother = (AliAODMCParticle*) mcTracks->At(mother->GetMother());
        if (mother->GetLabel()>=nPrimGen) continue;
        part = particle;
        // check for primary
        isPrimary = particle->IsPhysicalPrimary();
      } else { // ESD analysis
        AliMCParticle* particle = (AliMCParticle*) mcEvent->GetTrack(ipart);
        if (!particle) continue;
        // skip injected signals
        AliMCParticle* mother = particle;
        while (mother->GetMother()>=0) mother = (AliMCParticle*) mcEvent->GetTrack(mother->GetMother());
        if (mother->GetLabel()>=nPrimGen) continue;
        part = particle;
        // check for primary
        isPrimary = mcEvent->IsPhysicalPrimary(ipart);
      }
      
      // store only primaries and all reconstructed (non-zero mask)
      Int_t mask = masks[ipart];
      if (isPrimary) mask |= (1 << 30);
      if (!mask) continue; 
      AddTrack(part,mask);
    }
    delete[] masks;
  }
  fTree->Fill();
  PostData(0,fListOfHistos);
  PostData(1,fTree);
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
UInt_t AliAnalysisTaskCFTree::GetFilterMap(AliVParticle* particle){
  UInt_t mask = 0;
  if (particle->InheritsFrom("AliAODTrack")) {
    AliAODTrack* part = (AliAODTrack*) particle;
    mask = part->GetFilterMap();
    Double_t nCrossedRaws      = part->GetTPCNCrossedRows();
    Double_t nFindableClusters = part->GetTPCNclsF();
    Double_t nSharedClusters   = part->GetTPCnclsS();
    Double_t nClusters         = part->GetTPCncls();
    Bool_t itsRefit            = part->GetStatus() & AliVTrack::kITSrefit;
    if (nCrossedRaws/nFindableClusters > fFoundFractionCut) mask |= (1 << 26);
    if (nCrossedRaws>fCrossedRowsCut)                       mask |= (1 << 27);
    if (itsRefit)                                           mask |= (1 << 28);
    if (nSharedClusters/nClusters<=fSharedClusterCut)       mask |= (1 << 29);
    if (part->GetLabel()<0)                                 mask |= (1 << 31);
  } else if (particle->InheritsFrom("AliESDtrack")){
    AliESDtrack* part = (AliESDtrack*) particle;
    if (!fTrackFilter) AliFatal("Track filter undefined");
    mask |= fTrackFilter->IsSelected(part);
  }
  
  return mask;
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
AliCFParticle* AliAnalysisTaskCFTree::AddTrack(AliVParticle* part, UInt_t mask, UInt_t flag){

  // skip neutral mc particles
  Char_t charge = part->Charge();
  if (charge==0) return NULL;
  
  // set pt,eta,phi
  Float_t pt=0,eta=0,phi=0;
  if (flag==0 || flag==1){ // AOD, MC or Global ESD tracks
    pt  = part->Pt();
    eta = part->Eta();
    phi = part->Phi();
    if  (flag==1) mask &= (~fHybridConstrainedMask) & (~fTPConlyConstrainedMask);
  } 
  else if (flag==2) { // Hybrid constrained tracks (ESD)
    AliESDtrack* track = (AliESDtrack*) part;
    const AliExternalTrackParam* param = track->GetConstrainedParam();
    pt  = param->Pt();
    eta = param->Eta();
    phi = param->Phi();
    mask &= fHybridConstrainedMask;
  } 
  else if (flag==3) { // TPC only constrained tracks (ESD)
    AliESDtrack* track = (AliESDtrack*) part;
    AliESDtrack tpcTrack;
    if (!track->FillTPCOnlyTrack(tpcTrack)) return NULL;
    AliExternalTrackParam param;
    const AliESDVertex* vtxSPD = ((AliESDEvent*) fInputHandler->GetEvent())->GetPrimaryVertexSPD();
    if (!tpcTrack.RelateToVertexTPC(vtxSPD,fField,1.e30,&param)) return NULL;
    pt  = param.Pt();
    eta = param.Eta();
    phi = param.Phi();
    mask &= fTPConlyConstrainedMask;
  }

  // kinematic cuts
  if (pt < fPtMin || TMath::Abs(eta) > fTrackEtaCut) return NULL;

  return new ((*fParticles)[fParticles->GetEntriesFast()]) AliCFParticle(pt,eta,phi,charge,mask);
}
Commit	Line	Data
5fe7e785	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	// Analysis task to produce trees of lightweight events
	17	// evgeny.kryshen@cern.ch
	18
	19	// aliroot
	20	#include "AliAnalysisTaskCFTree.h"
	21	#include "AliCFParticle.h"
	22	#include "AliAnalysisManager.h"
	23	#include "AliInputEventHandler.h"
	24	#include "AliESDEvent.h"
	25	#include "AliAODEvent.h"
	26	#include "AliMCEvent.h"
	27	#include "AliHeader.h"
	28	#include "AliESDHeader.h"
	29	#include "AliAODHeader.h"
	30	#include "AliAODMCHeader.h"
	31	#include "AliGenEventHeader.h"
	32	#include "AliGenCocktailEventHeader.h"
	33	#include "AliAODTrack.h"
	34	#include "AliESDtrack.h"
	35	#include "AliMCParticle.h"
	36	#include "AliAODMCParticle.h"
	37	#include "AliExternalTrackParam.h"
	38	#include "AliCentrality.h"
	39	#include "AliAnalysisFilter.h"
	40	#include "AliPIDResponse.h"
	41	#include "AliTPCPIDResponse.h"
	42	// root
	43	#include "TMath.h"
	44	#include "TFile.h"
	45	#include "TList.h"
	46	#include "TH1I.h"
	47	#include "TChain.h"
	48	#include "TTree.h"
	49	#include "TClonesArray.h"
	50	ClassImp(AliAnalysisTaskCFTree)
	51
	52	//-----------------------------------------------------------------------------
	53	AliAnalysisTaskCFTree::AliAnalysisTaskCFTree(const char* name) :
	54	AliAnalysisTask(name,""),
	55	fDebug(0),
	56	fMode(0),
	57	fInputHandler(0x0),
	58	fMcHandler(0x0),
	59	fTrackFilter(0x0),
	60	fHybridConstrainedMask(0),
	61	fTPConlyConstrainedMask(0),
	62	fPIDResponse(0x0),
	63	fListOfHistos(0x0),
	64	fEventStatistics(0x0),
65	fTree(0x0),
66	fParticles(0x0),
67	fField(0),
68	fCentrality(0),
69	fZvtx(0),
70	fRunNumber(0),
71	fPeriod(0),
72	fOrbit(),
73	fBc(),
74	fSelectMask(0),
75	fSelectBit(AliVEvent::kMB),
76	fZVertexCut(10.),
77	fnTracksVertex(1),
78	fCentralityMethod("V0M"),
79	fTrackFilterBit(0xffffffff),
80	fTrackEtaCut(1.0),
81	fPtMin(0.15),
82	fSharedClusterCut(0.4),
83	fCrossedRowsCut(100),
84	fFoundFractionCut(0.8)
85	{
86	Info("AliAnalysisTaskCFTree","Calling Constructor");
87	DefineInput(0,TChain::Class());
88	DefineOutput(0,TList::Class());
89	DefineOutput(1,TTree::Class());
90	}
91	//-----------------------------------------------------------------------------
92
93
94	//-----------------------------------------------------------------------------
95	void AliAnalysisTaskCFTree::ConnectInputData(Option_t* /option/){
96	fInputHandler = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
97	fMcHandler = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
98	}
99	//-----------------------------------------------------------------------------
100
101
102	//-----------------------------------------------------------------------------
103	void AliAnalysisTaskCFTree::CreateOutputObjects(){
104	fListOfHistos = new TList();
105	fListOfHistos->SetOwner();
106	fEventStatistics = new TH1I("fEventStatistics","",10,0,10);
107	fEventStatistics->SetBit(TH1::kCanRebin);
108
109	fListOfHistos->Add(fEventStatistics);
110
111	fParticles = new TClonesArray("AliCFParticle",2000);
112	// create file-resident tree
113	TDirectory *owd = gDirectory;
114	OpenFile(1);
115	fTree = new TTree("events","events");
116	owd->cd();
117	fTree->Branch("cent",&fCentrality);
118	fTree->Branch("zvtx",&fZvtx);
119	fTree->Branch("field",&fField);
120	fTree->Branch("run",&fRunNumber);
121	fTree->Branch("period",&fPeriod);
122	fTree->Branch("orbit",&fOrbit);
123	fTree->Branch("bc",&fBc);
124	fTree->Branch("mask",&fSelectMask);
125	fTree->Branch("particles",&fParticles);
126	PostData(0,fListOfHistos);
127	PostData(1,fTree);
128	}
129	//-----------------------------------------------------------------------------
130
131
132	//-----------------------------------------------------------------------------
133	void AliAnalysisTaskCFTree::Exec(Option_t *){
134	fParticles->Clear();
135
136	fEventStatistics->Fill("before cuts",1);
137
138	AliVEvent* event = fInputHandler->GetEvent();
139	if (!event) return;
140	fEventStatistics->Fill("after event check",1);
141
142	fPIDResponse = fInputHandler->GetPIDResponse();
143
144	fSelectMask = fInputHandler->IsEventSelected();
145	if (!(fSelectMask & fSelectBit)) return;
146	fEventStatistics->Fill("after trigger selection",1);
147	fRunNumber = event->GetRunNumber();
148	fPeriod = event->GetPeriodNumber();
149	fOrbit = event->GetOrbitNumber();
150	fBc = event->GetBunchCrossNumber();
151	fField = event->GetMagneticField();
152	fCentrality = event->GetCentrality()->GetCentralityPercentile(fCentralityMethod);
153
154	if (fMode==0){ // Data analysis
155	const AliVVertex* vertex = event->GetPrimaryVertex();
156	if (!vertex) return;
157	fZvtx = vertex->GetZ();
158	TString name(vertex->GetName());
159	if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex")) return; // Reject TPC only vertex
160	fEventStatistics->Fill("after rejection of TPC only vertex",1);
161	if (TMath::Abs(fZvtx) >= fZVertexCut \|\| vertex->GetNContributors()<fnTracksVertex) return;
162	fEventStatistics->Fill("after vertex selection",1);
163
164	for (Int_t ipart=0;ipart<event->GetNumberOfTracks();ipart++){
165	AliVTrack* track = (AliVTrack*) event->GetTrack(ipart);
166	if (!track) continue;
167	UInt_t mask = GetFilterMap(track);
168	if (!(mask & fTrackFilterBit)) continue;
169
170	if (track->InheritsFrom("AliAODTrack")) AddTrack(track,mask,0);
171	else if (track->InheritsFrom("AliESDtrack")) {
172	if (mask) AddTrack(track,mask,1);
173	if (mask & fHybridConstrainedMask) AddTrack(track,mask,2);
174	if (mask & fTPConlyConstrainedMask) AddTrack(track,mask,3);
175	}
176	}
177	}
178	else { // MC analysis
179	AliMCEvent* mcEvent = 0;
180	TClonesArray* mcTracks = 0;
181	Int_t nPrimGen = 0;
182	Int_t nProduced = 0;
183
184	mcEvent = fMcHandler ? fMcHandler->MCEvent() : 0;
185	mcTracks = (TClonesArray*) event->FindListObject(AliAODMCParticle::StdBranchName());
186
187	if (!mcEvent && !mcTracks) { printf("No mc object found\n"); return; }
188	fEventStatistics->Fill("after mc objects check",1);
189
190	if (mcEvent) {
191	AliHeader* header = (AliHeader*) mcEvent->Header();
192	AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
193	AliGenEventHeader* mcHeader = dynamic_cast<AliGenEventHeader*> (cocktailHeader ? cocktailHeader->GetHeaders()->First() : header->GenEventHeader());
194	nProduced = mcEvent->GetNumberOfTracks();
195	nPrimGen = mcHeader->NProduced();
196	fZvtx = mcEvent->GetPrimaryVertex()->GetZ();
197	} else if (mcTracks) {
198	AliAODMCHeader* mcHeader = (AliAODMCHeader*) event->FindListObject(AliAODMCHeader::StdBranchName());
5a502471	199	if (!mcHeader) { printf("AliAODMCHeader not found\n"); return; }
	200	if (mcHeader->GetCocktailHeaders()) {
	201	AliGenEventHeader* header0 = mcHeader->GetCocktailHeader(0);
	202	if (!header0) { printf("first header expected but not found\n"); return; }
	203	nPrimGen = header0->NProduced();
	204	} else nPrimGen = nProduced;
5fe7e785	205	fZvtx = mcHeader->GetVtxZ();
	206	} else return;
	207	fEventStatistics->Fill("after mc header check",1);
	208
	209	if (TMath::Abs(fZvtx)>fZVertexCut) return;
	210	fEventStatistics->Fill("after MC vertex cut",1);
	211
	212	const AliVVertex* vertex = event->GetPrimaryVertex();
	213	if (!vertex) return;
	214	fEventStatistics->Fill("after check for vertex object",1);
	215	TString name(vertex->GetName());
	216	if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex")) fZvtx=1000; // Reject TPC only vertex
	217	fEventStatistics->Fill("after rejection of TPC only vertex",1);
	218	if (vertex->GetNContributors()<fnTracksVertex) fZvtx=-1000;
	219	fEventStatistics->Fill("after check for vertex contributors",1);
	220
	221	UInt_t* masks = new UInt_t[nProduced];
	222	for (Int_t i=0;i<nProduced;i++) masks[i]=0;
	223
	224	// Loop over reconstructed tracks to set masks
	225	for (Int_t ipart=0;ipart<event->GetNumberOfTracks();ipart++){
	226	AliVTrack* part = (AliVTrack*) event->GetTrack(ipart);
	227	if (!part) continue;
	228	Int_t label = TMath::Abs(part->GetLabel());
	229	if (label>=nProduced) continue;
	230	masks[label] \|= GetFilterMap(part);
	231	}
	232
	233	// Loop over mc tracks to be stored
	234	for (Int_t ipart=0;ipart<nProduced;ipart++){
	235	AliVParticle* part = 0;
	236	Bool_t isPrimary = 0;
	237	if (mcTracks) { // AOD analysis
	238	AliAODMCParticle* particle = (AliAODMCParticle*) mcTracks->At(ipart);
	239	if (!particle) continue;
	240	// skip injected signals
	241	AliAODMCParticle* mother = particle;
	242	while (!mother->IsPrimary()) mother = (AliAODMCParticle*) mcTracks->At(mother->GetMother());
	243	if (mother->GetLabel()>=nPrimGen) continue;
	244	part = particle;
	245	// check for primary
	246	isPrimary = particle->IsPhysicalPrimary();
	247	} else { // ESD analysis
	248	AliMCParticle* particle = (AliMCParticle*) mcEvent->GetTrack(ipart);
	249	if (!particle) continue;
	250	// skip injected signals
	251	AliMCParticle* mother = particle;
	252	while (mother->GetMother()>=0) mother = (AliMCParticle*) mcEvent->GetTrack(mother->GetMother());
	253	if (mother->GetLabel()>=nPrimGen) continue;
	254	part = particle;
	255	// check for primary
	256	isPrimary = mcEvent->IsPhysicalPrimary(ipart);
	257	}
	258
	259	// store only primaries and all reconstructed (non-zero mask)
	260	Int_t mask = masks[ipart];
	261	if (isPrimary) mask \|= (1 << 30);
	262	if (!mask) continue;
	263	AddTrack(part,mask);
	264	}
	265	delete[] masks;
	266	}
	267	fTree->Fill();
	268	PostData(0,fListOfHistos);
269	PostData(1,fTree);
270	}
271	//-----------------------------------------------------------------------------
272
273
274	//-----------------------------------------------------------------------------
275	UInt_t AliAnalysisTaskCFTree::GetFilterMap(AliVParticle* particle){
276	UInt_t mask = 0;
277	if (particle->InheritsFrom("AliAODTrack")) {
278	AliAODTrack* part = (AliAODTrack*) particle;
279	mask = part->GetFilterMap();
280	Double_t nCrossedRaws = part->GetTPCNCrossedRows();
281	Double_t nFindableClusters = part->GetTPCNclsF();
282	Double_t nSharedClusters = part->GetTPCnclsS();
283	Double_t nClusters = part->GetTPCncls();
284	Bool_t itsRefit = part->GetStatus() & AliVTrack::kITSrefit;
285	if (nCrossedRaws/nFindableClusters > fFoundFractionCut) mask \|= (1 << 26);
286	if (nCrossedRaws>fCrossedRowsCut) mask \|= (1 << 27);
287	if (itsRefit) mask \|= (1 << 28);
288	if (nSharedClusters/nClusters<=fSharedClusterCut) mask \|= (1 << 29);
289	if (part->GetLabel()<0) mask \|= (1 << 31);
290	} else if (particle->InheritsFrom("AliESDtrack")){
291	AliESDtrack* part = (AliESDtrack*) particle;
292	if (!fTrackFilter) AliFatal("Track filter undefined");
293	mask \|= fTrackFilter->IsSelected(part);
294	}
295
296	return mask;
297	}
298	//-----------------------------------------------------------------------------
299
300
301	//-----------------------------------------------------------------------------
302	AliCFParticle* AliAnalysisTaskCFTree::AddTrack(AliVParticle* part, UInt_t mask, UInt_t flag){
303
304	// skip neutral mc particles
305	Char_t charge = part->Charge();
306	if (charge==0) return NULL;
307
308	// set pt,eta,phi
309	Float_t pt=0,eta=0,phi=0;
310	if (flag==0 \|\| flag==1){ // AOD, MC or Global ESD tracks
311	pt = part->Pt();
312	eta = part->Eta();
313	phi = part->Phi();
314	if (flag==1) mask &= (~fHybridConstrainedMask) & (~fTPConlyConstrainedMask);
315	}
316	else if (flag==2) { // Hybrid constrained tracks (ESD)
317	AliESDtrack* track = (AliESDtrack*) part;
318	const AliExternalTrackParam* param = track->GetConstrainedParam();
319	pt = param->Pt();
320	eta = param->Eta();
321	phi = param->Phi();
322	mask &= fHybridConstrainedMask;
323	}
324	else if (flag==3) { // TPC only constrained tracks (ESD)
325	AliESDtrack* track = (AliESDtrack*) part;
326	AliESDtrack tpcTrack;
327	if (!track->FillTPCOnlyTrack(tpcTrack)) return NULL;
328	AliExternalTrackParam param;
329	const AliESDVertex* vtxSPD = ((AliESDEvent*) fInputHandler->GetEvent())->GetPrimaryVertexSPD();
330	if (!tpcTrack.RelateToVertexTPC(vtxSPD,fField,1.e30,&param)) return NULL;
331	pt = param.Pt();
332	eta = param.Eta();
333	phi = param.Phi();
334	mask &= fTPConlyConstrainedMask;
335	}
336
337	// kinematic cuts
338	if (pt < fPtMin \|\| TMath::Abs(eta) > fTrackEtaCut) return NULL;
339
340	return new ((*fParticles)[fParticles->GetEntriesFast()]) AliCFParticle(pt,eta,phi,charge,mask);
341	}