[u/mrichter/AliRoot.git] / ITS / macrosSDD / CheckSDDInESD.C

#if !defined(__CINT__) || defined(__MAKECINT__)
#include <TH1F.h>
#include <TH2F.h>
#include <TF1.h>
#include <TPad.h>
#include <TGraphErrors.h>
#include <TROOT.h>
#include <TFile.h>
#include <TTree.h>
#include <TMath.h>
#include <TCanvas.h>
#include <TStyle.h>
#include <TLatex.h>
#include "AliITSgeomTGeo.h"
#include "AliESDEvent.h"
#endif

enum {kAll, kTPCITS, kITSsa, kITSpureSA};

void CheckSDDInESD(TString filename="AliESDs.root", Int_t optTracks=kAll){


  TFile* esdFile = TFile::Open(filename.Data());
  if (!esdFile || !esdFile->IsOpen()) {
    printf("Error in opening ESD file");
    return;
  }

  AliESDEvent * esd = new AliESDEvent;
  TTree* tree = (TTree*) esdFile->Get("esdTree");
  if (!tree) {
    printf("Error: no ESD tree found");
    return;
  }
  esd->ReadFromTree(tree);
  TH1F* hpt=new TH1F("hpt","",100,0.,10.);
  TH1F* hphi=new TH1F("hphi","",100,-1,1);
  TH1F* hlam=new TH1F("hlam","",100,-2.,2.);
  TH1F* halpha=new TH1F("halpha","",100,-7,7);
  TH1F* hitscl=new TH1F("hitscl","",7,-0.5,6.5);
  TH1F* htpccl=new TH1F("htpccl","",200,-0.5,199.5);
  TH1F* hitsmap=new TH1F("hitsmap","",64,-0.5,63.5);
  TH1F* hclulay=new TH1F("hclulay","",7,-1.5,5.5);

  TH1F* hvx=new TH1F("hvx","",100,-1.,1.);
  TH1F* hvy=new TH1F("hvy","",100,-1.,1.);
  TH1F* hvz=new TH1F("hvz","",100,-20.,20.);
  TH1F* hdedx3=new TH1F("hdedx3","",100,0.,300.);
  TH1F* hdedx4=new TH1F("hdedx4","",100,0.,300.);
  TH1F* hdedx5=new TH1F("hdedx5","",100,0.,300.);
  TH1F* hdedx6=new TH1F("hdedx6","",100,0.,300.);
  TH1F* hStatus=new TH1F("hStatus","",11,-1.5,9.5);


  // -- Local coordinates


  // -- Module histos

  TH1F* hAllPMod  = new TH1F("hAllPmod","Crossing Tracks vs. Module",260,239.5,499.5);
  TH1F* hGoodPMod  = new TH1F("hGoodPmod","PointsAssocToTrack per Module",260,239.5,499.5);
  TH1F* hBadRegMod  = new TH1F("hBadRegmod","Tracks in BadRegion per Module",260,239.5,499.5);
  TH1F* hMissPMod  = new TH1F("hMissPmod","Missing Points per Module",260,239.5,499.5);
  TH1F* hSkippedMod  = new TH1F("hSkippedmod","Tracks in Skipped Module",260,239.5,499.5);
  TH1F* hOutAccMod  = new TH1F("hOutAccmod","Tracks outside zAcc per Module",260,239.5,499.5);
  TH1F* hNoRefitMod  = new TH1F("hNoRefitmod","Points rejected in refit per Module",260,239.5,499.5);

  TH1F* hAllPXloc  = new TH1F("hAllPxloc","Crossing Tracks vs. Xloc",75, -3.75, 3.75);
  TH1F* hGoodPXloc  = new TH1F("hGoodPxloc","PointsAssocToTrack vs. Xloc",75, -3.75, 3.75);
  TH1F* hBadRegXloc  = new TH1F("hBadRegxloc","Tracks in BadRegion vs. Xloc",75, -3.75, 3.75);
  TH1F* hMissPXloc  = new TH1F("hMissPxloc","Missing Points vs. Xloc",75, -3.75, 3.75);
  TH1F* hAllPZloc  = new TH1F("hAllPzloc","Crossing Tracks vs. Zloc",77, -3.85, 3.85);
  TH1F* hGoodPZloc  = new TH1F("hGoodPzloc","PointsAssocToTrack vs. Zloc",77, -3.85, 3.85);
  TH1F* hBadRegZloc  = new TH1F("hBadRegzloc","Tracks in BadRegion vs. Zloc",77, -3.85, 3.85);
  TH1F* hMissPZloc  = new TH1F("hMissPzloc","Missing Points vs. Zloc",77, -3.85, 3.85);
  TH2F* hdEdxVsMod=new TH2F("hdEdxVsMod","dE/dx vs. mod",260,239.5,499.5,100,0.,500.);

  gStyle->SetPalette(1);
  

  for (Int_t iEvent = 0; iEvent < tree->GetEntries(); iEvent++) {
    tree->GetEvent(iEvent);
    if (!esd) {
      printf("Error: no ESD object found for event %d", iEvent);
      return;
    }
    cout<<"-------- Event "<<iEvent<<endl;
    printf(" Tracks # = %d\n",esd->GetNumberOfTracks());
    const AliESDVertex *spdv=esd->GetVertex();
    printf(" SPD Primary Vertex in %f %f %f with %d contributors\n",spdv->GetX(),spdv->GetY(),spdv->GetZ(),spdv->GetNContributors());
    const AliESDVertex *trkv=esd->GetPrimaryVertex();
    printf(" Track Primary Vertex with %d contributors\n",trkv->GetNContributors());
    if(spdv->IsFromVertexer3D()){
      hvx->Fill(spdv->GetX());
      hvy->Fill(spdv->GetY());
      hvz->Fill(spdv->GetZ());
    }
    Double_t itss[4];
    for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
      AliESDtrack* track = esd->GetTrack(iTrack);
      Int_t nITSclus=track->GetNcls(0);
      UChar_t clumap=track->GetITSClusterMap();
      Int_t nPointsForPid=0;
      for(Int_t i=2; i<6; i++){
	if(clumap&(1<<i)) ++nPointsForPid;
      }
      //      track->PropagateTo(4.,5.);
      htpccl->Fill(track->GetNcls(1));
      Int_t status=track->GetStatus();
      Bool_t tpcin=0;
      hStatus->Fill(-1.);
      if(status & AliESDtrack::kTPCin){
	tpcin=1;
	hStatus->Fill(0.);
      }
      if(status & AliESDtrack::kTPCout){
	hStatus->Fill(1.);
      }
      if(status & AliESDtrack::kTPCrefit){
	hStatus->Fill(2.);
      }
      Bool_t itsin=0;
      if(status & AliESDtrack::kITSin){
	itsin=1;
	hStatus->Fill(3.);
      }
      if(status & AliESDtrack::kITSout){
	hStatus->Fill(4.);
      }
      if(status & AliESDtrack::kITSrefit){
	hStatus->Fill(5.);
      }
      if(!tpcin && itsin){
	hStatus->Fill(6.);
      }
      if(status & AliESDtrack::kITSpureSA){
	hStatus->Fill(7.);
      }

      if(status & AliESDtrack::kITSrefit){
	if((optTracks==kTPCITS) && !(status & AliESDtrack::kTPCin)) continue;
	if((optTracks==kITSsa) && (status & AliESDtrack::kTPCin)) continue;
	if((optTracks==kITSsa) && (status & AliESDtrack::kITSpureSA)) continue;
	if((optTracks==kITSpureSA) && (status & AliESDtrack::kITSpureSA)) continue;

	 track->GetITSdEdxSamples(itss);
	//	printf("Track %d (label %d) in ITS with %d clusters clumap %d pointspid= %d\n",iTrack,track->GetLabel(),nITSclus,clumap,nPointsForPid);
	//printf("   dedx=%f %f %f %f\n",itss[0],itss[1],itss[2],itss[3]);
	hitscl->Fill(nITSclus);
	hdedx3->Fill(itss[0]);
	hdedx4->Fill(itss[1]);
	hdedx5->Fill(itss[2]);
	hdedx6->Fill(itss[3]);
	hitsmap->Fill(clumap);
	hclulay->Fill(-1.);
	for(Int_t iLay=0;iLay<6;iLay++){
	  if(clumap&1<<iLay) hclulay->Fill(iLay);
	}
	hpt->Fill(track->Pt());
	hphi->Fill(TMath::ASin(track->GetSnp()));
	hlam->Fill(TMath::ATan(track->GetTgl()));
	halpha->Fill(track->GetAlpha());
	Int_t iMod,status;
	Float_t xloc,zloc;
	for(Int_t iLay=2; iLay<=3; iLay++){
	  Bool_t ok=track->GetITSModuleIndexInfo(iLay,iMod,status,xloc,zloc);
	  if(ok){
	    iMod+=240;
	    hAllPMod->Fill(iMod);
	    hAllPXloc->Fill(xloc);
	    hAllPZloc->Fill(zloc);
	    if(status==1){
	      hGoodPMod->Fill(iMod);
	      hGoodPXloc->Fill(xloc);
	      hGoodPZloc->Fill(zloc);
	      if(track->Pt()>1.) hdEdxVsMod->Fill(iMod,itss[iLay-2]);
	    }
	    else if(status==2){ 
	      hBadRegMod->Fill(iMod);
	      hBadRegXloc->Fill(xloc);
	      hBadRegZloc->Fill(zloc);
	    }
	    else if(status==3) hSkippedMod->Fill(iMod);
	    else if(status==4) hOutAccMod->Fill(iMod);
	    else if(status==5){
	      hMissPMod->Fill(iMod);
	      hMissPXloc->Fill(xloc);
	      hMissPZloc->Fill(zloc);
	    }
	    else if(status==6) hNoRefitMod->Fill(iMod);
	  }
	}
      }
    }
  }
  Float_t norm=hclulay->GetBinContent(1);
  if(norm<1.) norm=1.;
  hclulay->Scale(1./norm);
  gStyle->SetLineWidth(2);

  TCanvas* c1=new TCanvas("c1","Track quantities",900,900);
  c1->Divide(2,2);
  c1->cd(1);
  htpccl->Draw();
  htpccl->GetXaxis()->SetTitle("Clusters in TPC ");
  c1->cd(2);
  hitscl->Draw();
  hitscl->GetXaxis()->SetTitle("Clusters in ITS ");
  c1->cd(3);
  hclulay->Draw();
  hclulay->GetXaxis()->SetRange(2,7);
  hclulay->GetXaxis()->SetTitle("# ITS Layer");
  hclulay->GetYaxis()->SetTitle("Fraction of tracks with point in Layer x");
  c1->cd(4);

  TCanvas* c2=new TCanvas("c2","dedx per Layer",900,900);
  c2->Divide(2,2);
  c2->cd(1);
  hdedx3->Draw();
  hdedx3->GetXaxis()->SetTitle("dE/dx Lay3");
  c2->cd(2);
  hdedx4->Draw();
  hdedx4->GetXaxis()->SetTitle("dE/dx Lay4");
  c2->cd(3);
  hdedx5->Draw();
  hdedx5->GetXaxis()->SetTitle("dE/dx Lay5");
  c2->cd(4);
  hdedx6->Draw();
  hdedx6->GetXaxis()->SetTitle("dE/dx Lay6");

  hdEdxVsMod->SetStats(0);
  TCanvas* cdedx=new TCanvas("cdedx","dedx SDD",1400,600);
  cdedx->SetLogz();
  hdEdxVsMod->Draw("col"); 
  hdEdxVsMod->GetXaxis()->SetTitle("SDD Module Id");
  hdEdxVsMod->GetYaxis()->SetTitle("dE/dx (keV/300 #mum)");
  hdEdxVsMod->GetYaxis()->SetTitleOffset(1.25);


  TCanvas* cv=new TCanvas("cv","Vertex",600,900);
  cv->Divide(1,3);
  cv->cd(1);
  hvx->Draw();
  hvx->GetXaxis()->SetTitle("Xv (cm)");
  cv->cd(2);
  hvy->Draw();
  hvy->GetXaxis()->SetTitle("Yv (cm)");
  cv->cd(3);
  hvz->Draw();
  hvz->GetXaxis()->SetTitle("Xv (cm)");

  hGoodPMod->SetStats(0);
  hGoodPMod->SetTitle("");
  TCanvas* ceff0=new TCanvas("ceff0","ModuleIndexInfo",1000,600);
  hGoodPMod->Draw("e");
  hGoodPMod->GetXaxis()->SetTitle("SDD Module Id");
  hGoodPMod->GetYaxis()->SetTitle("Number of tracks");
  hMissPMod->SetLineColor(2);
  hMissPMod->SetMarkerColor(2);
  hMissPMod->SetMarkerStyle(22);
  hMissPMod->SetMarkerSize(0.5);
  hMissPMod->Draw("psame");
  hBadRegMod->SetLineColor(kGreen+1);
  hBadRegMod->SetMarkerColor(kGreen+1);
  hBadRegMod->SetMarkerStyle(20);
  hBadRegMod->SetMarkerSize(0.5);
  hBadRegMod->Draw("esame");
  hSkippedMod->SetLineColor(kYellow);
  hSkippedMod->Draw("esame");
  hOutAccMod->SetLineColor(4);
  hOutAccMod->Draw("esame");
  hNoRefitMod->SetLineColor(6);
  hNoRefitMod->Draw("esame");
  TLatex* t1=new TLatex(0.7,0.85,"Good Point");
  t1->SetNDC();
  t1->SetTextColor(1);
  t1->Draw();
  TLatex* t2=new TLatex(0.7,0.8,"Missing Point");
  t2->SetNDC();
  t2->SetTextColor(2);
  t2->Draw();
  TLatex* t3=new TLatex(0.7,0.75,"Bad Region");
  t3->SetNDC();
  t3->SetTextColor(kGreen+1);
  t3->Draw();
  ceff0->Update();

  TH1F* heff=new TH1F("heff","",260,239.5,499.5);
  for(Int_t imod=0; imod<260;imod++){
    Float_t numer=hGoodPMod->GetBinContent(imod+1)+hBadRegMod->GetBinContent(imod+1)+hOutAccMod->GetBinContent(imod+1)+hNoRefitMod->GetBinContent(imod+1);
    Float_t denom=hAllPMod->GetBinContent(imod+1);
    Float_t eff=0.;
    Float_t erreff=0.;
    if(denom>0){
      eff=numer/denom;
      erreff=TMath::Sqrt(eff*(1-eff)/denom);
    }
    heff->SetBinContent(imod+1,eff);
    heff->SetBinError(imod+1,erreff);
  }

  printf("---- Modules with efficiency < 90%% ----\n");
  heff->SetStats(0);
  TCanvas* ceff1=new TCanvas("ceff1","Efficiency",1000,600);
  heff->Draw();
  heff->GetXaxis()->SetTitle("SDD Module Id");
  heff->GetYaxis()->SetTitle("Fraction of tracks with point in good region");
  for(Int_t ibin=1; ibin<=heff->GetNbinsX(); ibin++){
    Float_t e=heff->GetBinContent(ibin);
    if(e<0.9){
      Int_t iMod=(Int_t)heff->GetBinCenter(ibin);
      Int_t lay,lad,det;
      AliITSgeomTGeo::GetModuleId(iMod,lay,lad,det);
      printf("Module %d - Layer %d Ladder %2d Det %d  -   Eff. %.3f\n",iMod,lay,lad,det,heff->GetBinContent(ibin));
    }
  }
  ceff1->Update();


  hGoodPXloc->SetTitle("");
  hGoodPZloc->SetTitle("");
  hGoodPXloc->SetStats(0);
  hGoodPZloc->SetStats(0);
  hGoodPXloc->SetMinimum(0);
  hGoodPZloc->SetMinimum(0);
  TCanvas* ceff2=new TCanvas("ceff2","LocalCoord",1000,600);
  ceff2->Divide(2,1);
  ceff2->cd(1);
  hGoodPXloc->Draw("e");
  hGoodPXloc->GetXaxis()->SetTitle("Xlocal (cm)");
  hGoodPXloc->GetYaxis()->SetTitle("Number of tracks");
  hMissPXloc->SetLineColor(2);
  hMissPXloc->SetMarkerColor(2);
  hMissPXloc->SetMarkerStyle(22);
  hMissPXloc->SetMarkerSize(0.5);
  hMissPXloc->Draw("psame");
  hBadRegXloc->SetLineColor(kGreen+1);
  hBadRegXloc->SetMarkerColor(kGreen+1);
  hBadRegXloc->SetMarkerStyle(20);
  hBadRegXloc->SetMarkerSize(0.5);
  hBadRegXloc->Draw("psame");
  t1->Draw();
  t2->Draw();
  t3->Draw();
  ceff2->cd(2);
  hGoodPZloc->Draw("e");
  hGoodPZloc->GetXaxis()->SetTitle("Zlocal (cm)");
  hGoodPZloc->GetYaxis()->SetTitle("Number of tracks");
  hMissPZloc->SetLineColor(2);
  hMissPZloc->SetMarkerColor(2);
  hMissPZloc->SetMarkerStyle(22);
  hMissPZloc->SetMarkerSize(0.5);
  hMissPZloc->Draw("psame");
  hBadRegZloc->SetLineColor(kGreen+1);
  hBadRegZloc->SetMarkerColor(kGreen+1);
  hBadRegZloc->SetMarkerStyle(20);
  hBadRegZloc->SetMarkerSize(0.5);
  hBadRegZloc->Draw("psame");
  t1->Draw();
  t2->Draw();
  t3->Draw();
  ceff2->Update();


}
Commit	Line	Data
db5391b9	1	#if !defined(__CINT__) \|\| defined(__MAKECINT__)
	2	#include <TH1F.h>
	3	#include <TH2F.h>
	4	#include <TF1.h>
	5	#include <TPad.h>
	6	#include <TGraphErrors.h>
	7	#include <TROOT.h>
	8	#include <TFile.h>
	9	#include <TTree.h>
	10	#include <TMath.h>
	11	#include <TCanvas.h>
	12	#include <TStyle.h>
	13	#include <TLatex.h>
	14	#include "AliITSgeomTGeo.h"
	15	#include "AliESDEvent.h"
	16	#endif
	17
	18	enum {kAll, kTPCITS, kITSsa, kITSpureSA};
	19
	20	void CheckSDDInESD(TString filename="AliESDs.root", Int_t optTracks=kAll){
	21
	22
	23	TFile* esdFile = TFile::Open(filename.Data());
	24	if (!esdFile \|\| !esdFile->IsOpen()) {
	25	printf("Error in opening ESD file");
	26	return;
	27	}
	28
	29	AliESDEvent * esd = new AliESDEvent;
	30	TTree* tree = (TTree*) esdFile->Get("esdTree");
	31	if (!tree) {
	32	printf("Error: no ESD tree found");
	33	return;
	34	}
	35	esd->ReadFromTree(tree);
	36	TH1F* hpt=new TH1F("hpt","",100,0.,10.);
	37	TH1F* hphi=new TH1F("hphi","",100,-1,1);
	38	TH1F* hlam=new TH1F("hlam","",100,-2.,2.);
	39	TH1F* halpha=new TH1F("halpha","",100,-7,7);
	40	TH1F* hitscl=new TH1F("hitscl","",7,-0.5,6.5);
	41	TH1F* htpccl=new TH1F("htpccl","",200,-0.5,199.5);
	42	TH1F* hitsmap=new TH1F("hitsmap","",64,-0.5,63.5);
	43	TH1F* hclulay=new TH1F("hclulay","",7,-1.5,5.5);
	44
	45	TH1F* hvx=new TH1F("hvx","",100,-1.,1.);
	46	TH1F* hvy=new TH1F("hvy","",100,-1.,1.);
	47	TH1F* hvz=new TH1F("hvz","",100,-20.,20.);
	48	TH1F* hdedx3=new TH1F("hdedx3","",100,0.,300.);
	49	TH1F* hdedx4=new TH1F("hdedx4","",100,0.,300.);
	50	TH1F* hdedx5=new TH1F("hdedx5","",100,0.,300.);
	51	TH1F* hdedx6=new TH1F("hdedx6","",100,0.,300.);
	52	TH1F* hStatus=new TH1F("hStatus","",11,-1.5,9.5);
	53
	54
	55	// -- Local coordinates
	56
	57
	58	// -- Module histos
	59
	60	TH1F* hAllPMod = new TH1F("hAllPmod","Crossing Tracks vs. Module",260,239.5,499.5);
	61	TH1F* hGoodPMod = new TH1F("hGoodPmod","PointsAssocToTrack per Module",260,239.5,499.5);
	62	TH1F* hBadRegMod = new TH1F("hBadRegmod","Tracks in BadRegion per Module",260,239.5,499.5);
	63	TH1F* hMissPMod = new TH1F("hMissPmod","Missing Points per Module",260,239.5,499.5);
	64	TH1F* hSkippedMod = new TH1F("hSkippedmod","Tracks in Skipped Module",260,239.5,499.5);
65	TH1F* hOutAccMod = new TH1F("hOutAccmod","Tracks outside zAcc per Module",260,239.5,499.5);
66	TH1F* hNoRefitMod = new TH1F("hNoRefitmod","Points rejected in refit per Module",260,239.5,499.5);
67
68	TH1F* hAllPXloc = new TH1F("hAllPxloc","Crossing Tracks vs. Xloc",75, -3.75, 3.75);
69	TH1F* hGoodPXloc = new TH1F("hGoodPxloc","PointsAssocToTrack vs. Xloc",75, -3.75, 3.75);
70	TH1F* hBadRegXloc = new TH1F("hBadRegxloc","Tracks in BadRegion vs. Xloc",75, -3.75, 3.75);
71	TH1F* hMissPXloc = new TH1F("hMissPxloc","Missing Points vs. Xloc",75, -3.75, 3.75);
72	TH1F* hAllPZloc = new TH1F("hAllPzloc","Crossing Tracks vs. Zloc",77, -3.85, 3.85);
73	TH1F* hGoodPZloc = new TH1F("hGoodPzloc","PointsAssocToTrack vs. Zloc",77, -3.85, 3.85);
74	TH1F* hBadRegZloc = new TH1F("hBadRegzloc","Tracks in BadRegion vs. Zloc",77, -3.85, 3.85);
75	TH1F* hMissPZloc = new TH1F("hMissPzloc","Missing Points vs. Zloc",77, -3.85, 3.85);
76	TH2F* hdEdxVsMod=new TH2F("hdEdxVsMod","dE/dx vs. mod",260,239.5,499.5,100,0.,500.);
77
78	gStyle->SetPalette(1);
79
80
81	for (Int_t iEvent = 0; iEvent < tree->GetEntries(); iEvent++) {
82	tree->GetEvent(iEvent);
83	if (!esd) {
84	printf("Error: no ESD object found for event %d", iEvent);
85	return;
86	}
87	cout<<"-------- Event "<<iEvent<<endl;
88	printf(" Tracks # = %d\n",esd->GetNumberOfTracks());
89	const AliESDVertex *spdv=esd->GetVertex();
e690d4d0	90	printf(" SPD Primary Vertex in %f %f %f with %d contributors\n",spdv->GetX(),spdv->GetY(),spdv->GetZ(),spdv->GetNContributors());
db5391b9	91	const AliESDVertex *trkv=esd->GetPrimaryVertex();
	92	printf(" Track Primary Vertex with %d contributors\n",trkv->GetNContributors());
	93	if(spdv->IsFromVertexer3D()){
e690d4d0	94	hvx->Fill(spdv->GetX());
	95	hvy->Fill(spdv->GetY());
	96	hvz->Fill(spdv->GetZ());
db5391b9	97	}
	98	Double_t itss[4];
	99	for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
	100	AliESDtrack* track = esd->GetTrack(iTrack);
	101	Int_t nITSclus=track->GetNcls(0);
	102	UChar_t clumap=track->GetITSClusterMap();
	103	Int_t nPointsForPid=0;
	104	for(Int_t i=2; i<6; i++){
	105	if(clumap&(1<<i)) ++nPointsForPid;
	106	}
	107	// track->PropagateTo(4.,5.);
90af9cac	108	htpccl->Fill(track->GetNcls(1));
db5391b9	109	Int_t status=track->GetStatus();
	110	Bool_t tpcin=0;
	111	hStatus->Fill(-1.);
	112	if(status & AliESDtrack::kTPCin){
	113	tpcin=1;
	114	hStatus->Fill(0.);
	115	}
	116	if(status & AliESDtrack::kTPCout){
	117	hStatus->Fill(1.);
	118	}
	119	if(status & AliESDtrack::kTPCrefit){
	120	hStatus->Fill(2.);
	121	}
	122	Bool_t itsin=0;
	123	if(status & AliESDtrack::kITSin){
	124	itsin=1;
	125	hStatus->Fill(3.);
	126	}
	127	if(status & AliESDtrack::kITSout){
	128	hStatus->Fill(4.);
	129	}
	130	if(status & AliESDtrack::kITSrefit){
	131	hStatus->Fill(5.);
	132	}
	133	if(!tpcin && itsin){
	134	hStatus->Fill(6.);
	135	}
	136	if(status & AliESDtrack::kITSpureSA){
	137	hStatus->Fill(7.);
	138	}
	139
	140	if(status & AliESDtrack::kITSrefit){
	141	if((optTracks==kTPCITS) && !(status & AliESDtrack::kTPCin)) continue;
	142	if((optTracks==kITSsa) && (status & AliESDtrack::kTPCin)) continue;
	143	if((optTracks==kITSsa) && (status & AliESDtrack::kITSpureSA)) continue;
	144	if((optTracks==kITSpureSA) && (status & AliESDtrack::kITSpureSA)) continue;
	145
	146	track->GetITSdEdxSamples(itss);
	147	// printf("Track %d (label %d) in ITS with %d clusters clumap %d pointspid= %d\n",iTrack,track->GetLabel(),nITSclus,clumap,nPointsForPid);
	148	//printf(" dedx=%f %f %f %f\n",itss[0],itss[1],itss[2],itss[3]);
	149	hitscl->Fill(nITSclus);
	150	hdedx3->Fill(itss[0]);
	151	hdedx4->Fill(itss[1]);
	152	hdedx5->Fill(itss[2]);
	153	hdedx6->Fill(itss[3]);
	154	hitsmap->Fill(clumap);
	155	hclulay->Fill(-1.);
	156	for(Int_t iLay=0;iLay<6;iLay++){
	157	if(clumap&1<<iLay) hclulay->Fill(iLay);
	158	}
	159	hpt->Fill(track->Pt());
	160	hphi->Fill(TMath::ASin(track->GetSnp()));
	161	hlam->Fill(TMath::ATan(track->GetTgl()));
	162	halpha->Fill(track->GetAlpha());
	163	Int_t iMod,status;
	164	Float_t xloc,zloc;
	165	for(Int_t iLay=2; iLay<=3; iLay++){
	166	Bool_t ok=track->GetITSModuleIndexInfo(iLay,iMod,status,xloc,zloc);
	167	if(ok){
	168	iMod+=240;
	169	hAllPMod->Fill(iMod);
	170	hAllPXloc->Fill(xloc);
	171	hAllPZloc->Fill(zloc);
	172	if(status==1){
173	hGoodPMod->Fill(iMod);
174	hGoodPXloc->Fill(xloc);
175	hGoodPZloc->Fill(zloc);
176	if(track->Pt()>1.) hdEdxVsMod->Fill(iMod,itss[iLay-2]);
177	}
178	else if(status==2){
179	hBadRegMod->Fill(iMod);
180	hBadRegXloc->Fill(xloc);
181	hBadRegZloc->Fill(zloc);
182	}
183	else if(status==3) hSkippedMod->Fill(iMod);
184	else if(status==4) hOutAccMod->Fill(iMod);
185	else if(status==5){
186	hMissPMod->Fill(iMod);
187	hMissPXloc->Fill(xloc);
188	hMissPZloc->Fill(zloc);
189	}
190	else if(status==6) hNoRefitMod->Fill(iMod);
191	}
192	}
193	}
194	}
195	}
196	Float_t norm=hclulay->GetBinContent(1);
197	if(norm<1.) norm=1.;
198	hclulay->Scale(1./norm);
199	gStyle->SetLineWidth(2);
200
201	TCanvas* c1=new TCanvas("c1","Track quantities",900,900);
202	c1->Divide(2,2);
203	c1->cd(1);
204	htpccl->Draw();
205	htpccl->GetXaxis()->SetTitle("Clusters in TPC ");
206	c1->cd(2);
207	hitscl->Draw();
208	hitscl->GetXaxis()->SetTitle("Clusters in ITS ");
209	c1->cd(3);
210	hclulay->Draw();
211	hclulay->GetXaxis()->SetRange(2,7);
212	hclulay->GetXaxis()->SetTitle("# ITS Layer");
213	hclulay->GetYaxis()->SetTitle("Fraction of tracks with point in Layer x");
214	c1->cd(4);
215
216	TCanvas* c2=new TCanvas("c2","dedx per Layer",900,900);
217	c2->Divide(2,2);
218	c2->cd(1);
219	hdedx3->Draw();
220	hdedx3->GetXaxis()->SetTitle("dE/dx Lay3");
221	c2->cd(2);
222	hdedx4->Draw();
223	hdedx4->GetXaxis()->SetTitle("dE/dx Lay4");
224	c2->cd(3);
225	hdedx5->Draw();
226	hdedx5->GetXaxis()->SetTitle("dE/dx Lay5");
227	c2->cd(4);
228	hdedx6->Draw();
229	hdedx6->GetXaxis()->SetTitle("dE/dx Lay6");
230
231	hdEdxVsMod->SetStats(0);
232	TCanvas* cdedx=new TCanvas("cdedx","dedx SDD",1400,600);
233	cdedx->SetLogz();
234	hdEdxVsMod->Draw("col");
235	hdEdxVsMod->GetXaxis()->SetTitle("SDD Module Id");
236	hdEdxVsMod->GetYaxis()->SetTitle("dE/dx (keV/300 #mum)");
237	hdEdxVsMod->GetYaxis()->SetTitleOffset(1.25);
238
239
240
241	TCanvas* cv=new TCanvas("cv","Vertex",600,900);
242	cv->Divide(1,3);
243	cv->cd(1);
244	hvx->Draw();
245	hvx->GetXaxis()->SetTitle("Xv (cm)");
246	cv->cd(2);
247	hvy->Draw();
248	hvy->GetXaxis()->SetTitle("Yv (cm)");
249	cv->cd(3);
250	hvz->Draw();
251	hvz->GetXaxis()->SetTitle("Xv (cm)");
252
253	hGoodPMod->SetStats(0);
254	hGoodPMod->SetTitle("");
255	TCanvas* ceff0=new TCanvas("ceff0","ModuleIndexInfo",1000,600);
256	hGoodPMod->Draw("e");
257	hGoodPMod->GetXaxis()->SetTitle("SDD Module Id");
258	hGoodPMod->GetYaxis()->SetTitle("Number of tracks");
259	hMissPMod->SetLineColor(2);
260	hMissPMod->SetMarkerColor(2);
261	hMissPMod->SetMarkerStyle(22);
262	hMissPMod->SetMarkerSize(0.5);
263	hMissPMod->Draw("psame");
264	hBadRegMod->SetLineColor(kGreen+1);
265	hBadRegMod->SetMarkerColor(kGreen+1);
266	hBadRegMod->SetMarkerStyle(20);
267	hBadRegMod->SetMarkerSize(0.5);
268	hBadRegMod->Draw("esame");
269	hSkippedMod->SetLineColor(kYellow);
270	hSkippedMod->Draw("esame");
271	hOutAccMod->SetLineColor(4);
272	hOutAccMod->Draw("esame");
273	hNoRefitMod->SetLineColor(6);
274	hNoRefitMod->Draw("esame");
275	TLatex* t1=new TLatex(0.7,0.85,"Good Point");
276	t1->SetNDC();
277	t1->SetTextColor(1);
278	t1->Draw();
279	TLatex* t2=new TLatex(0.7,0.8,"Missing Point");
280	t2->SetNDC();
281	t2->SetTextColor(2);
282	t2->Draw();
283	TLatex* t3=new TLatex(0.7,0.75,"Bad Region");
284	t3->SetNDC();
285	t3->SetTextColor(kGreen+1);
286	t3->Draw();
287	ceff0->Update();
288
289	TH1F* heff=new TH1F("heff","",260,239.5,499.5);
290	for(Int_t imod=0; imod<260;imod++){
291	Float_t numer=hGoodPMod->GetBinContent(imod+1)+hBadRegMod->GetBinContent(imod+1)+hOutAccMod->GetBinContent(imod+1)+hNoRefitMod->GetBinContent(imod+1);
292	Float_t denom=hAllPMod->GetBinContent(imod+1);
293	Float_t eff=0.;
294	Float_t erreff=0.;
295	if(denom>0){
296	eff=numer/denom;
297	erreff=TMath::Sqrt(eff*(1-eff)/denom);
298	}
299	heff->SetBinContent(imod+1,eff);
300	heff->SetBinError(imod+1,erreff);
301	}
302
303	printf("---- Modules with efficiency < 90%% ----\n");
304	heff->SetStats(0);
305	TCanvas* ceff1=new TCanvas("ceff1","Efficiency",1000,600);
306	heff->Draw();
307	heff->GetXaxis()->SetTitle("SDD Module Id");
308	heff->GetYaxis()->SetTitle("Fraction of tracks with point in good region");
309	for(Int_t ibin=1; ibin<=heff->GetNbinsX(); ibin++){
310	Float_t e=heff->GetBinContent(ibin);
311	if(e<0.9){
312	Int_t iMod=(Int_t)heff->GetBinCenter(ibin);
313	Int_t lay,lad,det;
314	AliITSgeomTGeo::GetModuleId(iMod,lay,lad,det);
315	printf("Module %d - Layer %d Ladder %2d Det %d - Eff. %.3f\n",iMod,lay,lad,det,heff->GetBinContent(ibin));
316	}
317	}
318	ceff1->Update();
319
320
321
322	hGoodPXloc->SetTitle("");
323	hGoodPZloc->SetTitle("");
324	hGoodPXloc->SetStats(0);
325	hGoodPZloc->SetStats(0);
326	hGoodPXloc->SetMinimum(0);
327	hGoodPZloc->SetMinimum(0);
328	TCanvas* ceff2=new TCanvas("ceff2","LocalCoord",1000,600);
329	ceff2->Divide(2,1);
330	ceff2->cd(1);
331	hGoodPXloc->Draw("e");
332	hGoodPXloc->GetXaxis()->SetTitle("Xlocal (cm)");
333	hGoodPXloc->GetYaxis()->SetTitle("Number of tracks");
334	hMissPXloc->SetLineColor(2);
335	hMissPXloc->SetMarkerColor(2);
336	hMissPXloc->SetMarkerStyle(22);
337	hMissPXloc->SetMarkerSize(0.5);
338	hMissPXloc->Draw("psame");
339	hBadRegXloc->SetLineColor(kGreen+1);
340	hBadRegXloc->SetMarkerColor(kGreen+1);
341	hBadRegXloc->SetMarkerStyle(20);
342	hBadRegXloc->SetMarkerSize(0.5);
343	hBadRegXloc->Draw("psame");
344	t1->Draw();
345	t2->Draw();
346	t3->Draw();
347	ceff2->cd(2);
348	hGoodPZloc->Draw("e");
349	hGoodPZloc->GetXaxis()->SetTitle("Zlocal (cm)");
350	hGoodPZloc->GetYaxis()->SetTitle("Number of tracks");
351	hMissPZloc->SetLineColor(2);
352	hMissPZloc->SetMarkerColor(2);
353	hMissPZloc->SetMarkerStyle(22);
354	hMissPZloc->SetMarkerSize(0.5);
355	hMissPZloc->Draw("psame");
356	hBadRegZloc->SetLineColor(kGreen+1);
357	hBadRegZloc->SetMarkerColor(kGreen+1);
358	hBadRegZloc->SetMarkerStyle(20);
359	hBadRegZloc->SetMarkerSize(0.5);
360	hBadRegZloc->Draw("psame");
361	t1->Draw();
362	t2->Draw();
363	t3->Draw();
364	ceff2->Update();
365
366
367
368	}
369
370
371
372
373