////////////////////////////////////////////////////////////////////////////
#include <TROOT.h>
+#include <TFile.h>
#include <TStyle.h>
#include <TClonesArray.h>
#include <TObjArray.h>
}
if(track) fkTrack = track;
- Double_t val[11]; memset(val, 0, 11*sizeof(Double_t));
+ Double_t val[7]; memset(val, 0, 7*sizeof(Double_t));
ULong_t status(fkESD->GetStatus());
val[0] =((status&AliESDtrack::kTRDin)?1:0) +
((status&AliESDtrack::kTRDStop)?1:0) +
((status&AliESDtrack::kTRDout)?2:0);
- val[1] = fkESD->Phi();
- val[2] = fkESD->Eta();
- val[3] = DebugLevel()>=1?GetPtBin(fkESD->Pt()):GetPtBinSignificant(fkESD->Pt());
+ val[1] = fPhi;//fkESD->Phi();
+ val[2] = fEta;//fkESD->Eta();
+ val[3] = GetPtBin(fPt/*fkESD->Pt()*/);
val[4] = 0.;
if(fkMC){
if(fkMC->GetLabel() == fkMC->GetTRDlabel()) val[4] = 0.;
else if(fkMC->GetLabel() == -fkMC->GetTRDlabel()) val[4] = 1.;
else val[4] = -1.;
}
- if(fkTrack){ // read track status in debug mode with friends
- //val[4] = fkTrack->GetStatusTRD(-1);
- for(Int_t ily(0); ily<AliTRDgeometry::kNlayer; ily++) val[5+ily]=fkTrack->GetStatusTRD(ily);
- }
+ val[5] = fkTrack?fkTrack->GetNumberOfTracklets():0;
+ // down scale PID resolution
+ Int_t spc(fSpecies); if(spc==-6) spc=0; if(spc==3) spc=2; if(spc==-3) spc=-2; if(spc>3) spc=3; if(spc<-3) spc=-3;
+ val[6] = spc;
+
+ if(fkTrack) AliDebug(2, Form("%3d[%s] tracklets[%d] label[%2d %+2d] species[%d %d] in[%c] out[%c] stop[%c]",
+ fkESD->GetId(), fkTrack->IsTPCseeded()?"TPC":"ITS", Int_t(val[5]), fkMC->GetLabel(), fkMC->GetTRDlabel(), fSpecies, spc,
+ (status&AliESDtrack::kTRDin)?'x':'o', (status&AliESDtrack::kTRDout)?'x':'o', (status&AliESDtrack::kTRDStop)?'x':'o'));
H->Fill(val);
return NULL;
}
//++++++++++++++++++++++
// cluster to detector
if(!(H = (THnSparseI*)gROOT->FindObject("hEFF"))){
- const Int_t mdim(11);
- Int_t npt=DebugLevel()>=1?20:3;
+ const Int_t mdim(7);
Int_t nlabel(1);
- const Char_t *eTitle[mdim] = {"label", "#phi [rad]", "eta", "p_{t} [bin]", "label", "status[0]", "status[1]", "status[2]", "status[3]", "status[4]", "status[5]"};
- const Int_t eNbins[mdim] = {5, 180, 50, npt, nlabel, 5, 5, 5, 5, 5, 5};
- const Double_t eMin[mdim] = {-0.5, -TMath::Pi(), -1., -0.5, -0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5},
- eMax[mdim] = {4.5, TMath::Pi(), 1., npt-.5, nlabel-0.5, 5.5, 5.5, 5.5, 5.5, 5.5, 5.5};
+ const Char_t *eTitle[mdim] = {"status", "#phi [rad]", "eta", "p_{t} [bin]", "label", "N_{trklt}", "chg*spec*rc"};
+ const Int_t eNbins[mdim] = { 5, 180, 50, fNpt, nlabel, AliTRDgeometry::kNlayer-2, 5};
+ const Double_t eMin[mdim] = { -0.5, -TMath::Pi(), -1., -0.5, -0.5, 1.5, -2.5},
+ eMax[mdim] = { 4.5, TMath::Pi(), 1., fNpt-.5, nlabel-0.5, 5.5, 2.5};
st = "basic efficiency;";
// define minimum info to be saved in non debug mode
Int_t ndim=DebugLevel()>=1?mdim:(HasMCdata()?5:4);
}
if(!fProj){
AliInfo("Building array of projections ...");
- fProj = new TObjArray(50); fProj->SetOwner(kTRUE);
+ fProj = new TObjArray(200); fProj->SetOwner(kTRUE);
}
+ // set pt/p segmentation. guess from data
+ THnSparse *H(NULL);
+ if(!(H = (THnSparse*)fContainer->FindObject("hEFF"))){
+ AliError("Missing/Wrong data @ hEFF.");
+ return kFALSE;
+ }
+ fNpt=H->GetAxis(3)->GetNbins()+1;
+ if(!MakeMomSegmentation()) return kFALSE;
+
if(!MakeProjectionBasicEff()) return kFALSE;
return kTRUE;
}
Int_t ndim(H->GetNdimensions()); //Bool_t debug(ndim>Int_t(kNdimCl));
TAxis *aa[11], *al(NULL); memset(aa, 0, sizeof(TAxis*) * 11);
for(Int_t id(0); id<ndim; id++) aa[id] = H->GetAxis(id);
- if(aa[3]->GetNbins()>3) SetDebugLevel(1);
if(H->GetNdimensions() > 4) al = H->GetAxis(4);
Int_t nlab=al?3:1;
//const Int_t nEtaPhi(4); Int_t rebinEtaPhiX[nEtaPhi] = {1, 2, 5, 1}, rebinEtaPhiY[nEtaPhi] = {2, 1, 1, 5};
AliTRDrecoProjection hp[15]; TObjArray php(15);
const Char_t *stat[] = {"!TRDin", "TRDin", "TRDin&TRDStop", "TRDin&TRDout", "TRDin&TRDout&TRDStop"};
- const Char_t *lab[] = {"Bad", "Good", "Accept"};
+ const Char_t *lab[] = {"MC-Bad", "MC-Good", "MC-Accept"};
Int_t ih(0);
for(Int_t ilab(0); ilab<nlab; ilab++){
for(Int_t istat(0); istat<5; istat++){
-// isel++; // new selection
hp[ih].Build(Form("HEff%d%d", ilab, istat),
- Form("Efficiency :: Lab[%s] Stat[#bf{%s}]", lab[ilab], stat[istat]),
- 2, 1, 3, aa);
+ Form("Efficiency :: Stat[#bf{%s}] %s", stat[istat], nlab>1?lab[ilab]:""), 2, 1, 3, aa);
//hp[ih].SetRebinStrategy(nEtaPhi, rebinEtaPhiX, rebinEtaPhiY);
php.AddLast(&hp[ih++]); //np[isel]++;
}
}
AliInfo(Form("Build %3d 3D projections.", ih));
+ AliTRDrecoProjection *pr0(NULL), *pr1(NULL);
Int_t istatus, ilab(0), coord[11]; memset(coord, 0, sizeof(Int_t) * 11); Double_t v = 0.;
for (Long64_t ib(0); ib < H->GetNbins(); ib++) {
v = H->GetBinContent(ib, coord); if(v<1.) continue;
istatus = coord[0]-1;
- if(al) ilab = coord[4];
+ ilab=al?0:coord[4];
Int_t isel = ilab*5+istatus;
+ if(isel>=ih){
+ AliError(Form("Wrong selection %d [%3d]", isel, ih));
+ return kFALSE;
+ }
+ if(!(pr0=(AliTRDrecoProjection*)php.At(isel))) {
+ AliError(Form("Missing projection %d", isel));
+ return kFALSE;
+ }
+ if(strcmp(pr0->H()->GetName(), Form("HEff%d%d", ilab, istatus))!=0){
+ AliError(Form("Projection mismatch :: request[HEff%d%d] found[%s]", ilab, istatus, pr0->H()->GetName()));
+ return kFALSE;
+ }
for(Int_t jh(0); jh<1/*np[isel]*/; jh++) ((AliTRDrecoProjection*)php.At(isel+jh))->Increment(coord, v);
}
- TH2 *h2(NULL); Int_t jh(0);
+ if(HasDump3D()){
+ TDirectory *cwd = gDirectory;
+ TFile::Open(Form("EffDump_%s.root", H->GetName()), "RECREATE");
+ for(Int_t ip(0); ip<php.GetEntriesFast(); ip++){
+ if(!(pr0 = (AliTRDrecoProjection*)php.At(ip))) continue;
+ if(!pr0->H()) continue;
+ TH3 *h3=(TH3*)pr0->H()->Clone();
+ h3->Write();
+ }
+ gFile->Close();
+ cwd->cd();
+ }
+
+ Int_t jh(0);
for(; ih--; ){
if(!hp[ih].H()) continue;
- hp[ih].Projection2D(1, 10, -1, kFALSE);
- if((h2 = (TH2*)gDirectory->Get(Form("%sEn", hp[ih].H()->GetName())))) fProj->AddAt(h2, jh++);
+ for(Int_t ipt(0); ipt<=fNpt; ipt++) fProj->AddAt(Projection2D(hp[ih].H(), ipt), jh++);
}
- AliTRDrecoProjection *pr0(NULL), *pr1(NULL);
- AliTRDrecoProjection prLab; TH2 *hLab[3] = {0}; TH1 *hpLab[3] = {0};
+/* AliTRDrecoProjection prLab; TH2 *hLab[3] = {0}; TH1 *hpLab[3] = {0};
for(ilab=0; ilab<nlab; ilab++){
if(!(pr0 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", ilab, 3)))) continue;
prLab=(*pr0);
prLab.Projection2D(1, 10, -1, kFALSE);
if((hLab[ilab] = (TH2*)gDirectory->Get(Form("%sEn", prLab.H()->GetName())))) fProj->AddAt(hLab[ilab], jh++);
if((hpLab[ilab] = prLab.H()->Project3D("z"))) fProj->AddAt(hpLab[ilab], jh++);
- }
+ }*/
for(Int_t istat(0); istat<5; istat++) {
if((pr0 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", 0, istat)))) {
+ // sum over MC labels if available
for(ilab=1; ilab<nlab; ilab++){
if(!(pr1 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", ilab, istat)))) continue;
(*pr0)+=(*pr1);
}
pr0->H()->SetNameTitle(Form("HEff%d", istat), Form("Efficiency :: Stat[#bf{%s}]", stat[istat]));
- pr0->Projection2D(1, 10, -1, kFALSE);
- if((h2 = (TH2*)gDirectory->Get(Form("%sEn", pr0->H()->GetName())))) fProj->AddAt(h2, jh++);
+ for(Int_t ipt(0); ipt<=fNpt; ipt++) fProj->AddAt(Projection2D(pr0->H(), ipt), jh++);
if(istat>1 && (pr1 = (AliTRDrecoProjection*)php.FindObject("HEff01"))) (*pr1)+=(*pr0);
if(istat>2 && (pr1 = (AliTRDrecoProjection*)php.FindObject("HEff02"))) (*pr1)+=(*pr0);
if(istat>3 && (pr1 = (AliTRDrecoProjection*)php.FindObject("HEff03"))) (*pr1)+=(*pr0);
}
}
- // All tracks
- TH2 *hEff[3] = {0};TH1 *hpEff[3] = {0};
- if((pr0 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", 0, 1)))) {
- pr0->H()->SetNameTitle("HEff", "Efficiency :: All Tracks");
- pr0->Projection2D(1, 10, -1, kFALSE);
- hEff[0] = (TH2*)gDirectory->Get(Form("%sEn", pr0->H()->GetName()));
- hpEff[0]= pr0->H()->Project3D("z");
- }
- // Tracked tracks
- if((pr0 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", 0, 2)))) {
- pr0->H()->SetNameTitle("H2EffT", "Efficiency :: Tracked Tracks");
- pr0->Projection2D(1, 10, -1, kFALSE);
- hEff[1] = (TH2*)gDirectory->Get(Form("%sEn", pr0->H()->GetName()));
- hpEff[1]= pr0->H()->Project3D("z");
- }
- // Propagated tracks
- if((pr0 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", 0, 3)))) {
- pr0->H()->SetNameTitle("HEffPrp", "Efficiency :: Propagated Tracks");
- pr0->Projection2D(1, 10, -1, kFALSE);
- hEff[2] = (TH2*)gDirectory->Get(Form("%sEn", pr0->H()->GetName()));
- hpEff[2]= pr0->H()->Project3D("z");
- }
- if(hEff[0]){
- if(hEff[1]){
- hEff[1]->Divide(hEff[0]);
- fProj->AddAt(hEff[1], jh++);
- }
- if(hEff[2]){
- TH2 *hEff1 = (TH2*)hEff[2]->Clone("H2EffPEn");
- hEff1->Divide(hEff[0]);
- fProj->AddAt(hEff1, jh++);
- }
+ // Project 2D tracks
+ const char suffix[] = {'A', 'T', 'P'};
+ const char *sname[] = {"All", "Trk", "Prp"};
+ for(Int_t istat(0); istat<3; istat++){
+ if(!(pr0 = (AliTRDrecoProjection*)php.FindObject(Form("HEff%d%d", 0, istat+1)))) continue;
+ pr0->H()->SetNameTitle(Form("HEff%c", suffix[istat]), Form("Efficiency :: %s Tracks", sname[istat]));
+ for(Int_t ipt(-1); ipt<=fNpt; ipt++) fProj->AddAt(Projection2D(pr0->H(), ipt), jh++);
}
- if(hpEff[0]){
- if(hpEff[1]){
- hpEff[1]->Divide(hpEff[0]);
- fProj->AddAt(hpEff[1], jh++);
- }
- if(hEff[2]){
- TH1 *hpEff1 = (TH1*)hpEff[2]->Clone("H2EffP_z");
- hpEff1->Divide(hpEff[0]);
- fProj->AddAt(hpEff1, jh++);
- }
+
+ // Efficiency
+ TH2F *h2T(NULL), *h2P(NULL);
+ for(Int_t ipt(-1); ipt<=fNpt; ipt++){
+ if(!(h2T = (TH2F*)fProj->FindObject(Form("HEffT%d_2D", ipt)))) continue;
+ if(!(h2P = (TH2F*)fProj->FindObject(Form("HEffP%d_2D", ipt)))) continue;
+ h2P->Divide(h2T);
+ PutTrendValue(ipt<0?"pt":(Form("pt%d", ipt)), GetMeanStat(h2P, 0.01, ">"));
}
- // process MC label
+/* // process MC label
if(hEff[2]){
for(ilab=0; ilab<nlab; ilab++){
if(!hLab[ilab]) continue;
- hLab[ilab]->Divide(hEff[2]);
- fProj->AddAt(hLab[ilab], jh++);
+ // remove fakes
+ TH2 *hEff1 = (TH2*)hLab[ilab]->Clone(Form("%sN", hLab[ilab]->GetName()));
+ for(Int_t ix(1); ix<=hLab[ilab]->GetNbinsX(); ix++){
+ for(Int_t iy(1); iy<=hLab[ilab]->GetNbinsY(); iy++){
+ if(hLab[ilab]->GetBinContent(ix, iy)<5) hEff1->SetBinContent(ix, iy, 0.);
+ }}
+ hEff1->Divide(hEff[2]);
+ fProj->AddAt(hEff1, jh++);
}
}
if(hpEff[2]){
for(ilab=0; ilab<nlab; ilab++){
if(!hpLab[ilab]) continue;
- hpLab[ilab]->Divide(hpEff[2]);
- fProj->AddAt(hpLab[ilab], jh++);
+ TH1 *hpEff1 = (TH1*)hpLab[ilab]->Clone(Form("%sN", hpLab[ilab]->GetName()));
+ hpEff1->Divide(hpEff[2]);
+ fProj->AddAt(hpEff1, jh++);
}
- }
+ }*/
AliInfo(Form("Done %3d 2D projections.", jh));
return kTRUE;
}
cOut = new TCanvas(Form("%s_Eff", GetName()), "TRD Efficiency", 1536, 1536); cOut->Divide(2,2,1.e-5,1.e-5);
// tracking eff :: eta-phi dependence
for(Int_t it(0); it<2; it++){
- if(!(h2 = (TH2*)fProj->FindObject(Form("H2Eff%cEn", cid[it])))) {
- AliError(Form("Missing \"H2Eff%c\".", cid[it]));
+ if(!(h2 = (TH2*)fProj->FindObject(Form("HEff%cEnN", cid[it])))) {
+ AliError(Form("Missing \"HEff%cEnN\".", cid[it]));
continue;
}
h2->SetContour(10); h2->Scale(1.e2); SetRangeZ(h2, 80, 100, 5);
h2->Draw("colz"); MakeDetectorPlot();
// tracking eff :: pt dependence
TH1 *h[2] = {0};
- if(!(h[0] = (TH1*)fProj->FindObject("H2EffP_z"))){
- AliError("Missing \"H2EffP_z\".");
+ if(!(h[0] = (TH1*)fProj->FindObject("HEffP_zN"))){
+ AliError("Missing \"HEffP_zN\".");
return;
}
- if(!(h[1] = (TH1*)fProj->FindObject("H2EffT_z"))){
- AliError("Missing \"H2EffT_z\".");
+ if(!(h[1] = (TH1*)fProj->FindObject("HEffT_zN"))){
+ AliError("Missing \"HEffT_zN\".");
return;
}
TH1 *h1[3] = {0};
h1[il]->SetFillColor(color[il]);
h1[il]->SetFillStyle(il==2?3002:1001);
h1[il]->SetLineColor(color[il]);
+ h1[il]->SetMarkerStyle(4);
+ h1[il]->SetMarkerColor(color[il]);
h1[il]->SetLineWidth(1);
}
for(Int_t ip(0);ip<=(nbins+1);ip++){
leg->SetBorderSize(0); leg->SetFillColor(kWhite); leg->SetFillStyle(1001);
for(Int_t ic(0); ic<3;ic++){ hs->Add(h1[ic]);leg->AddEntry(h1[ic], labEff[ic], "f");}
p=cOut->cd(4); p->SetLeftMargin(0.08266129); p->SetRightMargin(0.0141129);p->SetTopMargin(0.006355932);p->SetLogx();
- hs->Draw(); leg->Draw();
+ hs->Draw(); //hs->Draw("same nostack,e1p");
+ leg->Draw();
hs->GetXaxis()->SetRangeUser(0.6,10.);
hs->GetXaxis()->SetMoreLogLabels();
hs->GetXaxis()->SetTitleOffset(1.2);
hs->GetYaxis()->SetNdivisions(513);
- hs->SetMinimum(80.);
+ hs->SetMinimum(75.);
hs->GetYaxis()->CenterTitle();
cOut->SaveAs(Form("%s.gif", cOut->GetName()));
cOut = new TCanvas(Form("%s_MC", GetName()), "TRD Label", 1536, 1536); cOut->Divide(2,2,1.e-5,1.e-5);
for(Int_t ipad(0); ipad<3; ipad++){
p=cOut->cd(ipad+1);p->SetRightMargin(0.1572581);p->SetTopMargin(0.08262712);
- if(!(h2 = (TH2*)fProj->FindObject(Form("HEffLb%dEn", ipad)))) continue;
+ if(!(h2 = (TH2*)fProj->FindObject(Form("HEffLb%dEnN", ipad)))) continue;
h2->SetContour(10);
- h2->Scale(1.e2); SetRangeZ(h2, ipad==1?80:0., ipad==1?100.:10., ipad==1?30:0.01);
+ h2->Scale(1.e2); SetRangeZ(h2, ipad==1?50:0., ipad==1?90.:50., ipad==1?0.01:0.01);
h2->GetZaxis()->SetTitle("Efficiency [%]"); h2->GetZaxis()->CenterTitle();
h2->Draw("colz");
MakeDetectorPlot();
}
color[0] = kRed; color[1] = kGreen; color[2] = kBlue;
for(Int_t il=0;il<3;il++){
- if(!(h[il] = (TH1D*)fProj->FindObject(Form("HEffLb%d_z", il)))) continue;
+ if(!(h[il] = (TH1D*)fProj->FindObject(Form("HEffLb%d_zN", il)))) continue;
h1[il]=new TH1F(Form("h1Lab%d", il), "", nbins+2, ptBins);
for(Int_t ip(0);ip<=(nbins+1);ip++) h1[il]->SetBinContent(ip+1, 1.e2*h[il]->GetBinContent(ip));
h1[il]->SetFillColor(color[il]);
hs->Add(h1[2]);leg->AddEntry(h1[2], labMC[2], "f"); // accept
hs->Add(h1[0]);leg->AddEntry(h1[0], labMC[0], "f"); // bad
p=cOut->cd(4); p->SetLeftMargin(0.08266129); p->SetRightMargin(0.0141129);p->SetTopMargin(0.006355932); p->SetLogx();
- hs->Draw(); leg->Draw();
+ hs->Draw(/*"nostack,e1p"*/); leg->Draw();
cOut->Modified();cOut->Update();
hs->GetXaxis()->SetRangeUser(0.6,10.);
hs->GetXaxis()->SetMoreLogLabels();
hs->GetXaxis()->SetTitleOffset(1.2);
hs->GetYaxis()->SetNdivisions(513);
- hs->SetMinimum(80.);
+ hs->SetMinimum(50.);
hs->GetYaxis()->CenterTitle();
cOut->SaveAs(Form("%s.gif", cOut->GetName()));
}
+
+
+//____________________________________________________________________
+TH2* AliTRDefficiency::Projection2D(TH3 *h3, Int_t ipt)
+{
+ TAxis *ax(h3->GetXaxis()), *ay(h3->GetYaxis());
+ TH2F *h2(NULL);
+ if(ipt<0){
+ h2 =(TH2F*)h3->Project3D("yx");
+ h2->SetNameTitle(Form("%s%d_2D", h3->GetName(), ipt), h3->GetTitle());
+ } else {
+ h2 = new TH2F(Form("%s%d_2D", h3->GetName(), ipt),
+ Form("%s | #it{%4.2f<=p_{t}[GeV/c]<%4.2f};%s;%s;Entries", h3->GetTitle(),
+ ipt?fgPt[ipt-1]:0., ipt==fNpt?9.99:fgPt[ipt], ax->GetTitle(), ay->GetTitle()),
+ ax->GetNbins(), ax->GetXmin(), ax->GetXmax(),
+ ay->GetNbins(), ay->GetXmin(), ay->GetXmax());
+ for(Int_t ix(1); ix<=ax->GetNbins(); ix++){
+ for(Int_t iy(1); iy<=ay->GetNbins(); iy++){
+ h2->SetBinContent(ix, iy, h3->GetBinContent(ix, iy, ipt));
+ }
+ }
+ }
+ return h2;
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff