[u/mrichter/AliRoot.git] / PWG1 / AliComparisonDCA.cxx

//------------------------------------------------------------------------------
// 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");

  // 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 <iostream>

#include "TFile.h"
#include "TCint.h"
#include "TH3F.h"
#include "TH2F.h"
#include "TF1.h"
#include "TProfile.h"
#include "TProfile2D.h"
#include "TGraph2D.h"
#include "TCanvas.h"
#include "TGraph.h"
// 
#include "AliTracker.h"   
#include "AliESDEvent.h"   
#include "AliESD.h"
#include "AliESDfriend.h"
#include "AliESDfriendTrack.h"
#include "AliRecInfoCuts.h" 
#include "AliMCInfoCuts.h" 
#include "AliLog.h" 
#include "AliESDVertex.h" 
//
#include "AliMathBase.h"
#include "AliTreeDraw.h" 

#include "AliMCInfo.h" 
#include "AliESDRecInfo.h" 
#include "AliComparisonDCA.h" 

using namespace std;

ClassImp(AliComparisonDCA)

//_____________________________________________________________________________
AliComparisonDCA::AliComparisonDCA():
//  TNamed("AliComparisonDCA","AliComparisonDCA"),
  AliComparisonObject("AliComparisonDCA"),

  // DCA histograms
  fD0TanSPtB1(0),
  fD1TanSPtB1(0),
  fD0TanSPtL1(0),
  fD1TanSPtL1(0),
  fD0TanSPtInTPC(0),
  fD1TanSPtInTPC(0),
  fVertex(0),

  // Cuts 
  fCutsRC(0), 
  fCutsMC(0),  

  // histogram folder 
  fAnalysisFolder(0)
{
  Init();
}

//_____________________________________________________________________________
AliComparisonDCA::~AliComparisonDCA()
{
  //
  if(fD0TanSPtB1) delete fD0TanSPtB1; fD0TanSPtB1=0;
  if(fD1TanSPtB1) delete fD1TanSPtB1; fD1TanSPtB1=0;
  if(fD0TanSPtL1) delete fD0TanSPtL1; fD0TanSPtL1=0;
  if(fD1TanSPtL1) delete fD1TanSPtL1; fD1TanSPtL1=0;
  if(fD0TanSPtInTPC) delete fD0TanSPtInTPC; fD0TanSPtInTPC=0;
  if(fD1TanSPtInTPC) delete fD1TanSPtInTPC; fD1TanSPtInTPC=0;
  if(fVertex) delete fVertex; fVertex=0;
  if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;

}

//_____________________________________________________________________________
void AliComparisonDCA::Init()
{
  // DCA histograms
  fD0TanSPtB1 = new TH3F("DCAyTanSPtB1","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
  fD0TanSPtB1->SetXTitle("tan(#theta)");
  fD0TanSPtB1->SetYTitle("#sqrt{p_{t}(GeV/c)}");
  fD0TanSPtB1->SetZTitle("DCA_{xy}");

  fD1TanSPtB1 = new TH3F("DCAzTanSPtB1","DCAzTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
  fD1TanSPtB1->SetXTitle("tan(#theta)");
  fD1TanSPtB1->SetYTitle("#sqrt(p_{t}(GeV/c))");
  fD1TanSPtB1->SetZTitle("DCA_{z}");

  fD0TanSPtL1 = new TH3F("DCAyTanSPtL1","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
  fD0TanSPtL1->SetXTitle("tan(#theta)");
  fD0TanSPtL1->SetYTitle("#sqrt{p_{t}(GeV/c)}");
  fD0TanSPtL1->SetZTitle("DCA_{xy}");

  fD1TanSPtL1 = new TH3F("DCAzTanSPtL1","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
  fD1TanSPtL1->SetXTitle("tan(#theta)");
  fD1TanSPtL1->SetYTitle("#sqrt{p_{t}(GeV/c)}");
  fD1TanSPtL1->SetZTitle("DCA_{z}");

  fD0TanSPtInTPC = new TH3F("DCAyTanSPtInTPC","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
  fD0TanSPtInTPC->SetXTitle("tan(#theta)");
  fD0TanSPtInTPC->SetYTitle("#sqrt{p_{t}(GeV/c)}");
  fD0TanSPtInTPC->SetZTitle("DCA_{xy}");

  fD1TanSPtInTPC = new TH3F("DCAzTanSPtInTPC","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
  fD1TanSPtInTPC->SetXTitle("tan(#theta)");
  fD1TanSPtInTPC->SetYTitle("#sqrt{p_{t}(GeV/c)}");
  fD1TanSPtInTPC->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");

  // vertex (0,0,0)
  fVertex = new AliESDVertex();
  fVertex->SetXv(0.0);
  fVertex->SetYv(0.0);
  fVertex->SetZv(0.0);
}

//_____________________________________________________________________________
void AliComparisonDCA::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC)
{
  // Fill DCA comparison information
  AliExternalTrackParam *track = 0;
  Double_t field      = AliTracker::GetBz(); // nominal Bz field [kG]
  Double_t kMaxD      = 123456.0; // max distance

  Int_t clusterITS[200];
  Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
  Float_t dca1[2], cov1[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z

  Float_t mcpt = infoMC->GetParticle().Pt();
  Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5);
  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
  fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );
  fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );
  fVertex->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(fVertex,field,kMaxD,dca,cov);

      if(bDCAStatus) {
        fD0TanSPtInTPC->Fill(tantheta,spt,dca[0]);
        fD1TanSPtInTPC->Fill(tantheta,spt,dca[1]);
	  }

	  delete track;
    }
  }
  
 // ITS + TPC
 infoRC->GetESDtrack()->GetImpactParameters(dca1,cov1);

 if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)==0){
    fD0TanSPtB1->Fill(tantheta,spt,dca1[0]);
    fD1TanSPtB1->Fill(tantheta,spt,dca1[1]);
  }
    fD0TanSPtL1->Fill(tantheta,spt,dca1[0]);
    fD1TanSPtL1->Fill(tantheta,spt,dca1[1]);
}

//_____________________________________________________________________________
Long64_t AliComparisonDCA::Merge(TCollection* 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; 
    

    fD0TanSPtB1->Add(entry->fD0TanSPtB1);
    fD1TanSPtB1->Add(entry->fD1TanSPtB1);
    fD0TanSPtL1->Add(entry->fD0TanSPtL1);
    fD1TanSPtL1->Add(entry->fD1TanSPtL1);
    fD0TanSPtInTPC->Add(entry->fD0TanSPtInTPC);
    fD1TanSPtInTPC->Add(entry->fD1TanSPtInTPC);

    count++;
  }

return count;
}

//_____________________________________________________________________________
void AliComparisonDCA::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC){
  // Process comparison information
  Process(infoMC,infoRC);
}

//_____________________________________________________________________________
void AliComparisonDCA::Analyse()
{
  //
  // Analyse comparison information and store output histograms
  // in the analysis folder "folderDCA" 
  //
  
  TH1::AddDirectory(kFALSE);

  TGraph * gr[4]= { 0,0,0,0 };
  TGraph2D *gr2[4]= { 0,0,0,0};
  AliComparisonDCA * comp=this;
  TObjArray *aFolderObj = new TObjArray;

  // 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->fD0TanSPtB1,2,5);
  gr[0]->GetXaxis()->SetTitle("Tan(#theta)");
  gr[0]->GetYaxis()->SetTitle("#sigmaDCA_xy (cm)");
  gr[0]->SetName("DCAXYResolTan");
  gr[0]->Draw("Al*");

  aFolderObj->Add(gr[0]);

  c->cd(2);
  gr[1] = AliMathBase::MakeStat1D(comp->fD1TanSPtB1,2,5);
  gr[1]->GetXaxis()->SetTitle("Tan(#theta)");
  gr[1]->GetYaxis()->SetTitle("#sigmaDCA_z (cm)");
  gr[1]->SetName("DCAZResolTan");
  gr[1]->Draw("Al*");

  aFolderObj->Add(gr[1]);

  //
  // DCA mean value
  //
  c->cd(3);
  gr[2] = AliMathBase::MakeStat1D(comp->fD0TanSPtB1,2,4);
  gr[2]->GetXaxis()->SetTitle("Tan(#theta)");
  gr[2]->GetYaxis()->SetTitle("mean DCA_xy (cm)");
  gr[2]->SetName("DCAXYMeanTan");
  gr[2]->Draw("Al*");

  aFolderObj->Add(gr[2]);

  c->cd(4);
  gr[3] = AliMathBase::MakeStat1D(comp->fD1TanSPtB1,2,4);
  gr[3]->GetXaxis()->SetTitle("Tan(#theta)");
  gr[3]->GetYaxis()->SetTitle("mean DCA_z (cm)");
  gr[3]->SetName("DCAZMeanTan");
  gr[3]->Draw("Al*");

  aFolderObj->Add(gr[3]);

  // 2D DCA resolution 
  c->cd(5);
  gr2[0] = AliMathBase::MakeStat2D(comp->fD0TanSPtB1,4,2,5); 
  gr2[0]->GetXaxis()->SetTitle("Tan(#theta)");
  gr2[0]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
  gr2[0]->GetZaxis()->SetTitle("#sigmaDCA_xy (cm)");
  gr2[0]->SetName("DCAXYResolSPTTan");
  gr2[0]->GetHistogram()->Draw("colz");

  gr2[0]->GetHistogram()->SetName("DCAXYResolSPTTan");
  aFolderObj->Add(gr2[0]->GetHistogram());

  c->cd(6);
  gr2[1] = AliMathBase::MakeStat2D(comp->fD1TanSPtB1,4,2,5); 
  gr2[1]->GetXaxis()->SetTitle("Tan(#theta)");
  gr2[1]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
  gr2[1]->GetZaxis()->SetTitle("#sigmaDCA_z (cm)");
  gr2[1]->SetName("DCAZResolSPTTan");
  gr2[1]->GetHistogram()->Draw("colz");

  gr2[1]->GetHistogram()->SetName("DCAZResolSPTTan");
  aFolderObj->Add(gr2[1]->GetHistogram());

  // 2D DCA mean value  
  c->cd(7);
  gr2[2] = AliMathBase::MakeStat2D(comp->fD0TanSPtB1,4,2,4); 
  gr2[2]->GetXaxis()->SetTitle("Tan(#theta)");
  gr2[2]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
  gr2[2]->GetZaxis()->SetTitle("mean DCA_xy (cm)");
  gr2[2]->SetName("DCAXYMeanSPTTan");
  gr2[2]->GetHistogram()->Draw("colz");

  gr2[2]->GetHistogram()->SetName("DCAXYMeanSPTTan");
  aFolderObj->Add(gr2[2]->GetHistogram());

  c->cd(8);
  gr2[3] = AliMathBase::MakeStat2D(comp->fD1TanSPtB1,4,2,4); 
  gr2[3]->GetXaxis()->SetTitle("Tan(#theta)");
  gr2[3]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
  gr2[3]->GetZaxis()->SetTitle("mean DCA_z (cm)");
  gr2[3]->SetName("DCAZMeanSPTTan");
  gr2[3]->GetHistogram()->Draw("colz");

  gr2[3]->GetHistogram()->SetName("DCAZMeanSPTTan");
  aFolderObj->Add(gr2[3]->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;
}
Commit	Line	Data
3baa4bfd	1	//------------------------------------------------------------------------------
	2	// Implementation of AliComparisonDCA class. It keeps information from
	3	// comparison of reconstructed and MC particle tracks. In addtion,
	4	// it keeps selection cuts used during comparison. The comparison
	5	// information is stored in the ROOT histograms. Analysis of these
	6	// histograms can be done by using Analyse() class function. The result of
	7	// the analysis (histograms/graphs) are stored in the folder
	8	// which is a data member of AliComparisonDCA.
	9	//
	10	// Author: J.Otwinowski 04/02/2008
	11	//------------------------------------------------------------------------------
	12
	13	/*
	14
	15	// after running comparison task, read the file, and get component
b4126c69	16	gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
b4126c69	17	LoadMyLibs();
3baa4bfd	18	TFile f("Output.root");
	19	AliComparisonDCA * compObj = (AliComparisonDCA*)f.Get("AliComparisonDCA");
	20
b4126c69	21	// Analyse comparison data
3baa4bfd	22	compObj->Analyse();
	23
	24	// the output histograms/graphs will be stored in the folder "folderDCA"
	25	compObj->GetAnalysisFolder()->ls("*");
	26
	27	// user can save whole comparison object (or only folder with anlysed histograms)
	28	// in the seperate output file (e.g.)
	29	TFile fout("Analysed_DCA.root","recreate");
	30	compObj->Write(); // compObj->GetAnalysisFolder()->Write();
	31	fout.Close();
	32
	33	*/
	34
	35	#include <iostream>
	36
	37	#include "TFile.h"
	38	#include "TCint.h"
	39	#include "TH3F.h"
	40	#include "TH2F.h"
	41	#include "TF1.h"
	42	#include "TProfile.h"
	43	#include "TProfile2D.h"
	44	#include "TGraph2D.h"
	45	#include "TCanvas.h"
	46	#include "TGraph.h"
	47	//
	48	#include "AliTracker.h"
	49	#include "AliESDEvent.h"
	50	#include "AliESD.h"
	51	#include "AliESDfriend.h"
	52	#include "AliESDfriendTrack.h"
	53	#include "AliRecInfoCuts.h"
	54	#include "AliMCInfoCuts.h"
	55	#include "AliLog.h"
	56	#include "AliESDVertex.h"
	57	//
	58	#include "AliMathBase.h"
	59	#include "AliTreeDraw.h"
	60
	61	#include "AliMCInfo.h"
	62	#include "AliESDRecInfo.h"
	63	#include "AliComparisonDCA.h"
	64
	65	using namespace std;
	66
	67	ClassImp(AliComparisonDCA)
	68
	69	//_____________________________________________________________________________
	70	AliComparisonDCA::AliComparisonDCA():
	71	// TNamed("AliComparisonDCA","AliComparisonDCA"),
	72	AliComparisonObject("AliComparisonDCA"),
	73
	74	// DCA histograms
	75	fD0TanSPtB1(0),
	76	fD1TanSPtB1(0),
	77	fD0TanSPtL1(0),
	78	fD1TanSPtL1(0),
	79	fD0TanSPtInTPC(0),
	80	fD1TanSPtInTPC(0),
	81	fVertex(0),
	82
	83	// Cuts
	84	fCutsRC(0),
	85	fCutsMC(0),
86
87	// histogram folder
88	fAnalysisFolder(0)
89	{
90	Init();
3baa4bfd	91	}
	92
	93	//_____________________________________________________________________________
	94	AliComparisonDCA::~AliComparisonDCA()
	95	{
	96	//
	97	if(fD0TanSPtB1) delete fD0TanSPtB1; fD0TanSPtB1=0;
	98	if(fD1TanSPtB1) delete fD1TanSPtB1; fD1TanSPtB1=0;
	99	if(fD0TanSPtL1) delete fD0TanSPtL1; fD0TanSPtL1=0;
	100	if(fD1TanSPtL1) delete fD1TanSPtL1; fD1TanSPtL1=0;
	101	if(fD0TanSPtInTPC) delete fD0TanSPtInTPC; fD0TanSPtInTPC=0;
	102	if(fD1TanSPtInTPC) delete fD1TanSPtInTPC; fD1TanSPtInTPC=0;
	103	if(fVertex) delete fVertex; fVertex=0;
	104	if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
	105
	106	}
	107
	108	//_____________________________________________________________________________
	109	void AliComparisonDCA::Init()
	110	{
	111	// DCA histograms
	112	fD0TanSPtB1 = new TH3F("DCAyTanSPtB1","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
	113	fD0TanSPtB1->SetXTitle("tan(#theta)");
	114	fD0TanSPtB1->SetYTitle("#sqrt{p_{t}(GeV/c)}");
	115	fD0TanSPtB1->SetZTitle("DCA_{xy}");
	116
	117	fD1TanSPtB1 = new TH3F("DCAzTanSPtB1","DCAzTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
	118	fD1TanSPtB1->SetXTitle("tan(#theta)");
	119	fD1TanSPtB1->SetYTitle("#sqrt(p_{t}(GeV/c))");
	120	fD1TanSPtB1->SetZTitle("DCA_{z}");
	121
	122	fD0TanSPtL1 = new TH3F("DCAyTanSPtL1","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
	123	fD0TanSPtL1->SetXTitle("tan(#theta)");
	124	fD0TanSPtL1->SetYTitle("#sqrt{p_{t}(GeV/c)}");
	125	fD0TanSPtL1->SetZTitle("DCA_{xy}");
	126
	127	fD1TanSPtL1 = new TH3F("DCAzTanSPtL1","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
	128	fD1TanSPtL1->SetXTitle("tan(#theta)");
	129	fD1TanSPtL1->SetYTitle("#sqrt{p_{t}(GeV/c)}");
	130	fD1TanSPtL1->SetZTitle("DCA_{z}");
	131
	132	fD0TanSPtInTPC = new TH3F("DCAyTanSPtInTPC","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
	133	fD0TanSPtInTPC->SetXTitle("tan(#theta)");
	134	fD0TanSPtInTPC->SetYTitle("#sqrt{p_{t}(GeV/c)}");
	135	fD0TanSPtInTPC->SetZTitle("DCA_{xy}");
	136
	137	fD1TanSPtInTPC = new TH3F("DCAzTanSPtInTPC","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
	138	fD1TanSPtInTPC->SetXTitle("tan(#theta)");
	139	fD1TanSPtInTPC->SetYTitle("#sqrt{p_{t}(GeV/c)}");
	140	fD1TanSPtInTPC->SetZTitle("DCA_{z}");
	141
	142	// init cuts
	143	if(!fCutsMC)
	144	AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
	145	if(!fCutsRC)
	146	AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
	147
	148	// init folder
	149	fAnalysisFolder = CreateFolder("folderDCA","Analysis DCA Folder");
b4126c69	150
	151	// vertex (0,0,0)
	152	fVertex = new AliESDVertex();
	153	fVertex->SetXv(0.0);
	154	fVertex->SetYv(0.0);
	155	fVertex->SetZv(0.0);
3baa4bfd	156	}
	157
	158	//_____________________________________________________________________________
	159	void AliComparisonDCA::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC)
	160	{
	161	// Fill DCA comparison information
	162	AliExternalTrackParam *track = 0;
	163	Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]
	164	Double_t kMaxD = 123456.0; // max distance
	165
	166	Int_t clusterITS[200];
	167	Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
decf0997	168	Float_t dca1[2], cov1[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
3baa4bfd	169
	170	Float_t mcpt = infoMC->GetParticle().Pt();
	171	Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5);
	172	Float_t spt = TMath::Sqrt(mcpt);
	173
	174	// distance to Prim. vertex
	175	const Double_t* dv = infoMC->GetVDist();
	176
	177	Bool_t isPrim = TMath::Sqrt(dv[0]dv[0] + dv[1]dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
	178
	179	// Check selection cuts
	180	if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
	181	if (!isPrim) return;
	182	if (infoRC->GetStatus(1)!=3) return;
	183	if (!infoRC->GetESDtrack()) return;
	184	if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
	185	if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;
	186
	187	// calculate and set prim. vertex
	188	fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );
	189	fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );
	190	fVertex->SetZv( infoMC->GetParticle().Vz() - dv[2] );
	191
	192	// calculate track parameters at vertex
	193	if (infoRC->GetESDtrack()->GetTPCInnerParam())
	194	{
	195	if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )
	196	{
	197	Bool_t bDCAStatus = track->PropagateToDCA(fVertex,field,kMaxD,dca,cov);
	198
	199	if(bDCAStatus) {
	200	fD0TanSPtInTPC->Fill(tantheta,spt,dca[0]);
	201	fD1TanSPtInTPC->Fill(tantheta,spt,dca[1]);
	202	}
	203
	204	delete track;
	205	}
	206	}
decf0997	207
	208	// ITS + TPC
	209	infoRC->GetESDtrack()->GetImpactParameters(dca1,cov1);
3baa4bfd	210
3baa4bfd	211	if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)==0){
decf0997	212	fD0TanSPtB1->Fill(tantheta,spt,dca1[0]);
decf0997	213	fD1TanSPtB1->Fill(tantheta,spt,dca1[1]);
3baa4bfd	214	}
decf0997	215	fD0TanSPtL1->Fill(tantheta,spt,dca1[0]);
decf0997	216	fD1TanSPtL1->Fill(tantheta,spt,dca1[1]);
3baa4bfd	217	}
	218
	219	//_____________________________________________________________________________
	220	Long64_t AliComparisonDCA::Merge(TCollection* list)
	221	{
	222	// Merge list of objects (needed by PROOF)
	223
	224	if (!list)
	225	return 0;
	226
	227	if (list->IsEmpty())
	228	return 1;
	229
	230	TIterator* iter = list->MakeIterator();
	231	TObject* obj = 0;
	232
	233	// collection of generated histograms
	234	Int_t count=0;
	235	while((obj = iter->Next()) != 0)
	236	{
	237	AliComparisonDCA* entry = dynamic_cast<AliComparisonDCA*>(obj);
	238	if (entry == 0) continue;
	239
	240
	241	fD0TanSPtB1->Add(entry->fD0TanSPtB1);
	242	fD1TanSPtB1->Add(entry->fD1TanSPtB1);
	243	fD0TanSPtL1->Add(entry->fD0TanSPtL1);
	244	fD1TanSPtL1->Add(entry->fD1TanSPtL1);
	245	fD0TanSPtInTPC->Add(entry->fD0TanSPtInTPC);
	246	fD1TanSPtInTPC->Add(entry->fD1TanSPtInTPC);
	247
	248	count++;
	249	}
	250
	251	return count;
	252	}
	253
	254	//_____________________________________________________________________________
	255	void AliComparisonDCA::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC){
	256	// Process comparison information
	257	Process(infoMC,infoRC);
	258	}
	259
	260	//_____________________________________________________________________________
	261	void AliComparisonDCA::Analyse()
	262	{
	263	//
	264	// Analyse comparison information and store output histograms
	265	// in the analysis folder "folderDCA"
	266	//
	267
	268	TH1::AddDirectory(kFALSE);
	269
decf0997	270	TGraph * gr[4]= { 0,0,0,0 };
decf0997	271	TGraph2D *gr2[4]= { 0,0,0,0};
3baa4bfd	272	AliComparisonDCA * comp=this;
b4126c69	273	TObjArray *aFolderObj = new TObjArray;
3baa4bfd	274
	275	// write results in the folder
	276	// Canvas to draw analysed histograms
	277	TCanvas * c = new TCanvas("canDCA","DCA resolution");
decf0997	278	c->Divide(2,4);
3baa4bfd	279	//
	280	// DCA resolution
	281	//
decf0997	282	c->cd(1);
	283	gr[0] = AliMathBase::MakeStat1D(comp->fD0TanSPtB1,2,5);
	284	gr[0]->GetXaxis()->SetTitle("Tan(#theta)");
	285	gr[0]->GetYaxis()->SetTitle("#sigmaDCA_xy (cm)");
	286	gr[0]->SetName("DCAXYResolTan");
	287	gr[0]->Draw("Al*");
	288
	289	aFolderObj->Add(gr[0]);
	290
3baa4bfd	291	c->cd(2);
decf0997	292	gr[1] = AliMathBase::MakeStat1D(comp->fD1TanSPtB1,2,5);
	293	gr[1]->GetXaxis()->SetTitle("Tan(#theta)");
	294	gr[1]->GetYaxis()->SetTitle("#sigmaDCA_z (cm)");
	295	gr[1]->SetName("DCAZResolTan");
	296	gr[1]->Draw("Al*");
	297
	298	aFolderObj->Add(gr[1]);
	299
	300	//
	301	// DCA mean value
	302	//
	303	c->cd(3);
	304	gr[2] = AliMathBase::MakeStat1D(comp->fD0TanSPtB1,2,4);
	305	gr[2]->GetXaxis()->SetTitle("Tan(#theta)");
	306	gr[2]->GetYaxis()->SetTitle("mean DCA_xy (cm)");
	307	gr[2]->SetName("DCAXYMeanTan");
	308	gr[2]->Draw("Al*");
	309
	310	aFolderObj->Add(gr[2]);
	311
	312	c->cd(4);
	313	gr[3] = AliMathBase::MakeStat1D(comp->fD1TanSPtB1,2,4);
	314	gr[3]->GetXaxis()->SetTitle("Tan(#theta)");
	315	gr[3]->GetYaxis()->SetTitle("mean DCA_z (cm)");
	316	gr[3]->SetName("DCAZMeanTan");
	317	gr[3]->Draw("Al*");
	318
	319	aFolderObj->Add(gr[3]);
	320
	321	// 2D DCA resolution
	322	c->cd(5);
	323	gr2[0] = AliMathBase::MakeStat2D(comp->fD0TanSPtB1,4,2,5);
	324	gr2[0]->GetXaxis()->SetTitle("Tan(#theta)");
	325	gr2[0]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
	326	gr2[0]->GetZaxis()->SetTitle("#sigmaDCA_xy (cm)");
	327	gr2[0]->SetName("DCAXYResolSPTTan");
	328	gr2[0]->GetHistogram()->Draw("colz");
	329
	330	gr2[0]->GetHistogram()->SetName("DCAXYResolSPTTan");
	331	aFolderObj->Add(gr2[0]->GetHistogram());
	332
	333	c->cd(6);
	334	gr2[1] = AliMathBase::MakeStat2D(comp->fD1TanSPtB1,4,2,5);
	335	gr2[1]->GetXaxis()->SetTitle("Tan(#theta)");
	336	gr2[1]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
	337	gr2[1]->GetZaxis()->SetTitle("#sigmaDCA_z (cm)");
	338	gr2[1]->SetName("DCAZResolSPTTan");
	339	gr2[1]->GetHistogram()->Draw("colz");
	340
	341	gr2[1]->GetHistogram()->SetName("DCAZResolSPTTan");
	342	aFolderObj->Add(gr2[1]->GetHistogram());
	343
	344	// 2D DCA mean value
	345	c->cd(7);
	346	gr2[2] = AliMathBase::MakeStat2D(comp->fD0TanSPtB1,4,2,4);
	347	gr2[2]->GetXaxis()->SetTitle("Tan(#theta)");
	348	gr2[2]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
	349	gr2[2]->GetZaxis()->SetTitle("mean DCA_xy (cm)");
	350	gr2[2]->SetName("DCAXYMeanSPTTan");
	351	gr2[2]->GetHistogram()->Draw("colz");
	352
	353	gr2[2]->GetHistogram()->SetName("DCAXYMeanSPTTan");
	354	aFolderObj->Add(gr2[2]->GetHistogram());
	355
356	c->cd(8);
357	gr2[3] = AliMathBase::MakeStat2D(comp->fD1TanSPtB1,4,2,4);
358	gr2[3]->GetXaxis()->SetTitle("Tan(#theta)");
359	gr2[3]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV/c)}");
360	gr2[3]->GetZaxis()->SetTitle("mean DCA_z (cm)");
361	gr2[3]->SetName("DCAZMeanSPTTan");
362	gr2[3]->GetHistogram()->Draw("colz");
363
364	gr2[3]->GetHistogram()->SetName("DCAZMeanSPTTan");
365	aFolderObj->Add(gr2[3]->GetHistogram());
b4126c69	366
	367	// export objects to analysis folder
	368	fAnalysisFolder = ExportToFolder(aFolderObj);
3baa4bfd	369
b4126c69	370	// delete only TObjArray
b4126c69	371	if(aFolderObj) delete aFolderObj;
3baa4bfd	372	}
3baa4bfd	373
b4126c69	374	//_____________________________________________________________________________
	375	TFolder* AliComparisonDCA::ExportToFolder(TObjArray * array)
	376	{
	377	// recreate folder avery time and export objects to new one
	378	//
	379	AliComparisonDCA * comp=this;
	380	TFolder *folder = comp->GetAnalysisFolder();
	381
	382	TString name, title;
	383	TFolder *newFolder = 0;
	384	Int_t i = 0;
	385	Int_t size = array->GetSize();
	386
	387	if(folder) {
	388	// get name and title from old folder
	389	name = folder->GetName();
	390	title = folder->GetTitle();
	391
	392	// delete old one
	393	delete folder;
	394
	395	// create new one
	396	newFolder = CreateFolder(name.Data(),title.Data());
	397	newFolder->SetOwner();
	398
	399	// add objects to folder
	400	while(i < size) {
	401	newFolder->Add(array->At(i));
	402	i++;
	403	}
	404	}
	405
	406	return newFolder;
	407	}
	408
	409
3baa4bfd	410	//_____________________________________________________________________________
	411	TFolder* AliComparisonDCA::CreateFolder(TString name,TString title) {
	412	// create folder for analysed histograms
	413	TFolder *folder = 0;
	414	folder = new TFolder(name.Data(),title.Data());
	415
	416	return folder;
	417	}