[u/mrichter/AliRoot.git] / ITS / macrosSDD / TrendQAtrainSDD.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 <TGrid.h>
#include <TGridResult.h>
#include <TMath.h>
#include <TSystem.h>
#include <TNtuple.h>
#include <TCanvas.h>
#include <TStyle.h>
#include <TLatex.h>
#include <TLegend.h>
#include <TLegendEntry.h>
#include "AliITSgeomTGeo.h"
#endif

Double_t LangausFun(Double_t *x, Double_t *par);
void MakePlots(TString ntupleFileName);


void TrendQAtrainSDD(TString period,
		     TString recoPass,
		     TString qaTrain1,
		     TString qaTrain2,
		     Int_t firstRun=0,
		     Int_t lastRun=999999999,
		     TString fileName="QAresults.root"){

  gStyle->SetOptStat(0);
  

  TGrid::Connect("alien:");
  Int_t year=0;
  if(period.Contains("LHC09")) year=2009;
  else if(period.Contains("LHC10")) year=2010;
  else if(period.Contains("LHC11")) year=2011;

  TString outFilNam=Form("TrendingSDD_%s_%s_%s.root",period.Data(),recoPass.Data(),qaTrain1.Data());


  const Int_t nVariables=27;
  TNtuple* ntsdd=new TNtuple("ntsdd","SDD trending","nrun:fracTrackWithClu1:errfracTrackWithClu1:fracTrackWithClu2:errfracTrackWithClu2:fracTrackWithClu3:errfracTrackWithClu3:fracTrackWithClu4:errfracTrackWithClu4:fracTrackWithClu5:errfracTrackWithClu5:fracTrackWithClu6:errfracTrackWithClu6:meanTrPts3:errmeanTrPts3:meanTrPts4:errmeanTrPts4:minDrTime:errminDrTime:meanDrTime:errmeanDrTime:fracExtra:errfracExtra:meandEdxTB0:errmeandEdxTB0:meandEdxTB5:errmeandEdxTB5");
  Float_t xnt[nVariables];

  TBits* readRun=new TBits(999999);
  readRun->ResetAllBits();
  if(!gSystem->Exec(Form("ls -l %s > /dev/null 2>&1",outFilNam.Data()))){
    TFile* oldfil=new TFile(outFilNam.Data());
    TNtuple* ntmp=(TNtuple*)oldfil->Get("ntsdd");
    Bool_t isOK=kFALSE;
    if(ntmp){
      if(ntmp->GetNvar()==ntsdd->GetNvar()){
	isOK=kTRUE;
	TObjArray* arr1=(TObjArray*)ntsdd->GetListOfBranches();
	TObjArray* arr2=(TObjArray*)ntmp->GetListOfBranches();
	for(Int_t iV=0; iV<ntmp->GetNvar(); iV++){
	  TString vnam1=arr1->At(iV)->GetName();
	  TString vnam2=arr2->At(iV)->GetName();
	  if(vnam1!=vnam2) isOK=kFALSE;
	  ntmp->SetBranchAddress(vnam2.Data(),&xnt[iV]);
	}
	if(isOK){
	  for(Int_t nE=0; nE<ntmp->GetEntries(); nE++){
	    ntmp->GetEvent(nE);
	    Int_t theRun=(Int_t)(xnt[0]+0.0001);
	    readRun->SetBitNumber(theRun);
	    ntsdd->Fill(xnt);
	  }
	}
      }
    }
    if(!isOK){
      printf("Ntuple in local file not OK -> will be recreated\n");
    }
    oldfil->Close();
    delete oldfil;
  }

  if(!gGrid||!gGrid->IsConnected()) {
    printf("gGrid not found! exit macro\n");
    return;
  }

  TString  path=Form("/alice/data/%d/%s/",year,period.Data(),recoPass.Data());
  TGridResult *gr = gGrid->Query(path,fileName);
  Int_t nFiles = gr->GetEntries();
  printf("================>%d files found\n", nFiles);
  if (nFiles < 1) return;


  if (nFiles > 1){
    for (Int_t iFil = 0; iFil <nFiles ; iFil++) { 
      TString fileNameLong=Form("%s",gr->GetKey(iFil,"turl"));
      if(!fileNameLong.Contains(recoPass.Data())) continue;
      if(!fileNameLong.Contains(qaTrain1.Data()) &&
	 !fileNameLong.Contains(qaTrain2.Data())) continue;
      if(!fileNameLong.Contains(Form("%s/%s",qaTrain1.Data(),fileName.Data())) &&
	 !fileNameLong.Contains(Form("%s/%s",qaTrain2.Data(),fileName.Data()))) continue;
      TString runNumber=fileNameLong;
      runNumber.ReplaceAll(Form("alien:///alice/data/%d/%s/",year,period.Data()),"");
      runNumber.Remove(9,runNumber.Sizeof());
   
      Int_t iRun=atoi(runNumber.Data());
      if(readRun->TestBitNumber(iRun)){ 
	printf("Run %d aleady in local ntuple -> skipping it\n",iRun);
	continue;
      }
      if(iRun<firstRun) continue;
      if(iRun>lastRun) continue;

      TFile* f=TFile::Open(fileNameLong.Data());  

      TDirectoryFile* df=(TDirectoryFile*)f->Get("SDD_Performance");
      if(!df){
	printf("Run %d SDD_Performance MISSING -> Exit\n",iRun);
	continue;
      }
      TList* l=(TList*)df->Get("coutputRP");
      if(!df){
	printf("Run %d coutputRP TList MISSING -> Exit\n",iRun);
	continue;
      }  
      TH1F* hcllay=(TH1F*)l->FindObject("hCluInLay");
      Float_t fracT[6]={0.,0.,0.,0.,0.,0.};
      Float_t efracT[6]={0.,0.,0.,0.,0.,0.};
      if(hcllay->GetBinContent(1)>0){
	for(Int_t iLay=0; iLay<6; iLay++){
	  fracT[iLay]=hcllay->GetBinContent(iLay+2)/hcllay->GetBinContent(1);
	  efracT[iLay]=TMath::Sqrt(fracT[iLay]*(1-fracT[iLay])/hcllay->GetBinContent(1));
	}
      }
      TH1F* hmodT=(TH1F*)l->FindObject("hTPMod");
      TH1F* hgamod=(TH1F*)l->FindObject("hGAMod");
      Int_t bestMod=0;
      for(Int_t iMod=0; iMod<260;iMod++){
	Int_t gda=(Int_t)hgamod->GetBinContent(iMod+1);
	if(gda>bestMod) bestMod=gda;
      }
      Int_t nChunks=1;
      if(bestMod>512){
	nChunks=(Int_t)(bestMod/512.+0.5);
      }
      hgamod->Scale(1./nChunks);

      TH1F* hev=(TH1F*)l->FindObject("hNEvents");
      Int_t nTotEvents=hev->GetBinContent(2);
      Int_t nTrigEvents=hev->GetBinContent(3);
      Int_t nEvents=nTotEvents;
      printf("Run %d Number of Events = %d Triggered=%d\n",iRun,nTotEvents,nTrigEvents);
      if(nTrigEvents>0){ 
	nEvents=nTrigEvents;
      }

      Int_t nModGood3=0;
      Int_t nModGood4=0;
      Int_t nModBadAn=0;
      Float_t sumtp3=0;
      Float_t sumtp4=0;
      Float_t sumEtp3=0;
      Float_t sumEtp4=0;
      for(Int_t iMod=0; iMod<260; iMod++){
	Float_t tps=hmodT->GetBinContent(iMod+1);
	Float_t ga=hgamod->GetBinContent(iMod+1);
	if(ga<500) nModBadAn++;
	Float_t tpsN=0.;
	Float_t etpsN=0.;
	if(ga>0){
	  tpsN=tps/ga/(Float_t)nEvents;
	  etpsN=TMath::Sqrt(tps)/ga/(Float_t)nEvents;
	  if(iMod<84){
	    sumtp3+=tpsN;
	    sumEtp3+=(etpsN*etpsN);
	    nModGood3++;
	  }
	  else{
	    sumtp4+=tpsN;
	    sumEtp4+=(etpsN*etpsN);
	    nModGood4++;
	  }
	}
      }

      TH1F* htimT=(TH1F*)l->FindObject("hDrTimTPAll");
      TH1F* htimTe=(TH1F*)l->FindObject("hDrTimTPExtra");
      
      Double_t fracExtra=0.;
      Double_t errFracExtra=0.;
      if(htimT->GetEntries()>0){
	fracExtra=htimTe->GetEntries()/htimT->GetEntries();
	errFracExtra=TMath::Sqrt(htimTe->GetEntries())/htimT->GetEntries();
      }
      Double_t averPoints=0.;
      Double_t cntBins=0.;
      for(Int_t iBin=1; iBin<=htimT->GetNbinsX(); iBin++){
	Float_t tim=htimT->GetBinCenter(iBin);
	if(tim>2000. && tim<4000.){
	  averPoints+=htimT->GetBinContent(iBin);
	  cntBins+=1;
	}
      }
      Double_t minTime=-999.;
      Double_t errMinTime=0.;
      if(cntBins>0){ 
	averPoints/=cntBins;
	for(Int_t iBin=1; iBin<=htimT->GetNbinsX(); iBin++){
	  if(htimT->GetBinContent(iBin)>0.5*averPoints){
	    minTime=htimT->GetBinCenter(iBin);
	    errMinTime=0.5*htimT->GetBinWidth(iBin);
	    break;
	  }
	}
      }
      TH1F* hSigTim0=(TH1F*)l->FindObject("hSigTimeInt0");
      TH1F* hSigTim5=(TH1F*)l->FindObject("hSigTimeInt5");

      Int_t index=0;
      xnt[index++]=(Float_t)iRun;
      xnt[index++]=fracT[0];
      xnt[index++]=efracT[0];
      xnt[index++]=fracT[1];
      xnt[index++]=efracT[1];
      xnt[index++]=fracT[2];
      xnt[index++]=efracT[2];
      xnt[index++]=fracT[3];
      xnt[index++]=efracT[3];
      xnt[index++]=fracT[4];
      xnt[index++]=efracT[4];
      xnt[index++]=fracT[5];
      xnt[index++]=efracT[5];
      xnt[index++]=sumtp3/nModGood3;
      xnt[index++]=TMath::Sqrt(sumEtp3)/nModGood3;
      xnt[index++]=sumtp4/nModGood4;
      xnt[index++]=TMath::Sqrt(sumEtp4)/nModGood4;
      xnt[index++]=minTime;
      xnt[index++]=errMinTime;
      xnt[index++]=htimT->GetMean();
      xnt[index++]=htimT->GetMeanError();
      xnt[index++]=fracExtra;
      xnt[index++]=errFracExtra;
      xnt[index++]=hSigTim0->GetMean();
      xnt[index++]=hSigTim0->GetMeanError();
      xnt[index++]=hSigTim5->GetMean();
      xnt[index++]=hSigTim5->GetMeanError();
      
      ntsdd->Fill(xnt);
    }
  }

  TFile* outfil=new TFile(outFilNam.Data(),"recreate");
  outfil->cd();
  ntsdd->Write();
  outfil->Close();

  MakePlots(outFilNam);

}

void MakePlots(TString ntupleFileName){
  TFile* fil=new TFile(ntupleFileName.Data(),"read");
  if(!fil){
    printf("File with ntuple does not exist\n");
    return;
  }
  TNtuple* ntsdd=(TNtuple*)fil->Get("ntsdd");

  Float_t nrun;
  Float_t meanTrPts3,errmeanTrPts3,meanTrPts4,errmeanTrPts4;
  Float_t minDrTime,errminDrTime,meanDrTime,errmeanDrTime;
  Float_t fracTrackWithClu1,fracTrackWithClu2,errfracTrackWithClu1,errfracTrackWithClu2;
  Float_t fracTrackWithClu3,fracTrackWithClu4,errfracTrackWithClu3,errfracTrackWithClu4;
  Float_t fracTrackWithClu5,fracTrackWithClu6,errfracTrackWithClu5,errfracTrackWithClu6;
  Float_t fracExtra,errfracExtra;
  Float_t meandEdxTB0,errmeandEdxTB0,meandEdxTB5,errmeandEdxTB5;

  ntsdd->SetBranchAddress("nrun",&nrun);
  ntsdd->SetBranchAddress("fracTrackWithClu1",&fracTrackWithClu1);
  ntsdd->SetBranchAddress("errfracTrackWithClu1",&errfracTrackWithClu1);
  ntsdd->SetBranchAddress("fracTrackWithClu2",&fracTrackWithClu2);
  ntsdd->SetBranchAddress("errfracTrackWithClu2",&errfracTrackWithClu2);
  ntsdd->SetBranchAddress("fracTrackWithClu3",&fracTrackWithClu3);
  ntsdd->SetBranchAddress("errfracTrackWithClu3",&errfracTrackWithClu3);
  ntsdd->SetBranchAddress("fracTrackWithClu4",&fracTrackWithClu4);
  ntsdd->SetBranchAddress("errfracTrackWithClu4",&errfracTrackWithClu4);
  ntsdd->SetBranchAddress("fracTrackWithClu5",&fracTrackWithClu5);
  ntsdd->SetBranchAddress("errfracTrackWithClu5",&errfracTrackWithClu5);
  ntsdd->SetBranchAddress("fracTrackWithClu6",&fracTrackWithClu6);
  ntsdd->SetBranchAddress("errfracTrackWithClu6",&errfracTrackWithClu6);

  ntsdd->SetBranchAddress("meanTrPts3",&meanTrPts3);
  ntsdd->SetBranchAddress("errmeanTrPts3",&errmeanTrPts3);
  ntsdd->SetBranchAddress("meanTrPts4",&meanTrPts4);
  ntsdd->SetBranchAddress("errmeanTrPts4",&errmeanTrPts4);
  ntsdd->SetBranchAddress("minDrTime",&minDrTime);
  ntsdd->SetBranchAddress("errminDrTime",&errminDrTime);
  ntsdd->SetBranchAddress("meanDrTime",&meanDrTime);
  ntsdd->SetBranchAddress("errmeanDrTime",&errmeanDrTime);
  ntsdd->SetBranchAddress("fracExtra",&fracExtra);
  ntsdd->SetBranchAddress("errfracExtra",&errfracExtra);
  ntsdd->SetBranchAddress("meandEdxTB0",&meandEdxTB0);
  ntsdd->SetBranchAddress("errmeandEdxTB0",&errmeandEdxTB0);
  ntsdd->SetBranchAddress("meandEdxTB5",&meandEdxTB5);
  ntsdd->SetBranchAddress("errmeandEdxTB5",&errmeandEdxTB5);

  TH1F* histotrp3=new TH1F("histotrp3","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histotrp4=new TH1F("histotrp4","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histominTime=new TH1F("histominTime","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histomeanTime=new TH1F("histomeanTime","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histofracExtra=new TH1F("histofracExtra","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histodEdxTB0=new TH1F("histodEdxTB0","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histodEdxTB5=new TH1F("histodEdxTB5","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histoTrackClu1=new TH1F("histoTrackClu1","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histoTrackClu2=new TH1F("histoTrackClu2","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histoTrackClu3=new TH1F("histoTrackClu3","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histoTrackClu4=new TH1F("histoTrackClu4","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histoTrackClu5=new TH1F("histoTrackClu5","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
  TH1F* histoTrackClu6=new TH1F("histoTrackClu6","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());


  for(Int_t i=0; i<ntsdd->GetEntries();i++){
    ntsdd->GetEvent(i);
    histoTrackClu1->SetBinContent(i+1,fracTrackWithClu1);
    histoTrackClu1->SetBinError(i+1,errfracTrackWithClu1);
    histoTrackClu1->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histoTrackClu2->SetBinContent(i+1,fracTrackWithClu2);
    histoTrackClu2->SetBinError(i+1,errfracTrackWithClu2);
    histoTrackClu2->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histoTrackClu3->SetBinContent(i+1,fracTrackWithClu3);
    histoTrackClu3->SetBinError(i+1,errfracTrackWithClu3);
    histoTrackClu3->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histoTrackClu4->SetBinContent(i+1,fracTrackWithClu4);
    histoTrackClu4->SetBinError(i+1,errfracTrackWithClu4);
    histoTrackClu4->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histoTrackClu5->SetBinContent(i+1,fracTrackWithClu5);
    histoTrackClu5->SetBinError(i+1,errfracTrackWithClu5);
    histoTrackClu5->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histoTrackClu6->SetBinContent(i+1,fracTrackWithClu6);
    histoTrackClu6->SetBinError(i+1,errfracTrackWithClu6);
    histoTrackClu6->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histominTime->SetBinContent(i+1,minDrTime);
    histominTime->SetBinError(i+1,errminDrTime);
    histominTime->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histomeanTime->SetBinContent(i+1,meanDrTime);
    histomeanTime->SetBinError(i+1,errmeanDrTime);
    histomeanTime->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histotrp3->SetBinContent(i+1,meanTrPts3);
    histotrp3->SetBinError(i+1,errmeanTrPts3);
    histotrp3->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histotrp4->SetBinContent(i+1,meanTrPts4);
    histotrp4->SetBinError(i+1,errmeanTrPts3);
    histotrp4->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histofracExtra->SetBinContent(i+1,fracExtra);
    histofracExtra->SetBinError(i+1,errfracExtra);
    histofracExtra->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histodEdxTB0->SetBinContent(i+1,meandEdxTB0);
    histodEdxTB0->SetBinError(i+1,errmeandEdxTB0);
    histodEdxTB0->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
    histodEdxTB5->SetBinContent(i+1,meandEdxTB5);
    histodEdxTB5->SetBinError(i+1,errmeandEdxTB5);
    histodEdxTB5->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
  }

  gStyle->SetOptStat(0);


  TCanvas* c1=new TCanvas("c1","Occupancy");
  histotrp3->SetLineColor(1);
  histotrp3->SetMarkerStyle(20);
  histotrp3->Draw();
  histotrp3->SetMinimum(0.);
  histotrp4->SetLineColor(2);
  histotrp4->SetMarkerColor(2);
  histotrp4->SetMarkerStyle(22);
  histotrp4->Draw("same");
  histotrp3->GetYaxis()->SetTitle("Track Point Occupancy");
  TLegend* leg=new TLegend(0.7,0.15,0.88,0.35);
  TLegendEntry* ent=leg->AddEntry(histotrp3,"Layer3","PL");
  ent=leg->AddEntry(histotrp4,"Layer4","PL");
  ent->SetTextColor(histotrp4->GetMarkerColor());
  leg->SetFillStyle(0);
  leg->Draw();
  c1->Update();

  TCanvas* c2=new TCanvas("c2","DriftTime",600,900);
  c2->Divide(1,2);
  c2->cd(1);
  histominTime->Draw();
  histominTime->SetMinimum(250);
  histominTime->SetMaximum(1000);
  histominTime->GetYaxis()->SetTitle("Minimum Drift Time (ns)");
  c2->cd(2);
  histomeanTime->Draw();
  histomeanTime->GetYaxis()->SetTitle("Average Drift Time (ns)");
  c2->Update();

  TCanvas* c3=new TCanvas("c3","ExtraClusters");
  histofracExtra->Draw();
  histofracExtra->GetYaxis()->SetTitle("Fraction of Extra Clusters");
  c3->Update();


  TCanvas* c4=new TCanvas("c4","Charge");
  histodEdxTB0->SetLineColor(1);
  histodEdxTB0->SetMarkerStyle(20);
  histodEdxTB0->Draw();
  histodEdxTB0->SetMinimum(0.);
  histodEdxTB5->SetLineColor(4);
  histodEdxTB5->SetMarkerColor(4);
  histodEdxTB5->SetMarkerStyle(23);
  histodEdxTB5->Draw("same");
  histodEdxTB0->GetYaxis()->SetTitle("<dE/dx> (keV/300 #mum)");
  TLegend* leg2=new TLegend(0.6,0.15,0.88,0.35);
  ent=leg2->AddEntry(histodEdxTB0,"Small drift time","PL");
  ent=leg2->AddEntry(histodEdxTB5,"Large drift time","PL");
  ent->SetTextColor(histodEdxTB5->GetMarkerColor());
  leg2->SetFillStyle(0);
  leg2->Draw();
  c4->Update();

  TCanvas* c5=new TCanvas("c5","TrackWithSDD");
  histoTrackClu3->Draw();
  histoTrackClu3->SetLineColor(1);
  histoTrackClu3->SetMarkerStyle(20);
  histoTrackClu3->Draw();
  histoTrackClu3->SetMinimum(0.);
  histoTrackClu4->SetLineColor(2);
  histoTrackClu4->SetMarkerColor(2);
  histoTrackClu4->SetMarkerStyle(22);
  histoTrackClu4->Draw("same");
  histoTrackClu1->SetLineColor(kGray+1);
  histoTrackClu1->SetMarkerColor(kGray+1);
  histoTrackClu1->SetMarkerStyle(24);
  histoTrackClu1->Draw("same");
  histoTrackClu2->SetLineColor(kGray+2);
  histoTrackClu2->SetMarkerColor(kGray+2);
  histoTrackClu2->SetMarkerStyle(26);
  histoTrackClu2->Draw("same");
  histoTrackClu5->SetLineColor(4);
  histoTrackClu5->SetMarkerColor(4);
  histoTrackClu5->SetMarkerStyle(29);
  histoTrackClu5->Draw("same");
  histoTrackClu6->SetLineColor(kBlue+1);
  histoTrackClu6->SetMarkerColor(kBlue+1);
  histoTrackClu6->SetMarkerStyle(30);
  histoTrackClu6->Draw("same");
  histoTrackClu3->GetYaxis()->SetTitle("Fraction of Tracks with Point in SDD Layers");
  TLegend* leg3=new TLegend(0.7,0.15,0.88,0.35);
  ent=leg3->AddEntry(histoTrackClu1,"Layer1","PL");
  ent->SetTextColor(histoTrackClu1->GetMarkerColor());
  ent=leg3->AddEntry(histoTrackClu2,"Layer2","PL");
  ent->SetTextColor(histoTrackClu2->GetMarkerColor());
  ent=leg3->AddEntry(histoTrackClu3,"Layer3","PL");
  ent->SetTextColor(histoTrackClu3->GetMarkerColor());
  ent=leg3->AddEntry(histoTrackClu4,"Layer4","PL");
  ent->SetTextColor(histoTrackClu4->GetMarkerColor());
  ent=leg3->AddEntry(histoTrackClu5,"Layer5","PL");
  ent->SetTextColor(histoTrackClu5->GetMarkerColor());
  ent=leg3->AddEntry(histoTrackClu6,"Layer6","PL");
  ent->SetTextColor(histoTrackClu6->GetMarkerColor());

  leg3->SetFillStyle(0);
  leg3->Draw();
  c5->Update();

}

Double_t LangausFun(Double_t *x, Double_t *par) {

  //Fit parameters:
  //par[0]=Width (scale) parameter of Landau density
  //par[1]=Most Probable (MP, location) parameter of Landau density
  //par[2]=Total area (integral -inf to inf, normalization constant)
  //par[3]=Width (sigma) of convoluted Gaussian function
  //
  //In the Landau distribution (represented by the CERNLIB approximation), 
  //the maximum is located at x=-0.22278298 with the location parameter=0.
  //This shift is corrected within this function, so that the actual
  //maximum is identical to the MP parameter.

  // Numeric constants
  Double_t invsq2pi = 0.3989422804014;   // (2 pi)^(-1/2)
  Double_t mpshift  = -0.22278298;       // Landau maximum location

  // Control constants
  Double_t np = 100.0;      // number of convolution steps
  Double_t sc =   5.0;      // convolution extends to +-sc Gaussian sigmas

  // Variables
  Double_t xx;
  Double_t mpc;
  Double_t fland;
  Double_t sum = 0.0;
  Double_t xlow,xupp;
  Double_t step;
  Double_t i;


  // MP shift correction
  mpc = par[1] - mpshift * par[0]; 

  // Range of convolution integral
  xlow = x[0] - sc * par[3];
  xupp = x[0] + sc * par[3];

  step = (xupp-xlow) / np;

  // Convolution integral of Landau and Gaussian by sum
  for(i=1.0; i<=np/2; i++) {
    xx = xlow + (i-.5) * step;
    fland = TMath::Landau(xx,mpc,par[0]) / par[0];
    sum += fland * TMath::Gaus(x[0],xx,par[3]);

    xx = xupp - (i-.5) * step;
    fland = TMath::Landau(xx,mpc,par[0]) / par[0];
    sum += fland * TMath::Gaus(x[0],xx,par[3]);
  }

  return (par[2] * step * sum * invsq2pi / par[3]);
}
Commit	Line	Data
341adf04	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 <TGrid.h>
	11	#include <TGridResult.h>
	12	#include <TMath.h>
	13	#include <TSystem.h>
	14	#include <TNtuple.h>
	15	#include <TCanvas.h>
	16	#include <TStyle.h>
	17	#include <TLatex.h>
	18	#include <TLegend.h>
	19	#include <TLegendEntry.h>
	20	#include "AliITSgeomTGeo.h"
	21	#endif
	22
	23	Double_t LangausFun(Double_t x, Double_t par);
	24	void MakePlots(TString ntupleFileName);
	25
	26
	27	void TrendQAtrainSDD(TString period,
	28	TString recoPass,
e1a5dc0a	29	TString qaTrain1,
e1a5dc0a	30	TString qaTrain2,
341adf04	31	Int_t firstRun=0,
	32	Int_t lastRun=999999999,
	33	TString fileName="QAresults.root"){
	34
	35	gStyle->SetOptStat(0);
	36
	37
	38	TGrid::Connect("alien:");
	39	Int_t year=0;
	40	if(period.Contains("LHC09")) year=2009;
	41	else if(period.Contains("LHC10")) year=2010;
4d8feb56	42	else if(period.Contains("LHC11")) year=2011;
341adf04	43
e1a5dc0a	44	TString outFilNam=Form("TrendingSDD_%s_%s_%s.root",period.Data(),recoPass.Data(),qaTrain1.Data());
341adf04	45
341adf04	46
e1a5dc0a	47	const Int_t nVariables=27;
e1a5dc0a	48	TNtuple* ntsdd=new TNtuple("ntsdd","SDD trending","nrun:fracTrackWithClu1:errfracTrackWithClu1:fracTrackWithClu2:errfracTrackWithClu2:fracTrackWithClu3:errfracTrackWithClu3:fracTrackWithClu4:errfracTrackWithClu4:fracTrackWithClu5:errfracTrackWithClu5:fracTrackWithClu6:errfracTrackWithClu6:meanTrPts3:errmeanTrPts3:meanTrPts4:errmeanTrPts4:minDrTime:errminDrTime:meanDrTime:errmeanDrTime:fracExtra:errfracExtra:meandEdxTB0:errmeandEdxTB0:meandEdxTB5:errmeandEdxTB5");
341adf04	49	Float_t xnt[nVariables];
	50
	51	TBits* readRun=new TBits(999999);
	52	readRun->ResetAllBits();
	53	if(!gSystem->Exec(Form("ls -l %s > /dev/null 2>&1",outFilNam.Data()))){
	54	TFile* oldfil=new TFile(outFilNam.Data());
	55	TNtuple* ntmp=(TNtuple*)oldfil->Get("ntsdd");
	56	Bool_t isOK=kFALSE;
	57	if(ntmp){
	58	if(ntmp->GetNvar()==ntsdd->GetNvar()){
	59	isOK=kTRUE;
	60	TObjArray* arr1=(TObjArray*)ntsdd->GetListOfBranches();
	61	TObjArray* arr2=(TObjArray*)ntmp->GetListOfBranches();
	62	for(Int_t iV=0; iV<ntmp->GetNvar(); iV++){
	63	TString vnam1=arr1->At(iV)->GetName();
	64	TString vnam2=arr2->At(iV)->GetName();
	65	if(vnam1!=vnam2) isOK=kFALSE;
	66	ntmp->SetBranchAddress(vnam2.Data(),&xnt[iV]);
	67	}
	68	if(isOK){
	69	for(Int_t nE=0; nE<ntmp->GetEntries(); nE++){
	70	ntmp->GetEvent(nE);
	71	Int_t theRun=(Int_t)(xnt[0]+0.0001);
	72	readRun->SetBitNumber(theRun);
	73	ntsdd->Fill(xnt);
	74	}
	75	}
	76	}
	77	}
	78	if(!isOK){
	79	printf("Ntuple in local file not OK -> will be recreated\n");
	80	}
	81	oldfil->Close();
	82	delete oldfil;
	83	}
	84
	85	if(!gGrid\|\|!gGrid->IsConnected()) {
	86	printf("gGrid not found! exit macro\n");
	87	return;
	88	}
	89
	90	TString path=Form("/alice/data/%d/%s/",year,period.Data(),recoPass.Data());
	91	TGridResult *gr = gGrid->Query(path,fileName);
	92	Int_t nFiles = gr->GetEntries();
	93	printf("================>%d files found\n", nFiles);
	94	if (nFiles < 1) return;
	95
	96
	97	if (nFiles > 1){
	98	for (Int_t iFil = 0; iFil <nFiles ; iFil++) {
	99	TString fileNameLong=Form("%s",gr->GetKey(iFil,"turl"));
	100	if(!fileNameLong.Contains(recoPass.Data())) continue;
e1a5dc0a	101	if(!fileNameLong.Contains(qaTrain1.Data()) &&
	102	!fileNameLong.Contains(qaTrain2.Data())) continue;
	103	if(!fileNameLong.Contains(Form("%s/%s",qaTrain1.Data(),fileName.Data())) &&
	104	!fileNameLong.Contains(Form("%s/%s",qaTrain2.Data(),fileName.Data()))) continue;
341adf04	105	TString runNumber=fileNameLong;
	106	runNumber.ReplaceAll(Form("alien:///alice/data/%d/%s/",year,period.Data()),"");
	107	runNumber.Remove(9,runNumber.Sizeof());
	108
	109	Int_t iRun=atoi(runNumber.Data());
	110	if(readRun->TestBitNumber(iRun)){
	111	printf("Run %d aleady in local ntuple -> skipping it\n",iRun);
	112	continue;
	113	}
	114	if(iRun<firstRun) continue;
	115	if(iRun>lastRun) continue;
	116
	117	TFile* f=TFile::Open(fileNameLong.Data());
	118
	119	TDirectoryFile* df=(TDirectoryFile*)f->Get("SDD_Performance");
	120	if(!df){
	121	printf("Run %d SDD_Performance MISSING -> Exit\n",iRun);
	122	continue;
	123	}
	124	TList* l=(TList*)df->Get("coutputRP");
	125	if(!df){
	126	printf("Run %d coutputRP TList MISSING -> Exit\n",iRun);
	127	continue;
e1a5dc0a	128	}
	129	TH1F* hcllay=(TH1F*)l->FindObject("hCluInLay");
	130	Float_t fracT[6]={0.,0.,0.,0.,0.,0.};
	131	Float_t efracT[6]={0.,0.,0.,0.,0.,0.};
	132	if(hcllay->GetBinContent(1)>0){
	133	for(Int_t iLay=0; iLay<6; iLay++){
	134	fracT[iLay]=hcllay->GetBinContent(iLay+2)/hcllay->GetBinContent(1);
	135	efracT[iLay]=TMath::Sqrt(fracT[iLay]*(1-fracT[iLay])/hcllay->GetBinContent(1));
	136	}
341adf04	137	}
	138	TH1F* hmodT=(TH1F*)l->FindObject("hTPMod");
	139	TH1F* hgamod=(TH1F*)l->FindObject("hGAMod");
	140	Int_t bestMod=0;
	141	for(Int_t iMod=0; iMod<260;iMod++){
	142	Int_t gda=(Int_t)hgamod->GetBinContent(iMod+1);
	143	if(gda>bestMod) bestMod=gda;
	144	}
	145	Int_t nChunks=1;
	146	if(bestMod>512){
	147	nChunks=(Int_t)(bestMod/512.+0.5);
	148	}
	149	hgamod->Scale(1./nChunks);
	150
	151	TH1F* hev=(TH1F*)l->FindObject("hNEvents");
	152	Int_t nTotEvents=hev->GetBinContent(2);
	153	Int_t nTrigEvents=hev->GetBinContent(3);
	154	Int_t nEvents=nTotEvents;
	155	printf("Run %d Number of Events = %d Triggered=%d\n",iRun,nTotEvents,nTrigEvents);
	156	if(nTrigEvents>0){
	157	nEvents=nTrigEvents;
	158	}
	159
	160	Int_t nModGood3=0;
	161	Int_t nModGood4=0;
	162	Int_t nModBadAn=0;
	163	Float_t sumtp3=0;
	164	Float_t sumtp4=0;
	165	Float_t sumEtp3=0;
	166	Float_t sumEtp4=0;
	167	for(Int_t iMod=0; iMod<260; iMod++){
	168	Float_t tps=hmodT->GetBinContent(iMod+1);
	169	Float_t ga=hgamod->GetBinContent(iMod+1);
	170	if(ga<500) nModBadAn++;
	171	Float_t tpsN=0.;
	172	Float_t etpsN=0.;
	173	if(ga>0){
	174	tpsN=tps/ga/(Float_t)nEvents;
	175	etpsN=TMath::Sqrt(tps)/ga/(Float_t)nEvents;
	176	if(iMod<84){
	177	sumtp3+=tpsN;
	178	sumEtp3+=(etpsN*etpsN);
	179	nModGood3++;
	180	}
	181	else{
	182	sumtp4+=tpsN;
	183	sumEtp4+=(etpsN*etpsN);
	184	nModGood4++;
	185	}
	186	}
	187	}
	188
	189	TH1F* htimT=(TH1F*)l->FindObject("hDrTimTPAll");
	190	TH1F* htimTe=(TH1F*)l->FindObject("hDrTimTPExtra");
	191
	192	Double_t fracExtra=0.;
	193	Double_t errFracExtra=0.;
	194	if(htimT->GetEntries()>0){
	195	fracExtra=htimTe->GetEntries()/htimT->GetEntries();
	196	errFracExtra=TMath::Sqrt(htimTe->GetEntries())/htimT->GetEntries();
	197	}
	198	Double_t averPoints=0.;
	199	Double_t cntBins=0.;
	200	for(Int_t iBin=1; iBin<=htimT->GetNbinsX(); iBin++){
201	Float_t tim=htimT->GetBinCenter(iBin);
202	if(tim>2000. && tim<4000.){
203	averPoints+=htimT->GetBinContent(iBin);
204	cntBins+=1;
205	}
206	}
207	Double_t minTime=-999.;
208	Double_t errMinTime=0.;
209	if(cntBins>0){
210	averPoints/=cntBins;
211	for(Int_t iBin=1; iBin<=htimT->GetNbinsX(); iBin++){
212	if(htimT->GetBinContent(iBin)>0.5*averPoints){
213	minTime=htimT->GetBinCenter(iBin);
214	errMinTime=0.5*htimT->GetBinWidth(iBin);
215	break;
216	}
217	}
218	}
219	TH1F* hSigTim0=(TH1F*)l->FindObject("hSigTimeInt0");
220	TH1F* hSigTim5=(TH1F*)l->FindObject("hSigTimeInt5");
221
e1a5dc0a	222	Int_t index=0;
	223	xnt[index++]=(Float_t)iRun;
	224	xnt[index++]=fracT[0];
	225	xnt[index++]=efracT[0];
	226	xnt[index++]=fracT[1];
	227	xnt[index++]=efracT[1];
	228	xnt[index++]=fracT[2];
	229	xnt[index++]=efracT[2];
	230	xnt[index++]=fracT[3];
	231	xnt[index++]=efracT[3];
	232	xnt[index++]=fracT[4];
	233	xnt[index++]=efracT[4];
	234	xnt[index++]=fracT[5];
	235	xnt[index++]=efracT[5];
	236	xnt[index++]=sumtp3/nModGood3;
	237	xnt[index++]=TMath::Sqrt(sumEtp3)/nModGood3;
	238	xnt[index++]=sumtp4/nModGood4;
	239	xnt[index++]=TMath::Sqrt(sumEtp4)/nModGood4;
	240	xnt[index++]=minTime;
	241	xnt[index++]=errMinTime;
	242	xnt[index++]=htimT->GetMean();
	243	xnt[index++]=htimT->GetMeanError();
	244	xnt[index++]=fracExtra;
	245	xnt[index++]=errFracExtra;
	246	xnt[index++]=hSigTim0->GetMean();
	247	xnt[index++]=hSigTim0->GetMeanError();
	248	xnt[index++]=hSigTim5->GetMean();
	249	xnt[index++]=hSigTim5->GetMeanError();
341adf04	250
	251	ntsdd->Fill(xnt);
	252	}
	253	}
	254
	255	TFile* outfil=new TFile(outFilNam.Data(),"recreate");
	256	outfil->cd();
	257	ntsdd->Write();
	258	outfil->Close();
	259
	260	MakePlots(outFilNam);
	261
	262	}
	263
	264	void MakePlots(TString ntupleFileName){
	265	TFile* fil=new TFile(ntupleFileName.Data(),"read");
	266	if(!fil){
	267	printf("File with ntuple does not exist\n");
	268	return;
	269	}
	270	TNtuple* ntsdd=(TNtuple*)fil->Get("ntsdd");
	271
	272	Float_t nrun;
	273	Float_t meanTrPts3,errmeanTrPts3,meanTrPts4,errmeanTrPts4;
	274	Float_t minDrTime,errminDrTime,meanDrTime,errmeanDrTime;
e1a5dc0a	275	Float_t fracTrackWithClu1,fracTrackWithClu2,errfracTrackWithClu1,errfracTrackWithClu2;
	276	Float_t fracTrackWithClu3,fracTrackWithClu4,errfracTrackWithClu3,errfracTrackWithClu4;
	277	Float_t fracTrackWithClu5,fracTrackWithClu6,errfracTrackWithClu5,errfracTrackWithClu6;
341adf04	278	Float_t fracExtra,errfracExtra;
	279	Float_t meandEdxTB0,errmeandEdxTB0,meandEdxTB5,errmeandEdxTB5;
	280
	281	ntsdd->SetBranchAddress("nrun",&nrun);
e1a5dc0a	282	ntsdd->SetBranchAddress("fracTrackWithClu1",&fracTrackWithClu1);
	283	ntsdd->SetBranchAddress("errfracTrackWithClu1",&errfracTrackWithClu1);
	284	ntsdd->SetBranchAddress("fracTrackWithClu2",&fracTrackWithClu2);
	285	ntsdd->SetBranchAddress("errfracTrackWithClu2",&errfracTrackWithClu2);
	286	ntsdd->SetBranchAddress("fracTrackWithClu3",&fracTrackWithClu3);
	287	ntsdd->SetBranchAddress("errfracTrackWithClu3",&errfracTrackWithClu3);
	288	ntsdd->SetBranchAddress("fracTrackWithClu4",&fracTrackWithClu4);
	289	ntsdd->SetBranchAddress("errfracTrackWithClu4",&errfracTrackWithClu4);
	290	ntsdd->SetBranchAddress("fracTrackWithClu5",&fracTrackWithClu5);
	291	ntsdd->SetBranchAddress("errfracTrackWithClu5",&errfracTrackWithClu5);
	292	ntsdd->SetBranchAddress("fracTrackWithClu6",&fracTrackWithClu6);
	293	ntsdd->SetBranchAddress("errfracTrackWithClu6",&errfracTrackWithClu6);
	294
341adf04	295	ntsdd->SetBranchAddress("meanTrPts3",&meanTrPts3);
	296	ntsdd->SetBranchAddress("errmeanTrPts3",&errmeanTrPts3);
	297	ntsdd->SetBranchAddress("meanTrPts4",&meanTrPts4);
	298	ntsdd->SetBranchAddress("errmeanTrPts4",&errmeanTrPts4);
	299	ntsdd->SetBranchAddress("minDrTime",&minDrTime);
	300	ntsdd->SetBranchAddress("errminDrTime",&errminDrTime);
	301	ntsdd->SetBranchAddress("meanDrTime",&meanDrTime);
	302	ntsdd->SetBranchAddress("errmeanDrTime",&errmeanDrTime);
	303	ntsdd->SetBranchAddress("fracExtra",&fracExtra);
	304	ntsdd->SetBranchAddress("errfracExtra",&errfracExtra);
	305	ntsdd->SetBranchAddress("meandEdxTB0",&meandEdxTB0);
	306	ntsdd->SetBranchAddress("errmeandEdxTB0",&errmeandEdxTB0);
	307	ntsdd->SetBranchAddress("meandEdxTB5",&meandEdxTB5);
	308	ntsdd->SetBranchAddress("errmeandEdxTB5",&errmeandEdxTB5);
	309
	310	TH1F* histotrp3=new TH1F("histotrp3","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	311	TH1F* histotrp4=new TH1F("histotrp4","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	312	TH1F* histominTime=new TH1F("histominTime","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	313	TH1F* histomeanTime=new TH1F("histomeanTime","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	314	TH1F* histofracExtra=new TH1F("histofracExtra","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	315	TH1F* histodEdxTB0=new TH1F("histodEdxTB0","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	316	TH1F* histodEdxTB5=new TH1F("histodEdxTB5","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
e1a5dc0a	317	TH1F* histoTrackClu1=new TH1F("histoTrackClu1","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	318	TH1F* histoTrackClu2=new TH1F("histoTrackClu2","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	319	TH1F* histoTrackClu3=new TH1F("histoTrackClu3","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	320	TH1F* histoTrackClu4=new TH1F("histoTrackClu4","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	321	TH1F* histoTrackClu5=new TH1F("histoTrackClu5","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	322	TH1F* histoTrackClu6=new TH1F("histoTrackClu6","",(Int_t)ntsdd->GetEntries(),0.,ntsdd->GetEntries());
	323
341adf04	324
	325	for(Int_t i=0; i<ntsdd->GetEntries();i++){
	326	ntsdd->GetEvent(i);
e1a5dc0a	327	histoTrackClu1->SetBinContent(i+1,fracTrackWithClu1);
	328	histoTrackClu1->SetBinError(i+1,errfracTrackWithClu1);
	329	histoTrackClu1->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	330	histoTrackClu2->SetBinContent(i+1,fracTrackWithClu2);
	331	histoTrackClu2->SetBinError(i+1,errfracTrackWithClu2);
	332	histoTrackClu2->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	333	histoTrackClu3->SetBinContent(i+1,fracTrackWithClu3);
	334	histoTrackClu3->SetBinError(i+1,errfracTrackWithClu3);
	335	histoTrackClu3->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	336	histoTrackClu4->SetBinContent(i+1,fracTrackWithClu4);
	337	histoTrackClu4->SetBinError(i+1,errfracTrackWithClu4);
	338	histoTrackClu4->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	339	histoTrackClu5->SetBinContent(i+1,fracTrackWithClu5);
	340	histoTrackClu5->SetBinError(i+1,errfracTrackWithClu5);
	341	histoTrackClu5->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	342	histoTrackClu6->SetBinContent(i+1,fracTrackWithClu6);
	343	histoTrackClu6->SetBinError(i+1,errfracTrackWithClu6);
	344	histoTrackClu6->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
341adf04	345	histominTime->SetBinContent(i+1,minDrTime);
	346	histominTime->SetBinError(i+1,errminDrTime);
	347	histominTime->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	348	histomeanTime->SetBinContent(i+1,meanDrTime);
	349	histomeanTime->SetBinError(i+1,errmeanDrTime);
	350	histomeanTime->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	351	histotrp3->SetBinContent(i+1,meanTrPts3);
	352	histotrp3->SetBinError(i+1,errmeanTrPts3);
	353	histotrp3->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	354	histotrp4->SetBinContent(i+1,meanTrPts4);
	355	histotrp4->SetBinError(i+1,errmeanTrPts3);
	356	histotrp4->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	357	histofracExtra->SetBinContent(i+1,fracExtra);
	358	histofracExtra->SetBinError(i+1,errfracExtra);
	359	histofracExtra->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	360	histodEdxTB0->SetBinContent(i+1,meandEdxTB0);
	361	histodEdxTB0->SetBinError(i+1,errmeandEdxTB0);
	362	histodEdxTB0->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	363	histodEdxTB5->SetBinContent(i+1,meandEdxTB5);
	364	histodEdxTB5->SetBinError(i+1,errmeandEdxTB5);
	365	histodEdxTB5->GetXaxis()->SetBinLabel(i+1,Form("%d",(Int_t)nrun));
	366	}
	367
	368	gStyle->SetOptStat(0);
	369
	370
	371	TCanvas* c1=new TCanvas("c1","Occupancy");
	372	histotrp3->SetLineColor(1);
	373	histotrp3->SetMarkerStyle(20);
	374	histotrp3->Draw();
	375	histotrp3->SetMinimum(0.);
	376	histotrp4->SetLineColor(2);
	377	histotrp4->SetMarkerColor(2);
	378	histotrp4->SetMarkerStyle(22);
	379	histotrp4->Draw("same");
	380	histotrp3->GetYaxis()->SetTitle("Track Point Occupancy");
	381	TLegend* leg=new TLegend(0.7,0.15,0.88,0.35);
	382	TLegendEntry* ent=leg->AddEntry(histotrp3,"Layer3","PL");
	383	ent=leg->AddEntry(histotrp4,"Layer4","PL");
	384	ent->SetTextColor(histotrp4->GetMarkerColor());
	385	leg->SetFillStyle(0);
	386	leg->Draw();
	387	c1->Update();
	388
	389	TCanvas* c2=new TCanvas("c2","DriftTime",600,900);
	390	c2->Divide(1,2);
	391	c2->cd(1);
	392	histominTime->Draw();
	393	histominTime->SetMinimum(250);
	394	histominTime->SetMaximum(1000);
	395	histominTime->GetYaxis()->SetTitle("Minimum Drift Time (ns)");
	396	c2->cd(2);
	397	histomeanTime->Draw();
	398	histomeanTime->GetYaxis()->SetTitle("Average Drift Time (ns)");
	399	c2->Update();
	400
	401	TCanvas* c3=new TCanvas("c3","ExtraClusters");
	402	histofracExtra->Draw();
	403	histofracExtra->GetYaxis()->SetTitle("Fraction of Extra Clusters");
	404	c3->Update();
	405
	406
	407	TCanvas* c4=new TCanvas("c4","Charge");
	408	histodEdxTB0->SetLineColor(1);
409	histodEdxTB0->SetMarkerStyle(20);
410	histodEdxTB0->Draw();
411	histodEdxTB0->SetMinimum(0.);
412	histodEdxTB5->SetLineColor(4);
413	histodEdxTB5->SetMarkerColor(4);
414	histodEdxTB5->SetMarkerStyle(23);
415	histodEdxTB5->Draw("same");
e1a5dc0a	416	histodEdxTB0->GetYaxis()->SetTitle("<dE/dx> (keV/300 #mum)");
341adf04	417	TLegend* leg2=new TLegend(0.6,0.15,0.88,0.35);
	418	ent=leg2->AddEntry(histodEdxTB0,"Small drift time","PL");
	419	ent=leg2->AddEntry(histodEdxTB5,"Large drift time","PL");
	420	ent->SetTextColor(histodEdxTB5->GetMarkerColor());
	421	leg2->SetFillStyle(0);
	422	leg2->Draw();
	423	c4->Update();
	424
e1a5dc0a	425	TCanvas* c5=new TCanvas("c5","TrackWithSDD");
	426	histoTrackClu3->Draw();
	427	histoTrackClu3->SetLineColor(1);
	428	histoTrackClu3->SetMarkerStyle(20);
	429	histoTrackClu3->Draw();
	430	histoTrackClu3->SetMinimum(0.);
	431	histoTrackClu4->SetLineColor(2);
	432	histoTrackClu4->SetMarkerColor(2);
	433	histoTrackClu4->SetMarkerStyle(22);
	434	histoTrackClu4->Draw("same");
	435	histoTrackClu1->SetLineColor(kGray+1);
	436	histoTrackClu1->SetMarkerColor(kGray+1);
	437	histoTrackClu1->SetMarkerStyle(24);
	438	histoTrackClu1->Draw("same");
	439	histoTrackClu2->SetLineColor(kGray+2);
	440	histoTrackClu2->SetMarkerColor(kGray+2);
	441	histoTrackClu2->SetMarkerStyle(26);
	442	histoTrackClu2->Draw("same");
	443	histoTrackClu5->SetLineColor(4);
	444	histoTrackClu5->SetMarkerColor(4);
	445	histoTrackClu5->SetMarkerStyle(29);
	446	histoTrackClu5->Draw("same");
	447	histoTrackClu6->SetLineColor(kBlue+1);
	448	histoTrackClu6->SetMarkerColor(kBlue+1);
	449	histoTrackClu6->SetMarkerStyle(30);
	450	histoTrackClu6->Draw("same");
	451	histoTrackClu3->GetYaxis()->SetTitle("Fraction of Tracks with Point in SDD Layers");
	452	TLegend* leg3=new TLegend(0.7,0.15,0.88,0.35);
	453	ent=leg3->AddEntry(histoTrackClu1,"Layer1","PL");
	454	ent->SetTextColor(histoTrackClu1->GetMarkerColor());
	455	ent=leg3->AddEntry(histoTrackClu2,"Layer2","PL");
	456	ent->SetTextColor(histoTrackClu2->GetMarkerColor());
	457	ent=leg3->AddEntry(histoTrackClu3,"Layer3","PL");
	458	ent->SetTextColor(histoTrackClu3->GetMarkerColor());
	459	ent=leg3->AddEntry(histoTrackClu4,"Layer4","PL");
	460	ent->SetTextColor(histoTrackClu4->GetMarkerColor());
	461	ent=leg3->AddEntry(histoTrackClu5,"Layer5","PL");
	462	ent->SetTextColor(histoTrackClu5->GetMarkerColor());
	463	ent=leg3->AddEntry(histoTrackClu6,"Layer6","PL");
	464	ent->SetTextColor(histoTrackClu6->GetMarkerColor());
	465
	466	leg3->SetFillStyle(0);
	467	leg3->Draw();
	468	c5->Update();
341adf04	469
	470	}
	471
	472	Double_t LangausFun(Double_t x, Double_t par) {
	473
	474	//Fit parameters:
	475	//par[0]=Width (scale) parameter of Landau density
	476	//par[1]=Most Probable (MP, location) parameter of Landau density
	477	//par[2]=Total area (integral -inf to inf, normalization constant)
	478	//par[3]=Width (sigma) of convoluted Gaussian function
	479	//
	480	//In the Landau distribution (represented by the CERNLIB approximation),
	481	//the maximum is located at x=-0.22278298 with the location parameter=0.
	482	//This shift is corrected within this function, so that the actual
	483	//maximum is identical to the MP parameter.
	484
	485	// Numeric constants
	486	Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
	487	Double_t mpshift = -0.22278298; // Landau maximum location
	488
	489	// Control constants
	490	Double_t np = 100.0; // number of convolution steps
	491	Double_t sc = 5.0; // convolution extends to +-sc Gaussian sigmas
	492
	493	// Variables
	494	Double_t xx;
	495	Double_t mpc;
	496	Double_t fland;
	497	Double_t sum = 0.0;
	498	Double_t xlow,xupp;
	499	Double_t step;
	500	Double_t i;
	501
	502
	503	// MP shift correction
	504	mpc = par[1] - mpshift * par[0];
	505
	506	// Range of convolution integral
	507	xlow = x[0] - sc * par[3];
	508	xupp = x[0] + sc * par[3];
	509
	510	step = (xupp-xlow) / np;
	511
	512	// Convolution integral of Landau and Gaussian by sum
	513	for(i=1.0; i<=np/2; i++) {
	514	xx = xlow + (i-.5) * step;
	515	fland = TMath::Landau(xx,mpc,par[0]) / par[0];
	516	sum += fland * TMath::Gaus(x[0],xx,par[3]);
	517
	518	xx = xupp - (i-.5) * step;
	519	fland = TMath::Landau(xx,mpc,par[0]) / par[0];
	520	sum += fland * TMath::Gaus(x[0],xx,par[3]);
	521	}
	522
	523	return (par[2] * step * sum * invsq2pi / par[3]);
	524	}