[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEtrackFilter.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$ */

//
// Track filter class 
// Apply cut steps to all tracks in one event and returns a list of
// filtered tracks
//
// Author:
//   Markus Fasel <M.Fasel@gsi.de>
//
#include <TMath.h>
#include <TString.h>

#include "AliAnalysisCuts.h"
#include "AliCFContainer.h"
#include "AliCFParticleGenCuts.h"
#include "AliCFTrackKineCuts.h"
#include "AliCFTrackIsPrimaryCuts.h"
#include "AliCFTrackQualityCuts.h"
#include "AliESDEvent.h"
#include "AliESDtrack.h"
#include "AliMCEvent.h"
#include "AliMCParticle.h"

#include "AliHFEcontainer.h"
#include "AliHFEcutStep.h"
#include "AliHFEextraCuts.h"
#include "AliHFEtrackFilter.h"
#include "AliHFEtools.h"

ClassImp(AliHFEtrackFilter)

//__________________________________________________________________
AliHFEtrackFilter::AliHFEtrackFilter(const Char_t *name) :
  TNamed(name, ""),
  fFilteredTracks(NULL),
  fCutSteps(NULL),
  fEfficiencyContainers(NULL),
  fMC(NULL),
  fMCsignal(NULL),
  fPtBins(0),
  fEtaBins(0),
  fPhiBins(0),
  fPtBinning(NULL),
  fEtaBinning(NULL),
  fPhiBinning(NULL)
{
  //
  // Constructor
  //
  fCutSteps = new TObjArray;
  fFilteredTracks = new TObjArray;
  fEfficiencyContainers = new TObjArray(4);
}

//__________________________________________________________________
AliHFEtrackFilter::AliHFEtrackFilter(const AliHFEtrackFilter &ref):
  TNamed(ref),
  fFilteredTracks(NULL),
  fCutSteps(NULL),
  fEfficiencyContainers(NULL),
  fMC(NULL),
  fMCsignal(NULL),
  fPtBins(0),
  fEtaBins(0),
  fPhiBins(0),
  fPtBinning(NULL),
  fEtaBinning(NULL),
  fPhiBinning(NULL)
{
  //
  // Copy constructor
  //
  ref.Copy(*this);
}

//__________________________________________________________________
AliHFEtrackFilter &AliHFEtrackFilter::operator=(const AliHFEtrackFilter &o){
  //
  // Assignment operator
  //
  if(this != &o){
    o.Copy(*this);
  }
  return *this;
}

//__________________________________________________________________
AliHFEtrackFilter::~AliHFEtrackFilter(){
  //
  // Destructor
  //

  // remove correction framework content
  if(fPtBinning) delete fPtBinning;
  if(fEtaBinning) delete fEtaBinning;
  if(fPhiBinning) delete fPhiBinning;
  if(TestBit(kOwnCFContainers)) 
    fEfficiencyContainers->Delete();
  delete fEfficiencyContainers;
  // remove default content
  delete fFilteredTracks;
  delete fCutSteps;
}

//__________________________________________________________________
void AliHFEtrackFilter::Copy(TObject &o) const{
  //
  // Copy content into the object o
  //
  TNamed::Copy(o);
  AliHFEtrackFilter &target = dynamic_cast<AliHFEtrackFilter &>(o);

  // Make copy
  target.fFilteredTracks = dynamic_cast<TObjArray *>(fFilteredTracks->Clone());
  target.fCutSteps = dynamic_cast<TObjArray *>(fCutSteps->Clone());
  target.fEfficiencyContainers = dynamic_cast<TObjArray *>(fEfficiencyContainers->Clone());
  target.fMC = fMC;
  target.fMCsignal = dynamic_cast<AliHFEcutStep *>(fMCsignal->Clone());
  
  target.fPtBins = fPtBins;
  target.fEtaBins = fEtaBins;
  target.fPhiBins = fPhiBins;
  target.fPtBinning = new Double_t[fPtBins];
  memcpy(target.fPtBinning, fPtBinning, sizeof(Double_t) * fPtBins);
  target.fEtaBinning = new Double_t[fEtaBins];
  memcpy(target.fEtaBinning, fEtaBinning, sizeof(Double_t) * fEtaBins);
  target.fPhiBinning = new Double_t[fPhiBins];
  memcpy(target.fPhiBinning, fPhiBinning, sizeof(Double_t) * fPhiBins);
}

//__________________________________________________________________
void AliHFEtrackFilter::InitCF(){
  //
  // Initialize correction framework container
  // No HFE container given
  // Only for testing purpose
  //
   const Char_t *cnames[4] = {"_container", "_container_signal", "_container_MC", "_container_signalMC"};
  const Char_t *ctitlesAppends[4] = {"(All Rec)", "(All Signals)", "(All Rec (MC))", "(All Signals(MC))"};

  // Create the binning if not done from outside
  if(!fPtBins){
    fPtBins = 40;
    fPtBinning = AliHFEtools::MakeLogarithmicBinning(fPtBins, 0.1, 10); 
  }
  if(!fEtaBins){
    fEtaBins = 8;
    fEtaBinning = AliHFEtools::MakeLinearBinning(fEtaBins, -0.9, 0.9); 
  }
  if(!fPhiBins){
    fPhiBins = 18;
    fPhiBinning = AliHFEtools::MakeLinearBinning(fPhiBins, 0, 2*TMath::Pi());
  }
  Double_t chargeBins[3] = {-1.1, 0., 1.1};
  Int_t nStep = fCutSteps->GetEntriesFast()+1;
  Int_t nBins[4] = {fPtBins, fEtaBins, fPhiBins + 2};
  AliCFContainer *ctmp = NULL;
  for(Int_t icont = 0; icont < 4; icont++){
    TString containername =  GetName() + TString(cnames[icont]);
    TString containertitle = TString("Container for filter ") + GetName() + TString(ctitlesAppends[icont]);
    ctmp = new AliCFContainer(containername.Data(), containertitle.Data(), nStep, 4, nBins);
    SetBit(kOwnCFContainers, kTRUE);

    // Set the binning
    ctmp->SetBinLimits(0, fPtBinning);
    ctmp->SetBinLimits(1, fEtaBinning);
    ctmp->SetBinLimits(2, fPhiBinning);
    ctmp->SetBinLimits(3, chargeBins);

    // Label variable names
    ctmp->SetVarTitle(0, "pt");
    ctmp->SetVarTitle(1, "eta");
    ctmp->SetVarTitle(2, "phi");
    ctmp->SetVarTitle(3, "charge");

    // Label step name
    ctmp->SetStepTitle(0, "No Cuts");
    AliHFEcutStep *cutStep = NULL;
    for(Int_t istep = 0; istep < fCutSteps->GetEntriesFast(); istep++){
      cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(istep)); 
      if(cutStep) ctmp->SetStepTitle(istep + 1, cutStep->GetName());
    }
    fEfficiencyContainers->AddAt(ctmp, 0);
  }
  OwnContainers();
}

//__________________________________________________________________
void AliHFEtrackFilter::InitCF(AliHFEcontainer *cont){
  //
  // Initialize Correction Framework container 
  // Use the given HFE container
  // HFE standard way to use the correction Framework
  //
  
  const Char_t *cnames[4] = {"_container", "_container_signal", "_container_MC", "_container_signalMC"};
  const Char_t *ctitlesAppends[4] = {"(All Rec)", "(All Signals)", "(All Rec (MC))", "(All Signals(MC))"};
  Int_t nStep = fCutSteps->GetEntriesFast()+1;
  AliCFContainer *ctmp = NULL;
  for(Int_t icont = 0; icont < 4; icont++){
    TString contname = GetName() + TString(cnames[icont]); 
    TString contTitle = TString("Container for filter ") + GetName() + TString(ctitlesAppends[icont]);
    //printf("Adding container %s: %s\n", contname.Data(), contTitle.Data());
    cont->CreateContainer(contname.Data(), contTitle.Data(), nStep);
    fEfficiencyContainers->AddAt((ctmp = cont->GetCFContainer(contname.Data())), icont);
    
    // Label step name
    ctmp->SetStepTitle(0, "No Cuts");
    AliHFEcutStep *cutStep = NULL;
    for(Int_t istep = 0; istep < fCutSteps->GetEntriesFast(); istep++){
      cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(istep)); 
      if(cutStep) ctmp->SetStepTitle(istep + 1, cutStep->GetName());
    }
  }
  ReleaseContainers();
}

//__________________________________________________________________
void AliHFEtrackFilter::FilterTracks(AliESDEvent * const event){
  //
  // Perform track filtering
  // Check each cut step one by one and select tracks which pass
  // all cuts. If the correction framework is initialized, a correction#
  // framework container will be filled for each cut step.
  //
  AliHFEcutStep *cutStep = NULL;
  AliESDtrack *track = NULL;
  AliMCParticle *mctrack = NULL;
  Double_t cont[4] = {0., 0., 0., 0.}, contMC[4] = {0., 0., 0., 0.};
  Bool_t goodTrack = kTRUE, signal = kFALSE;
  Int_t nStep = fCutSteps->GetEntriesFast();
  AliCFContainer *call = dynamic_cast<AliCFContainer *>(fEfficiencyContainers->At(0));
  AliCFContainer *cSignal = dynamic_cast<AliCFContainer *>(fEfficiencyContainers->At(1));
  AliCFContainer *callMC = dynamic_cast<AliCFContainer *>(fEfficiencyContainers->At(2));
  AliCFContainer *cSignalMC = dynamic_cast<AliCFContainer *>(fEfficiencyContainers->At(3));
  for(Int_t itrack = 0; itrack < event->GetNumberOfTracks(); itrack++){
    signal = kFALSE;
    track = event->GetTrack(itrack);
    // check Monte-Carlo Information if available 
    if(fMC){
      mctrack = dynamic_cast<AliMCParticle *>(fMC->GetTrack(TMath::Abs(track->GetLabel())));
      if(mctrack){
        //AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMC->GetTrack(mctrack->Particle()->GetFirstMother()));
        //AliInfo(Form("Label %d, Mother %d", track->GetLabel(), mother->Particle()->GetPdgCode()));
        contMC[0] = mctrack->Pt();
        contMC[1] = mctrack->Eta();
        contMC[2] = mctrack->Phi();
        contMC[3] = mctrack->Charge();
        if(fMCsignal->IsSelected(mctrack)) signal = kTRUE;
        //if(TMath::Abs(mother->Particle()->GetPdgCode()) != 443) signal = kFALSE;
        //AliInfo(Form("Signal? %s", signal ? "Yes": "No"));
      }
    }

    // Fill Array without cut
    cont[0] = track->Pt();
    cont[1] = track->Eta();
    cont[2] = track->Phi();
    cont[3] = track->Charge();
    if(call) call->Fill(cont, 0);
    if(callMC) callMC->Fill(contMC, 0);
    if(signal){
      if(cSignal) cSignal->Fill(cont, 0);
      if(cSignalMC) cSignalMC->Fill(contMC, 0);
    }
    // cut the track
    goodTrack = kTRUE;
    for(Int_t icut = 0; icut < nStep; icut++){
      cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut));
      if(cutStep && (!cutStep->IsSelected(track))){
        // track cut away
        goodTrack = kFALSE;
        break;
      }

      // Track survived cut step: Fill container
      if(call) call->Fill(cont, icut + 1);
      if(callMC) callMC->Fill(contMC, icut + 1);
      if(signal){
        if(cSignal) cSignal->Fill(cont, icut + 1);
        if(cSignalMC) cSignalMC->Fill(contMC, icut + 1);
      }
    }

    // Append track to the list of filtered tracks
    if(goodTrack) fFilteredTracks->Add(track);
  }
}

//__________________________________________________________________
void AliHFEtrackFilter::GenerateCutSteps(){
  //
  // Make the default cuts
  //
  MakeCutStepRecKineITSTPC();
  MakeCutStepPrimary();
  MakeCutStepHFEITS();
  MakeCutStepHFETRD();

  MakeMCSignalCuts();
}

//__________________________________________________________________
void AliHFEtrackFilter::AddCutStep(AliHFEcutStep *step){
  //
  // Add new cut step to the filter
  //
  if(!fCutSteps->FindObject(step->GetName())) fCutSteps->Add(step);
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::GetCutStep(Int_t istep){
  //
  // Getter for single cut step
  //
  return dynamic_cast<AliHFEcutStep *>(fCutSteps->At(istep));
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::GetCutStep(const Char_t *name){
  //
  // Getter for single cut step (by step name)
  //
  return dynamic_cast<AliHFEcutStep *>(fCutSteps->FindObject(name));
}

//__________________________________________________________________
void AliHFEtrackFilter::Flush(){
  //
  // Empty track container
  //
  fFilteredTracks->Clear();
}

//__________________________________________________________________
void AliHFEtrackFilter::SetPtBins(Int_t nBins, Double_t *binning){
  //
  // User defined pt binning
  //
  fPtBins = nBins;
  fPtBinning = new Double_t[fPtBins + 1];
  memcpy(fPtBinning, binning, sizeof(Double_t) * nBins);
}

//__________________________________________________________________
void AliHFEtrackFilter::SetEtaBins(Int_t nBins, Double_t *binning){
  //
  // User defined eta binning
  //
  fEtaBins = nBins;
  fEtaBinning = new Double_t[fEtaBins + 1];
  memcpy(fEtaBinning, binning, sizeof(Double_t) * nBins);
}

//__________________________________________________________________
void AliHFEtrackFilter::SetPhiBins(Int_t nBins, Double_t *binning){
  //
  // User defined phi binning
  //
  fPhiBins = nBins;
  fPhiBinning = new Double_t[fPhiBins + 1];
  memcpy(fPhiBinning, binning, sizeof(Double_t) * nBins);
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::MakeCutStepRecKineITSTPC(){
  //
  // Make the cut step for Rec Kine
  // Cut step is already included in the filter
  //
  AliHFEcutStep *fCutStep = new AliHFEcutStep("RecKineITSTPC");

  AliCFTrackQualityCuts *trackQuality = new AliCFTrackQualityCuts((Char_t *)"QualityRec", (Char_t *)"REC Track Quality Cuts");
  trackQuality->SetMinNClusterTPC(80);
  trackQuality->SetMaxChi2PerClusterTPC(3.5);
  trackQuality->SetStatus(AliESDtrack::kTPCrefit | AliESDtrack::kITSrefit);
  trackQuality->SetMaxCovDiagonalElements(2., 2., 0.5, 0.5, 2); 
  fCutStep->AddCut(trackQuality);

  AliHFEextraCuts *hfecuts = new AliHFEextraCuts("HFETPC","Extra cuts from the HFE group");
  hfecuts->SetClusterRatioTPC(0.6, AliHFEextraCuts::kFoundOverCR);
  fCutStep->AddCut(hfecuts);
  
  AliCFTrackKineCuts *kineCuts = new AliCFTrackKineCuts((Char_t *)"RecKine", (Char_t *)"REC Kine Cuts");
  kineCuts->SetPtRange(0.1, 20);
  kineCuts->SetEtaRange(-0.8, 0.8);
  fCutStep->AddCut(kineCuts); 

  AddCutStep(fCutStep);
  return fCutStep;
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::MakeCutStepPrimary(){
  // 
  // Make cut on primaries
  // Cut step is already included in the filter
  //
  AliHFEcutStep *fCutStep = new AliHFEcutStep("Primary");

  AliCFTrackIsPrimaryCuts *primaryCut = new AliCFTrackIsPrimaryCuts((Char_t *)"PrimaryCuts", (Char_t *)"REC Primary Cuts");
  primaryCut->SetMaxDCAToVertexXY(0.5);
  primaryCut->SetMaxDCAToVertexZ(5.);
  primaryCut->SetAcceptKinkDaughters(kFALSE);
  fCutStep->AddCut(primaryCut);

  AddCutStep(fCutStep);
  return fCutStep;
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::MakeCutStepHFEITS(){
  //
  // Add special ITS cuts
  // Cut step is already included in the filter
  //
  AliHFEcutStep *fCutStep = new AliHFEcutStep("HFEITS");

  AliHFEextraCuts *hfecuts = new AliHFEextraCuts((Char_t *)"HFEPixelsCuts",(Char_t *)"Extra cuts from the HFE group");
  hfecuts->SetRequireITSpixel(AliHFEextraCuts::kFirst);
  //hfecuts->SetCheckITSstatus(kTRUE);
  fCutStep->AddCut(hfecuts);

  AddCutStep(fCutStep);
  return fCutStep;
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::MakeCutStepHFETRD(){
  //
  // Add special TRD cut
  // Cut step is already included in the filter
  //
  AliHFEcutStep *fCutStep = new AliHFEcutStep("HFETRD");

  AliHFEextraCuts *hfecuts = new AliHFEextraCuts("HFETRDCuts","Extra cuts from the HFE group");
  hfecuts->SetMinTrackletsTRD(0);
  fCutStep->AddCut(hfecuts);

  AddCutStep(fCutStep);
  return fCutStep;
}

//__________________________________________________________________
AliHFEcutStep *AliHFEtrackFilter::MakeMCSignalCuts(){
  //
  // Define MC Signal
  // Cut step is already included in the filter
  //
  fMCsignal = new AliHFEcutStep("MCSignal");
  AliCFParticleGenCuts *genCuts = new AliCFParticleGenCuts((Char_t *)"fCutsGenMC", (Char_t *)"Particle Generation Cuts");
  genCuts->SetRequireIsCharged();
  genCuts->SetRequireIsPrimary();
  genCuts->SetProdVtxRange2D();
  genCuts->SetProdVtxRangeX(0, 1);
  genCuts->SetProdVtxRangeY(0, 1);
  genCuts->SetRequirePdgCode(11, kTRUE);
  fMCsignal->AddCut(genCuts);
  
  AliCFTrackKineCuts *kineMCcuts = new AliCFTrackKineCuts((Char_t *)"fCutsKineMC",(Char_t *)"MC Kine Cuts");
  kineMCcuts->SetPtRange(0.1, 20.);
  kineMCcuts->SetEtaRange(-0.8, 0.8);
  fMCsignal->AddCut(kineMCcuts);

  return fMCsignal;
}

//__________________________________________________________________
void AliHFEtrackFilter::SetMC(AliMCEvent * const mc){
  //
  // Publish MC event to the single cut steps
  //
  fMC = mc;
  fMCsignal->SetMC(fMC);
  AliHFEcutStep *cs = NULL;
  for(Int_t icut = 0; icut < fCutSteps->GetEntriesFast(); icut++)
    if((cs = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))) cs->SetMC(fMC);
}

//__________________________________________________________________
void AliHFEtrackFilter::SetRecEvent(AliVEvent *rec){
  //
  // Publish MC event to the single cut steps
  //
  AliHFEcutStep *cs = NULL;
  for(Int_t icut = 0; icut < fCutSteps->GetEntriesFast(); icut++)
    if((cs = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))) cs->SetRecEvent(rec);
}


//__________________________________________________________________
AliCFContainer *AliHFEtrackFilter::GetEfficiencyContainer(Int_t icont){
  //
  // return EfficiencyContainer
  //
  if(icont >= 4 || icont < 0) return NULL;
  return dynamic_cast<AliCFContainer *>(fEfficiencyContainers->At(icont));
}
Commit	Line	Data
70da6c5a	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	**************************************************************************/
27de2dfb	15
	16	/* $Id$ */
	17
70da6c5a	18	//
	19	// Track filter class
	20	// Apply cut steps to all tracks in one event and returns a list of
	21	// filtered tracks
	22	//
	23	// Author:
	24	// Markus Fasel <M.Fasel@gsi.de>
	25	//
	26	#include <TMath.h>
	27	#include <TString.h>
	28
	29	#include "AliAnalysisCuts.h"
	30	#include "AliCFContainer.h"
	31	#include "AliCFParticleGenCuts.h"
	32	#include "AliCFTrackKineCuts.h"
	33	#include "AliCFTrackIsPrimaryCuts.h"
	34	#include "AliCFTrackQualityCuts.h"
	35	#include "AliESDEvent.h"
	36	#include "AliESDtrack.h"
	37	#include "AliMCEvent.h"
	38	#include "AliMCParticle.h"
	39
	40	#include "AliHFEcontainer.h"
	41	#include "AliHFEcutStep.h"
	42	#include "AliHFEextraCuts.h"
	43	#include "AliHFEtrackFilter.h"
	44	#include "AliHFEtools.h"
	45
	46	ClassImp(AliHFEtrackFilter)
	47
	48	//__________________________________________________________________
	49	AliHFEtrackFilter::AliHFEtrackFilter(const Char_t *name) :
	50	TNamed(name, ""),
	51	fFilteredTracks(NULL),
	52	fCutSteps(NULL),
	53	fEfficiencyContainers(NULL),
	54	fMC(NULL),
	55	fMCsignal(NULL),
	56	fPtBins(0),
	57	fEtaBins(0),
	58	fPhiBins(0),
	59	fPtBinning(NULL),
	60	fEtaBinning(NULL),
	61	fPhiBinning(NULL)
	62	{
	63	//
	64	// Constructor
	65	//
	66	fCutSteps = new TObjArray;
	67	fFilteredTracks = new TObjArray;
	68	fEfficiencyContainers = new TObjArray(4);
	69	}
	70
	71	//__________________________________________________________________
	72	AliHFEtrackFilter::AliHFEtrackFilter(const AliHFEtrackFilter &ref):
	73	TNamed(ref),
	74	fFilteredTracks(NULL),
	75	fCutSteps(NULL),
	76	fEfficiencyContainers(NULL),
	77	fMC(NULL),
	78	fMCsignal(NULL),
	79	fPtBins(0),
	80	fEtaBins(0),
	81	fPhiBins(0),
82	fPtBinning(NULL),
83	fEtaBinning(NULL),
84	fPhiBinning(NULL)
85	{
86	//
87	// Copy constructor
88	//
89	ref.Copy(*this);
90	}
91
92	//__________________________________________________________________
93	AliHFEtrackFilter &AliHFEtrackFilter::operator=(const AliHFEtrackFilter &o){
94	//
95	// Assignment operator
96	//
97	if(this != &o){
98	o.Copy(*this);
99	}
100	return *this;
101	}
102
103	//__________________________________________________________________
104	AliHFEtrackFilter::~AliHFEtrackFilter(){
105	//
106	// Destructor
107	//
108
109	// remove correction framework content
110	if(fPtBinning) delete fPtBinning;
111	if(fEtaBinning) delete fEtaBinning;
112	if(fPhiBinning) delete fPhiBinning;
113	if(TestBit(kOwnCFContainers))
114	fEfficiencyContainers->Delete();
115	delete fEfficiencyContainers;
116	// remove default content
117	delete fFilteredTracks;
118	delete fCutSteps;
119	}
120
121	//__________________________________________________________________
122	void AliHFEtrackFilter::Copy(TObject &o) const{
123	//
124	// Copy content into the object o
125	//
126	TNamed::Copy(o);
127	AliHFEtrackFilter &target = dynamic_cast<AliHFEtrackFilter &>(o);
128
129	// Make copy
130	target.fFilteredTracks = dynamic_cast<TObjArray *>(fFilteredTracks->Clone());
131	target.fCutSteps = dynamic_cast<TObjArray *>(fCutSteps->Clone());
132	target.fEfficiencyContainers = dynamic_cast<TObjArray *>(fEfficiencyContainers->Clone());
133	target.fMC = fMC;
134	target.fMCsignal = dynamic_cast<AliHFEcutStep *>(fMCsignal->Clone());
135
136	target.fPtBins = fPtBins;
137	target.fEtaBins = fEtaBins;
138	target.fPhiBins = fPhiBins;
139	target.fPtBinning = new Double_t[fPtBins];
140	memcpy(target.fPtBinning, fPtBinning, sizeof(Double_t) * fPtBins);
141	target.fEtaBinning = new Double_t[fEtaBins];
142	memcpy(target.fEtaBinning, fEtaBinning, sizeof(Double_t) * fEtaBins);
143	target.fPhiBinning = new Double_t[fPhiBins];
144	memcpy(target.fPhiBinning, fPhiBinning, sizeof(Double_t) * fPhiBins);
145	}
146
147	//__________________________________________________________________
148	void AliHFEtrackFilter::InitCF(){
149	//
150	// Initialize correction framework container
151	// No HFE container given
152	// Only for testing purpose
153	//
154	const Char_t *cnames[4] = {"_container", "_container_signal", "_container_MC", "_container_signalMC"};
155	const Char_t *ctitlesAppends[4] = {"(All Rec)", "(All Signals)", "(All Rec (MC))", "(All Signals(MC))"};
156
157	// Create the binning if not done from outside
158	if(!fPtBins){
159	fPtBins = 40;
160	fPtBinning = AliHFEtools::MakeLogarithmicBinning(fPtBins, 0.1, 10);
161	}
162	if(!fEtaBins){
163	fEtaBins = 8;
164	fEtaBinning = AliHFEtools::MakeLinearBinning(fEtaBins, -0.9, 0.9);
165	}
166	if(!fPhiBins){
167	fPhiBins = 18;
168	fPhiBinning = AliHFEtools::MakeLinearBinning(fPhiBins, 0, 2*TMath::Pi());
169	}
170	Double_t chargeBins[3] = {-1.1, 0., 1.1};
171	Int_t nStep = fCutSteps->GetEntriesFast()+1;
172	Int_t nBins[4] = {fPtBins, fEtaBins, fPhiBins + 2};
173	AliCFContainer *ctmp = NULL;
174	for(Int_t icont = 0; icont < 4; icont++){
175	TString containername = GetName() + TString(cnames[icont]);
176	TString containertitle = TString("Container for filter ") + GetName() + TString(ctitlesAppends[icont]);
177	ctmp = new AliCFContainer(containername.Data(), containertitle.Data(), nStep, 4, nBins);
178	SetBit(kOwnCFContainers, kTRUE);
179
180	// Set the binning
181	ctmp->SetBinLimits(0, fPtBinning);
182	ctmp->SetBinLimits(1, fEtaBinning);
183	ctmp->SetBinLimits(2, fPhiBinning);
184	ctmp->SetBinLimits(3, chargeBins);
185
186	// Label variable names
187	ctmp->SetVarTitle(0, "pt");
188	ctmp->SetVarTitle(1, "eta");
189	ctmp->SetVarTitle(2, "phi");
190	ctmp->SetVarTitle(3, "charge");
191
192	// Label step name
193	ctmp->SetStepTitle(0, "No Cuts");
194	AliHFEcutStep *cutStep = NULL;
195	for(Int_t istep = 0; istep < fCutSteps->GetEntriesFast(); istep++){
196	cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(istep));
bf892a6a	197	if(cutStep) ctmp->SetStepTitle(istep + 1, cutStep->GetName());
70da6c5a	198	}
	199	fEfficiencyContainers->AddAt(ctmp, 0);
	200	}
	201	OwnContainers();
	202	}
	203
	204	//__________________________________________________________________
	205	void AliHFEtrackFilter::InitCF(AliHFEcontainer *cont){
	206	//
	207	// Initialize Correction Framework container
	208	// Use the given HFE container
	209	// HFE standard way to use the correction Framework
	210	//
	211
	212	const Char_t *cnames[4] = {"_container", "_container_signal", "_container_MC", "_container_signalMC"};
	213	const Char_t *ctitlesAppends[4] = {"(All Rec)", "(All Signals)", "(All Rec (MC))", "(All Signals(MC))"};
	214	Int_t nStep = fCutSteps->GetEntriesFast()+1;
	215	AliCFContainer *ctmp = NULL;
	216	for(Int_t icont = 0; icont < 4; icont++){
	217	TString contname = GetName() + TString(cnames[icont]);
	218	TString contTitle = TString("Container for filter ") + GetName() + TString(ctitlesAppends[icont]);
	219	//printf("Adding container %s: %s\n", contname.Data(), contTitle.Data());
	220	cont->CreateContainer(contname.Data(), contTitle.Data(), nStep);
	221	fEfficiencyContainers->AddAt((ctmp = cont->GetCFContainer(contname.Data())), icont);
	222
	223	// Label step name
	224	ctmp->SetStepTitle(0, "No Cuts");
	225	AliHFEcutStep *cutStep = NULL;
	226	for(Int_t istep = 0; istep < fCutSteps->GetEntriesFast(); istep++){
	227	cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(istep));
bf892a6a	228	if(cutStep) ctmp->SetStepTitle(istep + 1, cutStep->GetName());
70da6c5a	229	}
	230	}
	231	ReleaseContainers();
	232	}
	233
	234	//__________________________________________________________________
	235	void AliHFEtrackFilter::FilterTracks(AliESDEvent * const event){
	236	//
	237	// Perform track filtering
	238	// Check each cut step one by one and select tracks which pass
	239	// all cuts. If the correction framework is initialized, a correction#
	240	// framework container will be filled for each cut step.
	241	//
	242	AliHFEcutStep *cutStep = NULL;
	243	AliESDtrack *track = NULL;
	244	AliMCParticle *mctrack = NULL;
	245	Double_t cont[4] = {0., 0., 0., 0.}, contMC[4] = {0., 0., 0., 0.};
	246	Bool_t goodTrack = kTRUE, signal = kFALSE;
	247	Int_t nStep = fCutSteps->GetEntriesFast();
	248	AliCFContainer call = dynamic_cast<AliCFContainer >(fEfficiencyContainers->At(0));
	249	AliCFContainer cSignal = dynamic_cast<AliCFContainer >(fEfficiencyContainers->At(1));
	250	AliCFContainer callMC = dynamic_cast<AliCFContainer >(fEfficiencyContainers->At(2));
	251	AliCFContainer cSignalMC = dynamic_cast<AliCFContainer >(fEfficiencyContainers->At(3));
	252	for(Int_t itrack = 0; itrack < event->GetNumberOfTracks(); itrack++){
	253	signal = kFALSE;
	254	track = event->GetTrack(itrack);
	255	// check Monte-Carlo Information if available
	256	if(fMC){
	257	mctrack = dynamic_cast<AliMCParticle *>(fMC->GetTrack(TMath::Abs(track->GetLabel())));
	258	if(mctrack){
	259	//AliMCParticle mother = dynamic_cast<AliMCParticle >(fMC->GetTrack(mctrack->Particle()->GetFirstMother()));
	260	//AliInfo(Form("Label %d, Mother %d", track->GetLabel(), mother->Particle()->GetPdgCode()));
	261	contMC[0] = mctrack->Pt();
	262	contMC[1] = mctrack->Eta();
	263	contMC[2] = mctrack->Phi();
	264	contMC[3] = mctrack->Charge();
	265	if(fMCsignal->IsSelected(mctrack)) signal = kTRUE;
	266	//if(TMath::Abs(mother->Particle()->GetPdgCode()) != 443) signal = kFALSE;
	267	//AliInfo(Form("Signal? %s", signal ? "Yes": "No"));
	268	}
	269	}
	270
	271	// Fill Array without cut
	272	cont[0] = track->Pt();
	273	cont[1] = track->Eta();
	274	cont[2] = track->Phi();
	275	cont[3] = track->Charge();
	276	if(call) call->Fill(cont, 0);
	277	if(callMC) callMC->Fill(contMC, 0);
	278	if(signal){
	279	if(cSignal) cSignal->Fill(cont, 0);
	280	if(cSignalMC) cSignalMC->Fill(contMC, 0);
	281	}
	282	// cut the track
	283	goodTrack = kTRUE;
	284	for(Int_t icut = 0; icut < nStep; icut++){
	285	cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut));
e3ae862b	286	if(cutStep && (!cutStep->IsSelected(track))){
70da6c5a	287	// track cut away
	288	goodTrack = kFALSE;
	289	break;
	290	}
	291
	292	// Track survived cut step: Fill container
	293	if(call) call->Fill(cont, icut + 1);
	294	if(callMC) callMC->Fill(contMC, icut + 1);
	295	if(signal){
	296	if(cSignal) cSignal->Fill(cont, icut + 1);
	297	if(cSignalMC) cSignalMC->Fill(contMC, icut + 1);
	298	}
	299	}
	300
	301	// Append track to the list of filtered tracks
	302	if(goodTrack) fFilteredTracks->Add(track);
	303	}
	304	}
	305
	306	//__________________________________________________________________
	307	void AliHFEtrackFilter::GenerateCutSteps(){
	308	//
	309	// Make the default cuts
	310	//
	311	MakeCutStepRecKineITSTPC();
	312	MakeCutStepPrimary();
	313	MakeCutStepHFEITS();
	314	MakeCutStepHFETRD();
	315
	316	MakeMCSignalCuts();
	317	}
	318
	319	//__________________________________________________________________
	320	void AliHFEtrackFilter::AddCutStep(AliHFEcutStep *step){
	321	//
	322	// Add new cut step to the filter
	323	//
	324	if(!fCutSteps->FindObject(step->GetName())) fCutSteps->Add(step);
	325	}
	326
	327	//__________________________________________________________________
	328	AliHFEcutStep *AliHFEtrackFilter::GetCutStep(Int_t istep){
	329	//
	330	// Getter for single cut step
	331	//
	332	return dynamic_cast<AliHFEcutStep *>(fCutSteps->At(istep));
	333	}
	334
	335	//__________________________________________________________________
	336	AliHFEcutStep AliHFEtrackFilter::GetCutStep(const Char_t name){
	337	//
	338	// Getter for single cut step (by step name)
	339	//
	340	return dynamic_cast<AliHFEcutStep *>(fCutSteps->FindObject(name));
	341	}
	342
	343	//__________________________________________________________________
	344	void AliHFEtrackFilter::Flush(){
	345	//
	346	// Empty track container
	347	//
	348	fFilteredTracks->Clear();
	349	}
	350
351	//__________________________________________________________________
352	void AliHFEtrackFilter::SetPtBins(Int_t nBins, Double_t *binning){
353	//
354	// User defined pt binning
355	//
356	fPtBins = nBins;
357	fPtBinning = new Double_t[fPtBins + 1];
358	memcpy(fPtBinning, binning, sizeof(Double_t) * nBins);
359	}
360
361	//__________________________________________________________________
362	void AliHFEtrackFilter::SetEtaBins(Int_t nBins, Double_t *binning){
363	//
364	// User defined eta binning
365	//
366	fEtaBins = nBins;
367	fEtaBinning = new Double_t[fEtaBins + 1];
368	memcpy(fEtaBinning, binning, sizeof(Double_t) * nBins);
369	}
370
371	//__________________________________________________________________
372	void AliHFEtrackFilter::SetPhiBins(Int_t nBins, Double_t *binning){
373	//
374	// User defined phi binning
375	//
376	fPhiBins = nBins;
377	fPhiBinning = new Double_t[fPhiBins + 1];
378	memcpy(fPhiBinning, binning, sizeof(Double_t) * nBins);
379	}
380
381	//__________________________________________________________________
382	AliHFEcutStep *AliHFEtrackFilter::MakeCutStepRecKineITSTPC(){
383	//
384	// Make the cut step for Rec Kine
385	// Cut step is already included in the filter
386	//
387	AliHFEcutStep *fCutStep = new AliHFEcutStep("RecKineITSTPC");
388
faee3b18	389	AliCFTrackQualityCuts trackQuality = new AliCFTrackQualityCuts((Char_t )"QualityRec", (Char_t *)"REC Track Quality Cuts");
faee3b18	390	trackQuality->SetMinNClusterTPC(80);
70da6c5a	391	trackQuality->SetMaxChi2PerClusterTPC(3.5);
	392	trackQuality->SetStatus(AliESDtrack::kTPCrefit \| AliESDtrack::kITSrefit);
	393	trackQuality->SetMaxCovDiagonalElements(2., 2., 0.5, 0.5, 2);
	394	fCutStep->AddCut(trackQuality);
	395
	396	AliHFEextraCuts *hfecuts = new AliHFEextraCuts("HFETPC","Extra cuts from the HFE group");
e3ae862b	397	hfecuts->SetClusterRatioTPC(0.6, AliHFEextraCuts::kFoundOverCR);
70da6c5a	398	fCutStep->AddCut(hfecuts);
70da6c5a	399
faee3b18	400	AliCFTrackKineCuts kineCuts = new AliCFTrackKineCuts((Char_t )"RecKine", (Char_t *)"REC Kine Cuts");
	401	kineCuts->SetPtRange(0.1, 20);
	402	kineCuts->SetEtaRange(-0.8, 0.8);
70da6c5a	403	fCutStep->AddCut(kineCuts);
	404
	405	AddCutStep(fCutStep);
	406	return fCutStep;
	407	}
	408
	409	//__________________________________________________________________
	410	AliHFEcutStep *AliHFEtrackFilter::MakeCutStepPrimary(){
	411	//
	412	// Make cut on primaries
	413	// Cut step is already included in the filter
	414	//
	415	AliHFEcutStep *fCutStep = new AliHFEcutStep("Primary");
	416
faee3b18	417	AliCFTrackIsPrimaryCuts primaryCut = new AliCFTrackIsPrimaryCuts((Char_t )"PrimaryCuts", (Char_t *)"REC Primary Cuts");
	418	primaryCut->SetMaxDCAToVertexXY(0.5);
	419	primaryCut->SetMaxDCAToVertexZ(5.);
70da6c5a	420	primaryCut->SetAcceptKinkDaughters(kFALSE);
	421	fCutStep->AddCut(primaryCut);
	422
	423	AddCutStep(fCutStep);
	424	return fCutStep;
	425	}
	426
	427	//__________________________________________________________________
	428	AliHFEcutStep *AliHFEtrackFilter::MakeCutStepHFEITS(){
	429	//
	430	// Add special ITS cuts
	431	// Cut step is already included in the filter
	432	//
	433	AliHFEcutStep *fCutStep = new AliHFEcutStep("HFEITS");
	434
faee3b18	435	AliHFEextraCuts hfecuts = new AliHFEextraCuts((Char_t )"HFEPixelsCuts",(Char_t *)"Extra cuts from the HFE group");
70da6c5a	436	hfecuts->SetRequireITSpixel(AliHFEextraCuts::kFirst);
	437	//hfecuts->SetCheckITSstatus(kTRUE);
	438	fCutStep->AddCut(hfecuts);
	439
	440	AddCutStep(fCutStep);
	441	return fCutStep;
	442	}
	443
	444	//__________________________________________________________________
	445	AliHFEcutStep *AliHFEtrackFilter::MakeCutStepHFETRD(){
	446	//
	447	// Add special TRD cut
	448	// Cut step is already included in the filter
	449	//
	450	AliHFEcutStep *fCutStep = new AliHFEcutStep("HFETRD");
	451
	452	AliHFEextraCuts *hfecuts = new AliHFEextraCuts("HFETRDCuts","Extra cuts from the HFE group");
	453	hfecuts->SetMinTrackletsTRD(0);
	454	fCutStep->AddCut(hfecuts);
	455
	456	AddCutStep(fCutStep);
	457	return fCutStep;
	458	}
	459
	460	//__________________________________________________________________
	461	AliHFEcutStep *AliHFEtrackFilter::MakeMCSignalCuts(){
	462	//
	463	// Define MC Signal
	464	// Cut step is already included in the filter
	465	//
	466	fMCsignal = new AliHFEcutStep("MCSignal");
faee3b18	467	AliCFParticleGenCuts genCuts = new AliCFParticleGenCuts((Char_t )"fCutsGenMC", (Char_t *)"Particle Generation Cuts");
70da6c5a	468	genCuts->SetRequireIsCharged();
70da6c5a	469	genCuts->SetRequireIsPrimary();
faee3b18	470	genCuts->SetProdVtxRange2D();
	471	genCuts->SetProdVtxRangeX(0, 1);
	472	genCuts->SetProdVtxRangeY(0, 1);
70da6c5a	473	genCuts->SetRequirePdgCode(11, kTRUE);
	474	fMCsignal->AddCut(genCuts);
	475
faee3b18	476	AliCFTrackKineCuts kineMCcuts = new AliCFTrackKineCuts((Char_t )"fCutsKineMC",(Char_t *)"MC Kine Cuts");
	477	kineMCcuts->SetPtRange(0.1, 20.);
	478	kineMCcuts->SetEtaRange(-0.8, 0.8);
70da6c5a	479	fMCsignal->AddCut(kineMCcuts);
	480
	481	return fMCsignal;
	482	}
	483
	484	//__________________________________________________________________
	485	void AliHFEtrackFilter::SetMC(AliMCEvent * const mc){
	486	//
	487	// Publish MC event to the single cut steps
	488	//
	489	fMC = mc;
	490	fMCsignal->SetMC(fMC);
bf892a6a	491	AliHFEcutStep *cs = NULL;
70da6c5a	492	for(Int_t icut = 0; icut < fCutSteps->GetEntriesFast(); icut++)
bf892a6a	493	if((cs = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))) cs->SetMC(fMC);
70da6c5a	494	}
	495
	496	//__________________________________________________________________
	497	void AliHFEtrackFilter::SetRecEvent(AliVEvent *rec){
	498	//
	499	// Publish MC event to the single cut steps
	500	//
bf892a6a	501	AliHFEcutStep *cs = NULL;
70da6c5a	502	for(Int_t icut = 0; icut < fCutSteps->GetEntriesFast(); icut++)
bf892a6a	503	if((cs = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))) cs->SetRecEvent(rec);
70da6c5a	504	}
	505
	506
	507	//__________________________________________________________________
	508	AliCFContainer *AliHFEtrackFilter::GetEfficiencyContainer(Int_t icont){
	509	//
	510	// return EfficiencyContainer
	511	//
	512	if(icont >= 4 \|\| icont < 0) return NULL;
	513	return dynamic_cast<AliCFContainer *>(fEfficiencyContainers->At(icont));
	514	}