Fix for FMD DA
[u/mrichter/AliRoot.git] / ITS / ShowITSRecPoints.C
index 6150079..4ef7b74 100644 (file)
@@ -5,18 +5,20 @@
 #include <TGraph2D.h>
 #include <TGeoManager.h>
 #include <TH1.h>
+#include <TH2.h>
+#include <TLatex.h>
+#include <TStyle.h>
 #include <TInterpreter.h>
-#include "AliCDBManager.h"
 #include "AliGeomManager.h"
 #include "AliHeader.h"
 #include "AliITS.h"
-#include "AliITSDetTypeRec.h"
-#include "AliITSgeom.h"
+#include "AliITSRecPointContainer.h"
+#include "AliITSgeomTGeo.h"
 #include "AliITSRecPoint.h"
 #include "AliRun.h"
 #endif
 
-Int_t ShowITSRecPoints(){
+Int_t ShowITSRecPoints(Int_t nevfordisp=0){
   ///////////////////////////////////////////////////////////////////////
   // Macro to check clusters in the 6 ITS layers                       //
   // Provides:                                                         //
@@ -27,18 +29,8 @@ Int_t ShowITSRecPoints(){
   if (gClassTable->GetID("AliRun") < 0) {
     gInterpreter->ExecuteMacro("loadlibs.C");
   }
-  // Set OCDB if needed
-  AliCDBManager* man = AliCDBManager::Instance();
-  if (!man->IsDefaultStorageSet()) {
-    printf("Setting a local default storage and run number 0\n");
-    man->SetDefaultStorage("local://$ALICE_ROOT");
-    man->SetRun(0);
-  }
-  else {
-    printf("Using deafult storage \n");
-  }
  
-  // retrives geometry 
+  // retrieves geometry 
   if(!gGeoManager){
     AliGeomManager::LoadGeometry("geometry.root");
   }
@@ -48,12 +40,6 @@ Int_t ShowITSRecPoints(){
     cerr<<"Can not open session RL=NULL"<< endl;
     return -1;
   }
-  Int_t retval = rl->LoadHeader();
-  if (retval){
-    cerr<<"LoadHeader returned error"<<endl;
-    return -1;
-  }
-
 
   AliITSLoader* ITSloader =  (AliITSLoader*) rl->GetLoader("ITSLoader");
   if(!ITSloader){
@@ -61,16 +47,15 @@ Int_t ShowITSRecPoints(){
     return -1;
   }
   ITSloader->LoadRecPoints("read");
-  AliITSgeom *geom = ITSloader->GetITSgeom();
+
+  AliITSRecPointContainer* rpcont=AliITSRecPointContainer::Instance();
 
   Float_t cluglo[3]={0.,0.,0.}; 
-  AliITSDetTypeRec* detTypeRec = new AliITSDetTypeRec();
-  detTypeRec->SetITSgeom(geom);
-  detTypeRec->SetDefaults();
 
-  Int_t modmin=geom->GetStartDet(0);
-  Int_t modmax=geom->GetLastDet(2);
-  Int_t totmod=modmax-modmin;
+  Int_t totmod=AliITSgeomTGeo::GetNModules();
+  Int_t modmin=AliITSgeomTGeo::GetModuleIndex(1,1,1);
+  Int_t modmax=AliITSgeomTGeo::GetNModules()-1;
+
   Float_t xlim[6]={4.5,7.5,16.,26.,40.,45.};
   Float_t zlim[6]={15.,15.,22.,30.,45.,55.};
 
@@ -83,8 +68,10 @@ Int_t ShowITSRecPoints(){
   TH1F** hzg=new TH1F*[6];
   TH1F** hr=new TH1F*[6];
   TH1F** hphi=new TH1F*[6];
+  TH2F** hzphi=new TH2F*[6];
   TH1F** hq=new TH1F*[6];
-
+  TH1F** hdrtime=new TH1F*[6];
+  
   Char_t name[10];
   for(Int_t iLay=0;iLay<6;iLay++){
     sprintf(name,"hmod%d",iLay+1);
@@ -123,67 +110,65 @@ Int_t ShowITSRecPoints(){
     hq[iLay]=new TH1F(name,"",100,0.,300.);    
     hq[iLay]->GetXaxis()->SetTitle("Charge (keV)");
     hq[iLay]->GetXaxis()->CenterTitle();
+    sprintf(name,"hdtim%d",iLay+1);
+    hdrtime[iLay]=new TH1F(name,"",100,0.,7000.);    
+    hdrtime[iLay]->GetXaxis()->SetTitle("Drift Time (ns)");
+    hdrtime[iLay]->GetXaxis()->CenterTitle();
+    sprintf(name,"hzphi%d",iLay+1);
+    hzphi[iLay]=new TH2F(name,Form("Layer %d",iLay+1),50,-TMath::Pi(),TMath::Pi(),50,-zlim[iLay],zlim[iLay]);
+    hzphi[iLay]->GetXaxis()->SetTitle("#varphi (rad)");
+    hzphi[iLay]->GetYaxis()->SetTitle("Zglob (cm)");
+    hzphi[iLay]->SetStats(0);
   }
 
-  TGraph **gpts=new TGraph*[3];
-  TGraph2D **gpts3d=new TGraph2D*[3];
-  for(Int_t i=0;i<3;i++){
-    gpts[i]=new TGraph(0);
-    gpts3d[i]=new TGraph2D(0);
-  }
+  TGraph *gptsXY=new TGraph(0);
+  TGraph *gptsRZ=new TGraph(0);
+
   Int_t totev=rl->GetNumberOfEvents();
   printf("Total Number of events = %d\n",totev);
 
   for(Int_t iev=0;iev<totev;iev++){
     rl->GetEvent(iev);
     TTree *TR = ITSloader->TreeR();
-    TClonesArray *ITSrec  = detTypeRec->RecPoints();
-    TBranch *branch = 0;
-    if(TR && ITSrec){
-      branch = ITSloader->TreeR()->GetBranch("ITSRecPoints");
-      if(branch)branch->SetAddress(&ITSrec);
-    }
-    Int_t nparticles = rl->GetHeader()->GetNtrack();
-    cout<<"Event #"<<iev<<"   #Particles="<<nparticles<<endl;
+    rpcont->FetchClusters(0,TR);
+    TClonesArray *ITSrec  = NULL;
+    if(iev%100==0) printf("Event #%d\n",iev);
 
 
     Int_t ipt=0;
-    for(Int_t subd=0;subd<3;subd++){
-
-      Int_t first = geom->GetStartDet(subd);
-      Int_t last = geom->GetLastDet(subd);
-
-      for (Int_t mod=first; mod<=last; mod++){
-       detTypeRec->ResetRecPoints();
-       branch->GetEvent(mod);
-       Int_t nrecp = ITSrec->GetEntries();
-       if(nrecp>0){
-         for(Int_t irec=0;irec<nrecp;irec++) {
-           AliITSRecPoint *recp = (AliITSRecPoint*)ITSrec->At(irec);
-           Int_t lay=recp->GetLayer();
-           hlayer->Fill(lay);
-           recp->GetGlobalXYZ(cluglo);
-           Float_t rad=TMath::Sqrt(cluglo[0]*cluglo[0]+cluglo[1]*cluglo[1]); 
-           Float_t phi=TMath::ATan2(cluglo[1],cluglo[0]);
-           if(iev<3){
-             gpts[iev]->SetPoint(ipt,cluglo[0],cluglo[1]);
-             gpts3d[iev]->SetPoint(ipt,cluglo[0],cluglo[1],cluglo[2]);
-             ipt++;
-           }
-           hmod[lay]->Fill(mod);
-           hzl[lay]->Fill(recp->GetDetLocalZ());
-           hxl[lay]->Fill(recp->GetDetLocalX());
-           hzg[lay]->Fill(cluglo[2]);
-           hyg[lay]->Fill(cluglo[1]);
-           hxg[lay]->Fill(cluglo[0]);
-           hr[lay]->Fill(rad);
-           hphi[lay]->Fill(phi);
-           hq[lay]->Fill(recp->GetQ());
+    for (Int_t mod=modmin; mod<=modmax; mod++){
+      ITSrec = rpcont->UncheckedGetClusters(mod);
+      Int_t nrecp = ITSrec->GetEntries();
+      if(nrecp>0){
+       for(Int_t irec=0;irec<nrecp;irec++) {
+         AliITSRecPoint *recp = (AliITSRecPoint*)ITSrec->At(irec);
+         Int_t lay=recp->GetLayer();
+         hlayer->Fill(lay);
+         recp->GetGlobalXYZ(cluglo);
+         Float_t rad=TMath::Sqrt(cluglo[0]*cluglo[0]+cluglo[1]*cluglo[1]); 
+         Float_t phi=TMath::ATan2(cluglo[1],cluglo[0]);
+         if(iev==nevfordisp){
+           gptsXY->SetPoint(ipt,cluglo[0],cluglo[1]);
+           if(cluglo[1]>0) gptsRZ->SetPoint(ipt,cluglo[2],rad);
+           else gptsRZ->SetPoint(ipt,cluglo[2],-rad);
+           ipt++;
          }
+         hmod[lay]->Fill(mod);
+         hzl[lay]->Fill(recp->GetDetLocalZ());
+         hxl[lay]->Fill(recp->GetDetLocalX());
+         hzg[lay]->Fill(cluglo[2]);
+         hyg[lay]->Fill(cluglo[1]);
+         hxg[lay]->Fill(cluglo[0]);
+         hr[lay]->Fill(rad);
+         hphi[lay]->Fill(phi);
+         hq[lay]->Fill(recp->GetQ());
+         hdrtime[lay]->Fill(recp->GetDriftTime());
+         hzphi[lay]->Fill(phi,cluglo[2]);
        }
       }
     }
   }
+  
   TCanvas **c=new TCanvas*[6];
   Char_t ctit[30];
   for(Int_t iLay=0;iLay<6;iLay++){
@@ -211,24 +196,64 @@ Int_t ShowITSRecPoints(){
     hq[iLay]->Draw();    
   }
 
-  TCanvas *cev0;
-  cev0=new TCanvas("cev0","Event 0",600,600);
-  gpts[0]->SetMarkerStyle(7);
-  gpts[0]->SetTitle(0);
-  gpts[0]->Draw("AP");
+  TCanvas* cdrtim=new TCanvas("cdrtim","SDD drift time",1000,600);
+  cdrtim->Divide(2,1);
+  cdrtim->cd(1);
+  hdrtime[2]->Draw();
+  cdrtim->cd(2);
+  hdrtime[3]->Draw();
 
-  TCanvas *cev1;
-  cev1=new TCanvas("cev1","Event 1",600,600);
-  gpts[1]->SetMarkerStyle(7);
-  gpts[1]->SetTitle(0);
-  gpts[1]->Draw("AP");
+  gStyle->SetPalette(1);
+  TLatex* tstat=new TLatex();
+  tstat->SetNDC();
+  TCanvas* cspd=new TCanvas("cspd","SPD",1000,600);
+  cspd->Divide(2,1);
+  cspd->cd(1);
+  hzphi[0]->Draw("colz");
+  tstat->DrawLatex(0.6,0.92,Form("# Clusters = %d",int(hzphi[0]->GetEntries())));
+  cspd->cd(2);
+  hzphi[1]->Draw("colz");
+  tstat->DrawLatex(0.6,0.92,Form("# Clusters = %d",int(hzphi[1]->GetEntries())));
+
+  TCanvas* csdd=new TCanvas("csdd","SDD",1000,600);
+  csdd->Divide(2,1);
+  csdd->cd(1);
+  hzphi[2]->Draw("colz");
+  tstat->DrawLatex(0.6,0.92,Form("# Clusters = %d",int(hzphi[2]->GetEntries())));  
+  csdd->cd(2);
+  hzphi[3]->Draw("colz");
+  tstat->DrawLatex(0.6,0.92,Form("# Clusters = %d",int(hzphi[3]->GetEntries())));  
 
-  TCanvas *cev2;
-  cev2=new TCanvas("cev2","Event 2",600,600);
-  gpts[2]->SetMarkerStyle(7);
-  gpts[2]->SetTitle(0);
-  gpts[2]->Draw("AP");
+  TCanvas* cssd=new TCanvas("cssd","SSD",1000,600);
+  cssd->Divide(2,1);
+  cssd->cd(1);
+  hzphi[4]->Draw("colz");
+  tstat->DrawLatex(0.6,0.92,Form("# Clusters = %d",int(hzphi[4]->GetEntries())));  
+  cssd->cd(2);
+  hzphi[5]->Draw("colz");
+  tstat->DrawLatex(0.6,0.92,Form("# Clusters = %d",int(hzphi[5]->GetEntries())));  
 
+  TCanvas *cev0;
+  sprintf(ctit,"Event %d XY",nevfordisp);
+  cev0=new TCanvas("cev0",ctit,600,600);
+  if(gptsXY->GetN()>0){
+    gptsXY->SetMarkerStyle(7);
+    gptsXY->SetTitle(0);
+    gptsXY->Draw("AP");
+    gptsXY->GetXaxis()->SetTitle("Xglob");
+    gptsXY->GetYaxis()->SetTitle("Yglob");
+   }
+
+  TCanvas *cev1;
+  sprintf(ctit,"Event %d Zr",nevfordisp);
+  cev1=new TCanvas("cev1",ctit,600,600);
+  if(gptsRZ->GetN()>0){
+    gptsRZ->SetMarkerStyle(7);
+    gptsRZ->SetTitle(0);
+    gptsRZ->Draw("AP");
+    gptsRZ->GetXaxis()->SetTitle("Zglob");
+    gptsRZ->GetYaxis()->SetTitle("Radius");
+  }
 
   return 0;
 }