--- /dev/null
+
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include <TH1F.h>
+#include <TH1I.h>
+#include <TH2D.h>
+#include <TFile.h>
+#include <Riostream.h>
+#include <TParticle.h>
+#include <TCanvas.h>
+#include <TLegend.h>
+#include <TSystem.h>
+#include <TROOT.h>
+#include <TClonesArray.h>
+#include <TTree.h>
+#include <TArrayF.h>
+#include <TStyle.h>
+#include "AliRun.h"
+#include "AliRunLoader.h"
+#include "AliStack.h"
+#include "AliESDEvent.h"
+#include "AliESDtrack.h"
+#include "AliTrackReference.h"
+#include "AliITSsegmentationUpgrade.h"
+#endif
+
+Bool_t IsTrackable(TClonesArray *trackRef);
+TArrayF radii;
+
+void efficiency(char *title="",Bool_t isPrimary=kTRUE){
+
+ gSystem->Load("libITSUpgradeBase.so");
+ gSystem->Load("libITSUpgradeRec.so");
+ gSystem->Load("libITSUpgradeSim.so");
+
+ // setting up the configuration to be considered
+ AliITSsegmentationUpgrade *seg = new AliITSsegmentationUpgrade();
+ if(!seg){
+ printf("no segmentation info available... Exiting");
+ return;
+ }
+ radii.Set(seg->GetNLayers());
+ for(Int_t l=0; l<seg->GetNLayers(); l++){
+ radii.AddAt(seg->GetRadius(l),l);
+ }
+ delete seg;
+
+ TString titlet(title);
+ titlet.ReplaceAll(" ","");
+
+ gROOT->SetStyle("Plain");
+ gStyle->SetPalette(1);
+ gStyle->SetOptStat(0);
+ gStyle->SetOptTitle(0);
+
+
+ //booking plots
+ // pt
+ Int_t nbins= 64;
+ Double_t ptmax = 1.6;
+
+ TH1F *hPtRec = new TH1F("hPtRec"," p_{T} distribution of reconstructed tracks",nbins,0,ptmax);
+ hPtRec->SetXTitle("p_{T}");
+ hPtRec->SetYTitle(Form("entries / %3.0f MeV/c",1000*(ptmax/nbins)));
+
+ TH1F *hPtRecMc = new TH1F("hPtRecMc"," MC p_{T} distribution of reconstructed tracks",nbins,0,ptmax);
+ hPtRecMc->SetXTitle("p_{T} GeV/c");
+ hPtRecMc->SetYTitle(Form("entries / %3.0f MeV/c",1000*(ptmax/nbins)));
+
+ TH1F *hPtRecGood = new TH1F("hPtRecGood"," P_{t} distribution of reconstructed tracks [positive labels]",nbins,0,ptmax);
+ hPtRecGood->SetXTitle("p_{T} GeV/c");
+ hPtRecGood->SetYTitle(Form("entries / %3.0f MeV/c",1000*(ptmax/nbins)));
+
+ TH1F *hPtMc = new TH1F("hPtMc", " p_{T} distribution of MC tracks",nbins,0,ptmax);
+ hPtMc->SetXTitle("p_{T} GeV/c");
+ hPtMc->SetYTitle(Form("entries / %3.0f MeV/c",1000*(ptmax/nbins)));
+
+ // Eta
+ Int_t nbinsEta = 24;
+ Double_t eta[2]={-1.2,1.2};
+
+ TH1F *hEtaRec = new TH1F("hEtaRec"," #eta distribution of reconstructed tracks",nbinsEta,eta[0],eta[1]);
+ hEtaRec->SetXTitle("#eta");
+ hEtaRec->SetYTitle(Form("entries / %3.3f ",((eta[1]-eta[0])/nbinsEta)));
+
+ TH1F *hEtaRecMc = new TH1F("hEtaRecMc"," #eta Mc distribution of reconstructed tracks",nbinsEta,eta[0],eta[1]);
+ hEtaRecMc->SetXTitle("#eta");
+ hEtaRecMc->SetYTitle(Form("entries / %3.3f ",((eta[1]-eta[0])/nbinsEta)));
+
+ TH1F *hEtaRecGood = new TH1F("hEtaRecGood"," #eta distribution of Good reconstructed tracks",nbinsEta,eta[0],eta[1]);
+ hEtaRecGood->SetXTitle("#eta");
+ hEtaRecGood->SetYTitle(Form("entries / %3.3f ",((eta[1]-eta[0])/nbinsEta)));
+
+ TH1F *hEtaMc = new TH1F("hEtaMc", " #eta distribution of MC tracks",nbinsEta,eta[0],eta[1]);
+ hEtaMc->SetXTitle("#eta");
+ hEtaMc->SetYTitle(Form("entries / %3.3f ",((eta[1]-eta[0])/nbinsEta)));
+
+
+ Double_t ptlimit=0.0;
+
+ Double_t cont=0;
+ Double_t ntrack=0;
+
+ // Init simulation
+ gAlice=NULL;
+ AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
+ if (!runLoader) {
+ Error("OpengAlice", "no files, please check the path");
+ return ;
+ }
+ runLoader->LoadgAlice();
+ gAlice = runLoader->GetAliRun();
+ runLoader->LoadHeader();
+ runLoader->LoadKinematics();
+ runLoader->LoadTrackRefs();
+
+ //Trackref
+ TTree *trackRefTree = 0x0;
+ TClonesArray *trackRef = new TClonesArray("AliTrackReference",1000);
+
+ // ESD
+ TFile* esdFile = TFile::Open("AliESDs.root");
+ if (!esdFile || !esdFile->IsOpen()) {
+ Error("CheckESD", "opening ESD file %s failed", "AliESDs.root");
+ return ;
+ }
+ AliESDEvent * esd = new AliESDEvent;
+ TTree* tree = (TTree*) esdFile->Get("esdTree");
+ if (!tree) {
+ Error("CheckESD", "no ESD tree found");
+ return ;
+ }
+ esd->ReadFromTree(tree);
+
+ printf("setting the low pt value to %f \n",ptlimit);
+ // Loop over events
+
+ for(Int_t i=0; i<runLoader->GetNumberOfEvents(); i++){
+ //for(Int_t i=0; i< 1; i++){
+ cout << "Event "<< i << endl;
+ runLoader->GetEvent(i);
+
+ AliStack *stack = runLoader->Stack();
+ trackRefTree=runLoader->TreeTR();
+ TBranch *br = trackRefTree->GetBranch("TrackReferences");
+ if(!br) {
+ printf("no TR branch available , exiting \n");
+ return;
+ }
+ br->SetAddress(&trackRef);
+
+ Int_t nParticles = stack->GetNtrack();
+ trackRefTree=runLoader->TreeTR();
+ trackRefTree->SetBranchAddress("TrackReferences",&trackRef);
+ for(Int_t p=0; p<nParticles; p++){
+ if(stack->Particle(p)->Pt()<ptlimit) continue;
+ if(isPrimary && !stack->IsPhysicalPrimary(p)) continue;
+ trackRefTree->GetEntry(stack->TreeKEntry(p));
+
+ if(IsTrackable(trackRef)) {
+ cont++;
+ hPtMc->Fill((stack->Particle(p))->Pt());
+ hEtaMc->Fill((stack->Particle(p))->Eta());
+ }
+ }//loop on kine particle
+ cout << " # trackable "<< cont;
+
+ // ESD
+ tree->GetEvent(i);
+ for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
+ AliESDtrack* track = esd->GetTrack(iTrack);
+ if(track->Pt()<ptlimit) continue;
+ if(isPrimary && !stack->IsPhysicalPrimary(TMath::Abs(track->GetLabel()))) continue;
+ ntrack++;
+ hPtRec->Fill(track->Pt());
+ hEtaRec->Fill(track->Eta());
+ if(track->GetLabel()>=0) {
+ hPtRecGood->Fill(stack->Particle(track->GetLabel())->Pt());
+ hEtaRecGood->Fill(stack->Particle(track->GetLabel())->Eta());
+ }
+ hPtRecMc->Fill(stack->Particle(TMath::Abs(track->GetLabel()))->Pt());
+ hEtaRecMc->Fill(stack->Particle(TMath::Abs(track->GetLabel()))->Eta());
+ }// loop on tracks
+ cout << " tracks " << ntrack << endl;
+
+ }//loop events
+
+ cout << "Global Efficiency "<< ntrack/cont << endl;
+
+ TCanvas *c = new TCanvas();
+ c->cd();
+ c->cd()->SetLogy();
+ c->cd()->SetGridy();
+ c->cd()->SetGridx();
+ c->cd()->SetObjectStat(0);
+ hPtMc->SetMarkerStyle(20);
+ hPtMc->Draw("err");
+ hPtRec->SetLineColor(kRed);
+ hPtRec->SetMarkerStyle(20);
+ hPtRec->SetMarkerColor(kRed);
+ hPtRec->Sumw2();
+ hPtRec->Draw("sameerr");
+ hPtRecMc->SetLineColor(kBlue);
+ hPtRecMc->SetMarkerStyle(20);
+ hPtRecMc->SetMarkerColor(kBlue);
+ hPtRecMc->Sumw2();
+ hPtRecMc->Draw("sameerr");
+ hPtRecGood->SetLineColor(kGreen);
+ hPtRecGood->SetMarkerStyle(20);
+ hPtRecGood->SetMarkerColor(kGreen);
+ hPtRecGood->Sumw2();
+ hPtRecGood->Draw("sameerr");
+
+ TLegend *leg = new TLegend(0.5,0.7,0.99,0.95);
+ leg->SetHeader(Form("%s",title));
+ leg->AddEntry(hPtRecGood,"Good tracks (positive labels)","p");
+ leg->AddEntry(hPtRec,"Reconstructed tracks (Rec p_{T})","p");
+ leg->AddEntry(hPtRecMc,"Reconstructed tracks ","p");
+ leg->AddEntry(hPtMc,"MC trackable particles","p");
+ leg->Draw();
+ c->SaveAs(Form("ptdist_%s.png",titlet.Data()));
+
+ TCanvas *cc = new TCanvas();
+ cc->cd();
+ cc->cd()->SetGridy();
+ cc->cd()->SetGridx();
+
+ Double_t max = 1.2;
+
+ TH1F *effrec=new TH1F("hEffRec","tracking efficiency of reconstructed tracks",nbins,0,ptmax);
+ effrec->Divide(hPtRecMc,hPtMc,1,1,"b");
+ effrec->SetXTitle("p_{T} GeV/c");
+ effrec->SetLineColor(kBlue);
+ effrec->SetMarkerStyle(20);
+ effrec->SetMarkerColor(kBlue);
+ effrec->SetMaximum(max);
+ effrec->SetMinimum(0);
+ effrec->DrawCopy("err");
+
+
+ TH1F *effrecgood=new TH1F("hEffRecGood","tracking efficiency of reconstructed tracks [positive labels]",nbins,0,ptmax);
+ effrecgood->SetXTitle("P_{T}");
+ effrecgood->Divide(hPtRecGood,hPtMc,1,1,"b");
+ effrecgood->SetLineColor(kGreen);
+ effrecgood->SetMarkerStyle(20);
+ effrecgood->SetMarkerColor(kGreen);
+ effrecgood->SetMaximum(max);
+ effrecgood->SetMinimum(0);
+ effrecgood->DrawCopy("sameerr");
+
+ TH1F *purity=new TH1F("purity","track purity [positive labels / all tracks]",nbins,0,ptmax);
+ purity->Divide(hPtRecGood,hPtRecMc,1,1,"b");
+ purity->SetLineColor(kGray);
+ purity->SetMarkerStyle(28);
+ purity->SetMarkerColor(kOrange);
+ purity->SetMaximum(max);
+ purity->SetMinimum(0);
+ purity->DrawCopy("sameerr");
+
+ TH1F *fakes=new TH1F("fakes","fake ratio [negative labels / reconstructable]",nbins,0,ptmax);
+
+ fakes->Add(hPtRecMc,hPtRecGood,1,-1);
+ fakes->Divide(hPtMc);
+ fakes->SetLineColor(kViolet);
+ fakes->SetLineStyle(3);
+ fakes->SetLineWidth(3);
+ fakes->SetMarkerColor(kViolet);
+ fakes->SetMaximum(max);
+ fakes->SetMinimum(0);
+ fakes->DrawCopy("samehist");
+
+ TLegend *legEff = new TLegend(0.1,0.83,0.65,0.99);
+ legEff->SetHeader(Form("%s",title));
+ legEff->AddEntry(effrec,"Overall Efficiency ","p");
+ legEff->AddEntry(effrecgood,"Efficiency good tracks (positive labels)","p");
+ legEff->AddEntry(purity,"Purity (positive labels / All tracks)","p");
+ legEff->AddEntry(fakes,"Fakes","l");
+ legEff->Draw();
+
+ cc->SaveAs(Form("efficiency_%s.png",titlet.Data()));
+
+
+ TCanvas *cEta = new TCanvas();
+ cEta->cd();
+ cEta->cd()->SetGridy();
+ cEta->cd()->SetGridx();
+ cEta->cd()->SetObjectStat(0);
+ hEtaMc->SetMarkerStyle(20);
+ hEtaMc->Draw("err");
+ hEtaMc->SetMinimum(0);
+ hEtaRec->SetLineColor(kRed);
+ hEtaRec->SetMarkerStyle(20);
+ hEtaRec->SetMarkerColor(kRed);
+ hEtaRec->Sumw2();
+ hEtaRec->Draw("sameerr");
+ hEtaRecMc->SetLineColor(kBlue);
+ hEtaRecMc->SetMarkerStyle(20);
+ hEtaRecMc->SetMarkerColor(kBlue);
+ hEtaRecMc->Sumw2();
+ hEtaRecMc->Draw("sameerr");
+ hEtaRecGood->SetLineColor(kGreen);
+ hEtaRecGood->SetMarkerStyle(20);
+ hEtaRecGood->SetMarkerColor(kGreen);
+ hEtaRecGood->Sumw2();
+ hEtaRecGood->Draw("sameerr");
+
+ TLegend *legEta = new TLegend(0.2,0.3,0.59,0.55);
+ legEta->SetHeader(Form("%s",title));
+ legEta->AddEntry(hEtaRecGood,"Good tracks (positive labels)","p");
+ legEta->AddEntry(hEtaRec,"Reconstructed tracks (Rec p_{T})","p");
+ legEta->AddEntry(hEtaRecMc,"Reconstructed tracks ","p");
+ legEta->AddEntry(hEtaMc,"MC trackable particles","p");
+ legEta->Draw();
+ cEta->SaveAs(Form("etadist_%s.png",title));
+
+ TCanvas *ccEta = new TCanvas();
+ ccEta->cd();
+ ccEta->cd()->SetGridy();
+ ccEta->cd()->SetGridx();
+
+ TH1F *effrecEta=new TH1F("hEffRecEta","#eta tracking efficiency of reconstructed tracks",nbinsEta,eta[0],eta[1]);
+ effrecEta->Divide(hEtaRecMc,hEtaMc,1,1,"b");
+ effrecEta->SetXTitle("#eta");
+ effrecEta->SetLineColor(kBlue);
+ effrecEta->SetMarkerStyle(20);
+ effrecEta->SetMarkerColor(kBlue);
+ effrecEta->SetMaximum(max);
+ effrecEta->SetMinimum(0);
+ effrecEta->DrawCopy("err");
+
+
+ TH1F *effrecgoodEta=new TH1F("hEffRecGoodEta","#eta tracking efficiency of reconstructed tracks [positive labels]",nbinsEta,eta[0],eta[1]);
+ effrecgoodEta->SetXTitle("eta");
+ effrecgoodEta->Divide(hEtaRecGood,hEtaMc,1,1,"b");
+ effrecgoodEta->SetLineColor(kGreen);
+ effrecgoodEta->SetMarkerStyle(20);
+ effrecgoodEta->SetMarkerColor(kGreen);
+ effrecgoodEta->SetMaximum(max);
+ effrecgoodEta->SetMinimum(0);
+ effrecgoodEta->DrawCopy("sameerr");
+
+ TH1F *purityEta=new TH1F("purityEta","#eta track purity [positive labels / all tracks]",nbinsEta,eta[0],eta[1]);
+ purityEta->Divide(hEtaRecGood,hEtaRecMc,1,1,"b");
+ purityEta->SetLineColor(kGray);
+ purityEta->SetMarkerStyle(28);
+ purityEta->SetMarkerColor(kOrange);
+ purityEta->SetMaximum(max);
+ purityEta->SetMinimum(0);
+ purityEta->DrawCopy("sameerr");
+
+ TH1F *fakesEta=new TH1F("fakesEta","#eta fake ratio [negative labels / reconstructable]",nbinsEta,eta[0],eta[1]);
+ fakesEta->Add(hEtaRecMc,hEtaRecGood,1,-1);
+ fakesEta->Divide(hEtaMc);
+ fakesEta->SetLineColor(kViolet);
+ fakesEta->SetLineStyle(3);
+ fakesEta->SetLineWidth(3);
+ fakesEta->SetMarkerColor(kViolet);
+ fakesEta->SetMaximum(max);
+ fakesEta->SetMinimum(0);
+ fakesEta->DrawCopy("samehist");
+
+ TLegend *legEffEta = new TLegend(0.3,0.3,0.75,0.55);
+ legEffEta->SetHeader(Form("%s",title));
+ legEffEta->AddEntry(effrecEta,"Overall Efficiency ","p");
+ legEffEta->AddEntry(effrecgoodEta,"Efficiency good tracks (positive labels)","p");
+ legEffEta->AddEntry(purityEta,"Purity (positive labels / All tracks)","p");
+ legEffEta->AddEntry(fakesEta,"Fakes","l");
+ legEffEta->Draw();
+
+ ccEta->SaveAs(Form("efficiencyEta_%s.png",titlet.Data()));
+
+
+}// main
+
+//___________________________________________________
+Bool_t IsTrackable(TClonesArray *trackRef){
+
+ Bool_t isOk=kFALSE;
+
+ Int_t nTrackRef =0;
+ TArrayF isInLayer(radii.GetSize());
+ for(Int_t l=0; l<radii.GetSize(); l++ ) isInLayer.AddAt(0,l);
+ for(Int_t nT =0; nT<trackRef->GetEntries(); nT++){
+ AliTrackReference *trR = (AliTrackReference*)trackRef->At(nT);
+ if(!trR) continue;
+ if(trR->DetectorId()!=0)continue;
+ Double_t rPart = TMath::Sqrt(trR->X()*trR->X()+trR->Y()*trR->Y());
+ for(Int_t iLay=0;iLay<8;iLay++){
+ if(TMath::Abs(rPart-radii.At(iLay))<0.01) isInLayer.AddAt(1,iLay);
+ }
+ }
+
+ for(Int_t iLayer=0;iLayer<radii.GetSize();iLayer++){
+ if(isInLayer.At(iLayer)) nTrackRef++;
+ }
+
+ if(nTrackRef>2) isOk=kTRUE;
+ return isOk;
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff