]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG1/AliComparisonDCA.cxx
short description of the package from Naomi's presentation
[u/mrichter/AliRoot.git] / PWG1 / AliComparisonDCA.cxx
index fe5cec011a1cd2198a94747e7663cd311cb3c1c5..24f3ac61d59a51d0b6fcc563638a11c5d8ebdb80 100644 (file)
-//------------------------------------------------------------------------------\r
-// Implementation of AliComparisonDCA class. It keeps information from \r
-// comparison of reconstructed and MC particle tracks. In addtion, \r
-// it keeps selection cuts used during comparison. The comparison \r
-// information is stored in the ROOT histograms. Analysis of these \r
-// histograms can be done by using Analyse() class function. The result of \r
-// the analysis (histograms) are stored in the output picture_dca.root file.\r
-//  \r
-// Author: J.Otwinowski 04/02/2008 \r
-//------------------------------------------------------------------------------\r
-\r
-/*\r
-  // after running analysis, read the file, and get component\r
-  gSystem->Load("libPWG1.so");\r
-  TFile f("Output.root");\r
-  AliComparisonDCA * comp = (AliComparisonDCA*)f.Get("AliComparisonDCA");\r
-\r
-  // analyse comparison data (output stored in pictures_dca.root)\r
-  comp->Analyse();\r
-\r
-  // TPC track length parameterisation\r
-  TF1 fl("fl","((min(250./(abs(x+0.000001)),250)-90))",0,2);  // length function\r
-  TF1 fl2("fl2","[0]/((min(250./(abs(x+0.000001)),250)-90))^[1]",0,2);\r
-  fl2.SetParameter(1,1);\r
-  fl2.SetParameter(0,1);\r
-*/\r
-\r
-#include <iostream>\r
-\r
-#include "TFile.h"\r
-#include "TCint.h"\r
-#include "TH3F.h"\r
-#include "TH2F.h"\r
-#include "TF1.h"\r
-#include "TProfile.h"\r
-#include "TProfile2D.h"\r
-#include "TGraph2D.h"\r
-#include "TCanvas.h"\r
-#include "TGraph.h"\r
-// \r
-#include "AliTracker.h"   \r
-#include "AliESDEvent.h"   \r
-#include "AliESD.h"\r
-#include "AliESDfriend.h"\r
-#include "AliESDfriendTrack.h"\r
-#include "AliRecInfoCuts.h" \r
-#include "AliMCInfoCuts.h" \r
-#include "AliLog.h" \r
-#include "AliESDVertex.h" \r
-//\r
-#include "AliMathBase.h"\r
-#include "AliTreeDraw.h" \r
-\r
-#include "AliMCInfo.h" \r
-#include "AliESDRecInfo.h" \r
-#include "AliComparisonDCA.h" \r
-\r
-using namespace std;\r
-\r
-ClassImp(AliComparisonDCA)\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonDCA::AliComparisonDCA():\r
-  TNamed("AliComparisonDCA","AliComparisonDCA"),\r
-\r
-  // DCA histograms\r
-  fD0TanSPtB1(0),\r
-  fD1TanSPtB1(0),\r
-  fD0TanSPtL1(0),\r
-  fD1TanSPtL1(0),\r
-  fD0TanSPtInTPC(0),\r
-  fD1TanSPtInTPC(0),\r
-  fVertex(0),\r
-\r
-  // Cuts \r
-  fCutsRC(0), \r
-  fCutsMC(0)  \r
-{\r
-  InitHisto();\r
-  InitCuts();\r
-\r
-  // vertex (0,0,0)\r
-  fVertex = new AliESDVertex();\r
-  fVertex->SetXv(0.0);\r
-  fVertex->SetYv(0.0);\r
-  fVertex->SetZv(0.0);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonDCA::~AliComparisonDCA()\r
-{\r
-  //\r
-  if(fD0TanSPtB1) delete fD0TanSPtB1; fD0TanSPtB1=0;\r
-  if(fD1TanSPtB1) delete fD1TanSPtB1; fD1TanSPtB1=0;\r
-  if(fD0TanSPtL1) delete fD0TanSPtL1; fD0TanSPtL1=0;\r
-  if(fD1TanSPtL1) delete fD1TanSPtL1; fD1TanSPtL1=0;\r
-  if(fD0TanSPtInTPC) delete fD0TanSPtInTPC; fD0TanSPtInTPC=0;\r
-  if(fD1TanSPtInTPC) delete fD1TanSPtInTPC; fD1TanSPtInTPC=0;\r
-  if(fVertex) delete fVertex; fVertex=0;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonDCA::InitHisto()\r
-{\r
-  // DCA histograms\r
-  fD0TanSPtB1 = new TH3F("DCAyTanSPtB1","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);\r
-  fD0TanSPtB1->SetXTitle("tan(#theta)");\r
-  fD0TanSPtB1->SetYTitle("#sqrt{p_{t}(GeV/c)}");\r
-  fD0TanSPtB1->SetZTitle("DCA_{xy}");\r
-\r
-  fD1TanSPtB1 = new TH3F("DCAzTanSPtB1","DCAzTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);\r
-  fD1TanSPtB1->SetXTitle("tan(#theta)");\r
-  fD1TanSPtB1->SetYTitle("#sqrt(p_{t}(GeV/c))");\r
-  fD1TanSPtB1->SetZTitle("DCA_{z}");\r
-\r
-  fD0TanSPtL1 = new TH3F("DCAyTanSPtL1","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);\r
-  fD0TanSPtL1->SetXTitle("tan(#theta)");\r
-  fD0TanSPtL1->SetYTitle("#sqrt{p_{t}(GeV/c)}");\r
-  fD0TanSPtL1->SetZTitle("DCA_{xy}");\r
-\r
-  fD1TanSPtL1 = new TH3F("DCAzTanSPtL1","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);\r
-  fD1TanSPtL1->SetXTitle("tan(#theta)");\r
-  fD1TanSPtL1->SetYTitle("#sqrt{p_{t}(GeV/c)}");\r
-  fD1TanSPtL1->SetZTitle("DCA_{z}");\r
-\r
-  fD0TanSPtInTPC = new TH3F("DCAyTanSPtInTPC","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);\r
-  fD0TanSPtInTPC->SetXTitle("tan(#theta)");\r
-  fD0TanSPtInTPC->SetYTitle("#sqrt{p_{t}(GeV/c)}");\r
-  fD0TanSPtInTPC->SetZTitle("DCA_{xy}");\r
-\r
-  fD1TanSPtInTPC = new TH3F("DCAzTanSPtInTPC","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);\r
-  fD1TanSPtInTPC->SetXTitle("tan(#theta)");\r
-  fD1TanSPtInTPC->SetYTitle("#sqrt{p_{t}(GeV/c)}");\r
-  fD1TanSPtInTPC->SetZTitle("DCA_{z}");\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonDCA::InitCuts()\r
-{\r
-  if(!fCutsMC) \r
-    AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");\r
-  if(!fCutsRC) \r
-    AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonDCA::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC)\r
-{\r
-  // Fill DCA comparison information\r
-  AliExternalTrackParam *track = 0;\r
-  Double_t kRadius    = 3.0;      // beam pipe radius\r
-  Double_t kMaxStep   = 5.0;      // max step\r
-  Double_t field      = AliTracker::GetBz(); // nominal Bz field [kG]\r
-  Double_t kMaxD      = 123456.0; // max distance\r
-\r
-  Int_t clusterITS[200];\r
-  Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
-\r
-  Float_t mcpt = infoMC->GetParticle().Pt();\r
-  Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5);\r
-  Float_t spt = TMath::Sqrt(mcpt);\r
-\r
-  // distance to Prim. vertex \r
-  const Double_t* dv = infoMC->GetVDist(); \r
-\r
-  Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
-\r
-  // Check selection cuts\r
-  if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
-  if (!isPrim) return;\r
-  if (infoRC->GetStatus(1)!=3) return;\r
-  if (!infoRC->GetESDtrack()) return;  \r
-  if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;\r
-  if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;\r
-\r
-  // calculate and set prim. vertex\r
-  fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );\r
-  fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );\r
-  fVertex->SetZv( infoMC->GetParticle().Vz() - dv[2] );\r
-\r
-  // calculate track parameters at vertex\r
-  if (infoRC->GetESDtrack()->GetTPCInnerParam())\r
-  {\r
-    if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )\r
-    {\r
-      Bool_t bStatus = AliTracker::PropagateTrackTo(track,kRadius,infoMC->GetMass(),kMaxStep,kTRUE);\r
-      Bool_t bDCAStatus = track->PropagateToDCA(fVertex,field,kMaxD,dca,cov);\r
-\r
-      if(bStatus && bDCAStatus) {\r
-        fD0TanSPtInTPC->Fill(tantheta,spt,dca[0]);\r
-        fD1TanSPtInTPC->Fill(tantheta,spt,dca[1]);\r
-         }\r
-\r
-         delete track;\r
-    }\r
-  }\r
-\r
- if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)==0){\r
-    fD0TanSPtB1->Fill(tantheta,spt,dca[0]);\r
-    fD1TanSPtB1->Fill(tantheta,spt,dca[1]);\r
-  }\r
-    fD0TanSPtL1->Fill(tantheta,spt,dca[0]);\r
-    fD1TanSPtL1->Fill(tantheta,spt,dca[1]);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-Long64_t AliComparisonDCA::Merge(TCollection* list) \r
-{\r
-  // Merge list of objects (needed by PROOF)\r
-\r
-  if (!list)\r
-  return 0;\r
-\r
-  if (list->IsEmpty())\r
-  return 1;\r
-\r
-  TIterator* iter = list->MakeIterator();\r
-  TObject* obj = 0;\r
-\r
-  // collection of generated histograms\r
-  Int_t count=0;\r
-  while((obj = iter->Next()) != 0) \r
-  {\r
-    AliComparisonDCA* entry = dynamic_cast<AliComparisonDCA*>(obj);\r
-    if (entry == 0) { \r
-      Error("Add","Attempt to add object of class: %s to a %s",\r
-         obj->ClassName(),this->ClassName());\r
-      return -1;\r
-    }\r
-\r
-    fD0TanSPtB1->Add(entry->fD0TanSPtB1);\r
-    fD1TanSPtB1->Add(entry->fD1TanSPtB1);\r
-    fD0TanSPtL1->Add(entry->fD0TanSPtL1);\r
-    fD1TanSPtL1->Add(entry->fD1TanSPtL1);\r
-    fD0TanSPtInTPC->Add(entry->fD0TanSPtInTPC);\r
-    fD1TanSPtInTPC->Add(entry->fD1TanSPtInTPC);\r
-\r
-    count++;\r
-  }\r
-\r
-return count;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonDCA::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC){\r
-  // Process comparison information\r
-  Process(infoMC,infoRC);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonDCA::Analyse()\r
-{\r
-  // Analyse output histograms\r
-  \r
-  AliComparisonDCA * comp=this;\r
-  TGraph2D * gr=0;\r
-  TGraph * gr0=0;\r
-\r
-  TFile *fp = new TFile("pictures_dca.root","recreate");\r
-  fp->cd();\r
-\r
-  TCanvas * c = new TCanvas("DCA","DCA resloution");\r
-  c->cd();\r
-\r
-  // DCA resolution\r
-  gr0 = AliMathBase::MakeStat1D(comp->fD0TanSPtB1,2,5);\r
-  gr0->GetXaxis()->SetTitle("Tan(#theta)");\r
-  gr0->GetYaxis()->SetTitle("#sigmaDCA (cm)");\r
-  gr0->Write("DCAResolTan");\r
-  //\r
-  gr = AliMathBase::MakeStat2D(comp->fD0TanSPtB1,4,2,5); \r
-  gr->GetXaxis()->SetTitle("Tan(#theta)");\r
-  gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");\r
-  gr->GetYaxis()->SetTitle("#sigmaDCA (cm)");\r
-  gr->GetHistogram()->Write("DCAResolSPTTan");\r
-\r
-  gPad->Clear();\r
-  gr0->Draw("al*");\r
-\r
-  fp->Close();\r
-}\r
+//------------------------------------------------------------------------------
+// Implementation of AliComparisonDCA class. It keeps information from 
+// comparison of reconstructed and MC particle tracks. In addtion, 
+// it keeps selection cuts used during comparison. The comparison 
+// information is stored in the ROOT histograms. Analysis of these 
+// histograms can be done by using Analyse() class function. The result of 
+// the analysis (histograms/graphs) are stored in the folder
+// which is a data member of AliComparisonDCA.
+//  
+// Author: J.Otwinowski 04/02/2008 
+//------------------------------------------------------------------------------
+
+/*
+  // after running comparison task, read the file, and get component
+  gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
+  LoadMyLibs();
+  TFile f("Output.root");
+  //AliComparisonDCA * compObj = (AliComparisonDCA*)f.Get("AliComparisonDCA");
+  AliComparisonDCA * compObj = (AliComparisonDCA*)cOutput->FindObject("AliComparisonDCA");
+
+  // Analyse comparison data
+  compObj->Analyse();
+
+  // the output histograms/graphs will be stored in the folder "folderDCA" 
+  compObj->GetAnalysisFolder()->ls("*");
+  // user can save whole comparison object (or only folder with anlysed histograms) 
+  // in the seperate output file (e.g.)
+  TFile fout("Analysed_DCA.root","recreate");
+  compObj->Write(); // compObj->GetAnalysisFolder()->Write();
+  fout.Close();
+
+*/
+
+#include <TAxis.h>
+#include <TCanvas.h>
+#include <TGraph.h>
+#include <TGraph2D.h>
+#include <TH1.h>
+
+#include "AliComparisonDCA.h" 
+#include "AliESDEvent.h"   
+#include "AliESDRecInfo.h" 
+#include "AliESDVertex.h" 
+#include "AliLog.h" 
+#include "AliMCInfo.h" 
+#include "AliMCInfoCuts.h" 
+#include "AliMathBase.h"
+#include "AliRecInfoCuts.h" 
+#include "AliTracker.h"   
+
+using namespace std;
+
+ClassImp(AliComparisonDCA)
+
+//_____________________________________________________________________________
+AliComparisonDCA::AliComparisonDCA():
+  AliComparisonObject("AliComparisonDCA"),
+
+  // DCA histograms
+  fDCAHisto(0),
+  /*
+  fD0TanSPtTPCITS(0),
+  fD1TanSPtTPCITS(0),
+  fD0TanSPt(0),
+  fD1TanSPt(0),
+  fD0TanSPtTPC(0),
+  fD1TanSPtTPC(0),
+  */
+
+  // Cuts 
+  fCutsRC(0), 
+  fCutsMC(0),  
+
+  // histogram folder 
+  fAnalysisFolder(0)
+{
+  // default constructor       
+}
+
+//_____________________________________________________________________________
+AliComparisonDCA::AliComparisonDCA(Char_t* name="AliComparisonDCA", Char_t* title="AliComparisonDCA",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):
+  AliComparisonObject(name,title),
+
+  // DCA histograms
+  fDCAHisto(0),
+  /*
+  fD0TanSPtTPCITS(0),
+  fD1TanSPtTPCITS(0),
+  fD0TanSPt(0),
+  fD1TanSPt(0),
+  fD0TanSPtTPC(0),
+  fD1TanSPtTPC(0),
+  */
+
+  // Cuts 
+  fCutsRC(0), 
+  fCutsMC(0),  
+
+  // histogram folder 
+  fAnalysisFolder(0)
+{
+  // named constructor  
+
+  SetAnalysisMode(analysisMode);
+  SetHptGenerator(hptGenerator);
+  Init();
+}
+
+
+//_____________________________________________________________________________
+AliComparisonDCA::~AliComparisonDCA()
+{
+  // destructor
+  if(fDCAHisto)  delete fDCAHisto; fDCAHisto=0; 
+  /*
+  if(fD0TanSPtTPCITS) delete fD0TanSPtTPCITS; fD0TanSPtTPCITS=0;
+  if(fD1TanSPtTPCITS) delete fD1TanSPtTPCITS; fD1TanSPtTPCITS=0;
+  if(fD0TanSPt) delete fD0TanSPt; fD0TanSPt=0;
+  if(fD1TanSPt) delete fD1TanSPt; fD1TanSPt=0;
+  if(fD0TanSPtTPC) delete fD0TanSPtTPC; fD0TanSPtTPC=0;
+  if(fD1TanSPtTPC) delete fD1TanSPtTPC; fD1TanSPtTPC=0;
+  */
+  if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
+
+}
+
+//_____________________________________________________________________________
+void AliComparisonDCA::Init()
+{
+  // DCA histograms
+
+ Int_t nPBins = 31;
+    Double_t binsP[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
+    Double_t pMin = 0., pMax = 10.;
+
+    if(IsHptGenerator() == kTRUE) {
+      nPBins = 100;
+      pMin = 0.; pMax = 100.;
+    }
+
+   //dca_r, dca_z, eta, pt
+   Int_t binsQA[4]    = {100,100,20,nPBins};
+   Double_t xminQA[4] = {-10.,-10.,-1., pMin};
+   Double_t xmaxQA[4] = {10.,10.,1., pMax};
+
+   fDCAHisto = new THnSparseF("fDCAHisto","dca_r:dca_z:eta:pt",4,binsQA,xminQA,xmaxQA);
+   if(!IsHptGenerator()) fDCAHisto->SetBinEdges(3,binsP);
+
+   fDCAHisto->GetAxis(0)->SetTitle("dca_r (cm)");
+   fDCAHisto->GetAxis(1)->SetTitle("dca_z (cm)");
+   fDCAHisto->GetAxis(2)->SetTitle("eta");
+   fDCAHisto->GetAxis(3)->SetTitle("pt (GeV/c)");
+   fDCAHisto->Sumw2();
+       
+  /*   
+  fD0TanSPtTPCITS = new TH3F("DCAyTanSPtTPCITS","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
+  fD0TanSPtTPCITS->SetXTitle("tan(#theta)");
+  fD0TanSPtTPCITS->SetYTitle("#sqrt{p_{t}(GeV)}");
+  fD0TanSPtTPCITS->SetZTitle("DCA_{xy}");
+
+  fD1TanSPtTPCITS = new TH3F("DCAzTanSPtTPCITS","DCAzTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
+  fD1TanSPtTPCITS->SetXTitle("tan(#theta)");
+  fD1TanSPtTPCITS->SetYTitle("#sqrt(p_{t}(GeV))");
+  fD1TanSPtTPCITS->SetZTitle("DCA_{z}");
+
+  fD0TanSPt = new TH3F("DCAyTanSPt","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
+  fD0TanSPt->SetXTitle("tan(#theta)");
+  fD0TanSPt->SetYTitle("#sqrt{p_{t}(GeV)}");
+  fD0TanSPt->SetZTitle("DCA_{xy}");
+
+  fD1TanSPt = new TH3F("DCAzTanSPt","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
+  fD1TanSPt->SetXTitle("tan(#theta)");
+  fD1TanSPt->SetYTitle("#sqrt{p_{t}(GeV)}");
+  fD1TanSPt->SetZTitle("DCA_{z}");
+
+  fD0TanSPtTPC = new TH3F("DCAyTanSPtTPC","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
+  fD0TanSPtTPC->SetXTitle("tan(#theta)");
+  fD0TanSPtTPC->SetYTitle("#sqrt{p_{t}(GeV)}");
+  fD0TanSPtTPC->SetZTitle("DCA_{xy}");
+
+  fD1TanSPtTPC = new TH3F("DCAzTanSPtTPC","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
+  fD1TanSPtTPC->SetXTitle("tan(#theta)");
+  fD1TanSPtTPC->SetYTitle("#sqrt{p_{t}(GeV)}");
+  fD1TanSPtTPC->SetZTitle("DCA_{z}");
+  */
+
+  // init cuts
+  if(!fCutsMC) 
+    AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
+  if(!fCutsRC) 
+    AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
+  // init folder
+  fAnalysisFolder = CreateFolder("folderDCA","Analysis DCA Folder");
+}
+
+//_____________________________________________________________________________
+void AliComparisonDCA::ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo * const infoRC)
+{
+  // Fill DCA comparison information
+  AliExternalTrackParam *track = 0;
+  //Double_t field      = AliTracker::GetBz(); // nominal Bz field [kG]
+  Double_t kMaxD      = 123456.0; // max distance
+
+  Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
+
+  Float_t mcpt = infoMC->GetParticle().Pt();
+  Float_t mceta = infoMC->GetParticle().Eta();
+  //Float_t spt = TMath::Sqrt(mcpt);
+
+  // distance to Prim. vertex 
+  const Double_t* dv = infoMC->GetVDist(); 
+
+  Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
+
+  // Check selection cuts
+  if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; 
+  if (!isPrim) return;
+  if (infoRC->GetStatus(1)!=3) return;
+  if (!infoRC->GetESDtrack()) return;  
+  if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
+  //if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;
+
+  // calculate and set prim. vertex
+  AliESDVertex vertexMC;
+  vertexMC.SetXv( infoMC->GetParticle().Vx() - dv[0] );
+  vertexMC.SetYv( infoMC->GetParticle().Vy() - dv[1] );
+  vertexMC.SetZv( infoMC->GetParticle().Vz() - dv[2] );
+
+  // calculate track parameters at vertex
+  if (infoRC->GetESDtrack()->GetTPCInnerParam())
+  {
+    if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )
+    {
+      
+      //Bool_t bDCAStatus = track->PropagateToDCA(&vertexMC,field,kMaxD,dca,cov);
+      Double_t field[3];  track->GetBxByBz(field); 
+      Bool_t bDCAStatus = track->PropagateToDCABxByBz(&vertexMC,field,kMaxD,dca,cov);
+
+      if(bDCAStatus) {
+        Double_t vDCAHisto[4]={dca[0],dca[1],mceta,mcpt};
+        fDCAHisto->Fill(vDCAHisto);
+      }
+    delete track;
+    }
+  }
+}
+
+//_____________________________________________________________________________
+void AliComparisonDCA::ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo * const infoRC)
+{
+  // Fill DCA comparison information
+  Int_t clusterITS[200];
+  Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
+
+  Float_t mcpt = infoMC->GetParticle().Pt();
+  Float_t mceta = infoMC->GetParticle().Eta();
+  //Float_t spt = TMath::Sqrt(mcpt);
+
+  // distance to Prim. vertex 
+  const Double_t* dv = infoMC->GetVDist(); 
+  Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
+
+  // Check selection cuts
+  if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; 
+  if (!isPrim) return;
+  if (infoRC->GetStatus(1)!=3) return;
+  if (!infoRC->GetESDtrack()) return;  
+  if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
+  //if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;
+
+  infoRC->GetESDtrack()->GetImpactParameters(dca,cov);
+
+  // ITS + TPC
+  if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)>fCutsRC->GetMinNClustersITS())
+  {
+    Double_t vDCAHisto[4]={dca[0],dca[1],mceta,mcpt};
+    fDCAHisto->Fill(vDCAHisto);
+  }
+}
+
+void AliComparisonDCA::ProcessConstrained(AliMCInfo* const /*infoMC*/, AliESDRecInfo * const /*infoRC*/)
+{
+  // Fill DCA comparison information
+  
+  AliDebug(AliLog::kWarning, "Warning: Not implemented");
+}
+
+//_____________________________________________________________________________
+Long64_t AliComparisonDCA::Merge(TCollection* const list) 
+{
+  // Merge list of objects (needed by PROOF)
+
+  if (!list)
+  return 0;
+
+  if (list->IsEmpty())
+  return 1;
+
+  TIterator* iter = list->MakeIterator();
+  TObject* obj = 0;
+
+  // collection of generated histograms
+  Int_t count=0;
+  while((obj = iter->Next()) != 0) 
+  {
+    AliComparisonDCA* entry = dynamic_cast<AliComparisonDCA*>(obj);
+    if (entry == 0) continue; 
+
+    fDCAHisto->Add(entry->fDCAHisto);
+    /*
+    fD0TanSPtTPCITS->Add(entry->fD0TanSPtTPCITS);
+    fD1TanSPtTPCITS->Add(entry->fD1TanSPtTPCITS);
+    fD0TanSPt->Add(entry->fD0TanSPt);
+    fD1TanSPt->Add(entry->fD1TanSPt);
+    fD0TanSPtTPC->Add(entry->fD0TanSPtTPC);
+    fD1TanSPtTPC->Add(entry->fD1TanSPtTPC);
+    */
+
+    count++;
+  }
+
+return count;
+}
+
+//_____________________________________________________________________________
+void AliComparisonDCA::Exec(AliMCInfo* const infoMC, AliESDRecInfo * const infoRC)
+{
+  // Process comparison information
+  if(GetAnalysisMode() == 0) ProcessTPC(infoMC,infoRC);
+  else if(GetAnalysisMode() == 1) ProcessTPCITS(infoMC,infoRC);
+  else if(GetAnalysisMode() == 2) ProcessConstrained(infoMC,infoRC);
+  else {
+    printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
+    return;
+  }
+}
+
+//_____________________________________________________________________________
+void AliComparisonDCA::Analyse()
+{
+  //
+  // Analyse comparison information and store output histograms
+  // in the analysis folder "folderDCA" 
+  //
+  
+  TH1::AddDirectory(kFALSE);
+  TObjArray *aFolderObj = new TObjArray;
+
+  /*
+  TGraph * gr[8]= { 0,0,0,0,0,0,0,0 };
+  TGraph2D *gr2[8]= { 0,0,0,0,0,0,0,0};
+  AliComparisonDCA * comp=this;
+
+  // write results in the folder 
+  // Canvas to draw analysed histograms
+  TCanvas * c = new TCanvas("canDCA","DCA resolution");
+  c->Divide(2,4);
+  //
+  // DCA resolution
+  //
+  c->cd(1);
+  gr[0] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPC,2,5);
+  gr[0]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[0]->GetYaxis()->SetTitle("#sigmaDCA_xy (cm)");
+  gr[0]->SetName("DCAXYResolTanTPC");
+  gr[0]->SetTitle("resol. DCA_xy (TPC only)");
+  gr[0]->Draw("Al*");
+
+  aFolderObj->Add(gr[0]);
+
+  c->cd(2);
+  gr[1] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPC,2,5);
+  gr[1]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[1]->GetYaxis()->SetTitle("#sigmaDCA_z (cm)");
+  gr[1]->SetName("DCAZResolTanTPC");
+  gr[1]->SetTitle("resol. DCA_z (TPC only)");
+  gr[1]->Draw("Al*");
+
+  aFolderObj->Add(gr[1]);
+
+  c->cd(3);
+  gr[2] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPCITS,2,5);
+  gr[2]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[2]->GetYaxis()->SetTitle("#sigmaDCA_xy (cm)");
+  gr[2]->SetName("DCAXYResolTanTPCITS");
+  gr[2]->SetTitle("resol. DCA_xy (TPC+ITS)");
+  gr[2]->Draw("Al*");
+
+  aFolderObj->Add(gr[2]);
+
+  c->cd(4);
+  gr[3] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPCITS,2,5);
+  gr[3]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[3]->GetYaxis()->SetTitle("#sigmaDCA_z (cm)");
+  gr[3]->SetName("DCAZResolTanTPCITS");
+  gr[3]->SetTitle("resol. DCA_z (TPC+ITS)");
+  gr[3]->Draw("Al*");
+
+  aFolderObj->Add(gr[3]);
+
+  //
+  // DCA mean value
+  //
+  c->cd(5);
+  gr[4] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPC,2,4);
+  gr[4]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[4]->GetYaxis()->SetTitle("mean DCA_xy (cm)");
+  gr[4]->SetName("DCAXYMeanTanTPC");
+  gr[4]->SetTitle("mean DCA_xy (TPC only)");
+  gr[4]->Draw("Al*");
+
+  aFolderObj->Add(gr[4]);
+
+  c->cd(6);
+  gr[5] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPC,2,4);
+  gr[5]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[5]->GetYaxis()->SetTitle("mean DCA_z (cm)");
+  gr[5]->SetName("DCAZMeanTanTPC");
+  gr[5]->SetTitle("mean DCA_z (TPC only)");
+  gr[5]->Draw("Al*");
+
+  aFolderObj->Add(gr[5]);
+
+  c->cd(7);
+  gr[6] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPCITS,2,4);
+  gr[6]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[6]->GetYaxis()->SetTitle("mean DCA_xy (cm)");
+  gr[6]->SetName("DCAXYMeanTanTPCITS");
+  gr[6]->SetTitle("mean DCA_xy (TPC+ITS)");
+  gr[6]->Draw("Al*");
+
+  aFolderObj->Add(gr[6]);
+
+  c->cd(8);
+  gr[7] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPCITS,2,4);
+  gr[7]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr[7]->GetYaxis()->SetTitle("mean DCA_z (cm)");
+  gr[7]->SetName("DCAZMeanTanTPCITS");
+  gr[7]->SetTitle("mean DCA_z (TPC+ITS)");
+  gr[7]->Draw("Al*");
+
+  aFolderObj->Add(gr[7]);
+
+  // 2D DCA resolution 
+  TCanvas * c1 = new TCanvas("canDCA1","2D DCA resolution");
+  c1->Divide(2,4);
+
+  // TPC only
+  c1->cd(1);
+  gr2[0] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPC,4,2,5); 
+  gr2[0]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[0]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[0]->GetZaxis()->SetTitle("#sigmaDCA_xy (cm)");
+  gr2[0]->SetName("DCAXYResolSPTTanTPC");
+  gr2[0]->SetTitle("#sigma DCA_xy (TPC only)");
+  gr2[0]->GetHistogram()->Draw("colz");
+
+  gr2[0]->GetHistogram()->SetName("DCAXYResolSPTTanTPC");
+  aFolderObj->Add(gr2[0]->GetHistogram());
+
+  c1->cd(2);
+  gr2[1] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPC,4,2,5); 
+  gr2[1]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[1]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[1]->GetZaxis()->SetTitle("#sigmaDCA_z (cm)");
+  gr2[1]->SetName("DCAZResolSPTTanTPC");
+  gr2[1]->SetTitle("#sigma DCA_z (TPC only)");
+  gr2[1]->GetHistogram()->Draw("colz");
+
+  gr2[1]->GetHistogram()->SetName("DCAZResolSPTTanTPC");
+  aFolderObj->Add(gr2[1]->GetHistogram());
+
+  // TPC+ITS
+  c1->cd(3);
+  gr2[2] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPCITS,4,2,5); 
+  gr2[2]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[2]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[2]->GetZaxis()->SetTitle("#sigmaDCA_xy (cm)");
+  gr2[2]->SetName("DCAXYResolSPTTanTPCITS");
+  gr2[2]->SetTitle("#sigma DCA_xy (TPC+ITS)");
+  gr2[2]->GetHistogram()->Draw("colz");
+
+  gr2[2]->GetHistogram()->SetName("DCAXYResolSPTTanTPCITS");
+  aFolderObj->Add(gr2[2]->GetHistogram());
+
+  c1->cd(4);
+  gr2[3] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPCITS,4,2,5); 
+  gr2[3]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[3]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[3]->GetZaxis()->SetTitle("#sigmaDCA_z (cm)");
+  gr2[3]->SetName("DCAZResolSPTTanTPCITS");
+  gr2[3]->SetTitle("#sigma DCA_z (TPC+ITS)");
+  gr2[3]->GetHistogram()->Draw("colz");
+
+  gr2[3]->GetHistogram()->SetName("DCAZResolSPTTanTPCITS");
+  aFolderObj->Add(gr2[3]->GetHistogram());
+
+  // 2D DCA mean value  
+  c1->cd(5);
+  gr2[4] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPC,4,2,4); 
+  gr2[4]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[4]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[4]->GetZaxis()->SetTitle("mean DCA_xy (cm)");
+  gr2[4]->SetName("DCAXYMeanSPTTanTPC");
+  gr2[4]->SetTitle("mean DCA_xy (TPC only)");
+  gr2[4]->GetHistogram()->Draw("colz");
+
+  gr2[4]->GetHistogram()->SetName("DCAXYMeanSPTTanTPC");
+  aFolderObj->Add(gr2[4]->GetHistogram());
+
+  c1->cd(6);
+  gr2[5] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPC,4,2,4); 
+  gr2[5]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[5]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[5]->GetZaxis()->SetTitle("mean DCA_z (cm)");
+  gr2[5]->SetName("DCAZMeanSPTTanTPC");
+  gr2[5]->SetTitle("mean DCA_z (TPC only)");
+  gr2[5]->GetHistogram()->Draw("colz");
+
+  gr2[5]->GetHistogram()->SetName("DCAZMeanSPTTanTPC");
+  aFolderObj->Add(gr2[5]->GetHistogram());
+
+  c1->cd(7);
+  gr2[6] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPCITS,4,2,4); 
+  gr2[6]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[6]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[6]->GetZaxis()->SetTitle("mean DCA_xy (cm)");
+  gr2[6]->SetName("DCAXYMeanSPTTanTPCITS");
+  gr2[6]->SetTitle("mean DCA_xy (TPC+ITS)");
+  gr2[6]->GetHistogram()->Draw("colz");
+
+  gr2[6]->GetHistogram()->SetName("DCAXYMeanSPTTanTPCITS");
+  aFolderObj->Add(gr2[6]->GetHistogram());
+
+  c1->cd(8);
+  gr2[7] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPCITS,4,2,4); 
+  gr2[7]->GetXaxis()->SetTitle("Tan(#theta)");
+  gr2[7]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
+  gr2[7]->GetZaxis()->SetTitle("mean DCA_z (cm)");
+  gr2[7]->SetName("DCAZMeanSPTTanTPCITS");
+  gr2[7]->SetTitle("mean DCA_z (TPC+ITS)");
+  gr2[7]->GetHistogram()->Draw("colz");
+
+  gr2[7]->GetHistogram()->SetName("DCAZMeanSPTTanTPCITS");
+  aFolderObj->Add(gr2[7]->GetHistogram());
+
+  */
+  // export objects to analysis folder
+  fAnalysisFolder = ExportToFolder(aFolderObj);
+
+  // delete only TObjArray
+  if(aFolderObj) delete aFolderObj;
+}
+
+//_____________________________________________________________________________
+TFolder* AliComparisonDCA::ExportToFolder(TObjArray * array) 
+{
+  // recreate folder avery time and export objects to new one
+  //
+  AliComparisonDCA * comp=this;
+  TFolder *folder = comp->GetAnalysisFolder();
+
+  TString name, title;
+  TFolder *newFolder = 0;
+  Int_t i = 0;
+  Int_t size = array->GetSize();
+
+  if(folder) { 
+     // get name and title from old folder
+     name = folder->GetName();  
+     title = folder->GetTitle();  
+
+        // delete old one
+     delete folder;
+
+        // create new one
+     newFolder = CreateFolder(name.Data(),title.Data());
+     newFolder->SetOwner();
+
+        // add objects to folder
+     while(i < size) {
+          newFolder->Add(array->At(i));
+          i++;
+        }
+  }
+
+return newFolder;
+}
+
+
+//_____________________________________________________________________________
+TFolder* AliComparisonDCA::CreateFolder(TString name,TString title) { 
+// create folder for analysed histograms
+TFolder *folder = 0;
+  folder = new TFolder(name.Data(),title.Data());
+
+  return folder;
+}