TRD/fastClusterAna.C

   1 void fastClusterAna() {
   2
   3 /////////////////////////////////////////////////////////////////////////
   4 //
   5 // Example macro for the analysis of the TRD cluster
   6 //
   7 /////////////////////////////////////////////////////////////////////////
   8
   9   // Dynamically link some shared libs
  10   if (gClassTable->GetID("AliRun") < 0) {
  11     gROOT->LoadMacro("loadlibs.C");
  12     loadlibs();
  13   }
  14
  15   // Input file name
  16   Char_t *alifile = "galice_r_v0.root";
  17
  18   // Event number
  19   Int_t   nEvent  = 0;
  20
  21   TH2F *HLocal  = new TH2F("HLocal" ,"rec. points local row/col-position"
  22                                     ,21,-0.5,20.5,81,-0.5,80.5);
  23   TH2F *HGlobal = new TH2F("HGlobal","rec. points global x/y-position"
  24                                     ,800,-400,400,800,-400,400);
  25
  26   // Connect the AliRoot file containing Geometry, Kine, Hits, and Digits
  27   TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
  28   if (!gafl) {
  29     cout << "Open the ALIROOT-file " << alifile << endl;
  30     gafl = new TFile(alifile);
  31   }
  32   else {
  33     cout << alifile << " is already open" << endl;
  34   }
  35
  36   // Get AliRun object from file or create it if not on file
  37   gAlice = (AliRun*) gafl->Get("gAlice");
  38   if (gAlice)
  39     cout << "AliRun object found on file" << endl;
  40   else
  41     gAlice = new AliRun("gAlice","Alice test program");
  42
  43   // Import the Trees for the event nEvent in the file
  44   Int_t nparticles = gAlice->GetEvent(nEvent);
  45   cout << "nparticles = " << nparticles << endl;
  46   if (nparticles <= 0) break;
  47
  48   // Get the pointer to the tree
  49   TTree          *RecTree       = gAlice->TreeR();
  50   RecTree->Print();
  51   // Get the pointer to the detector classes
  52   AliTRDv0       *TRD           = (AliTRDv0*) gAlice->GetDetector("TRD");
  53   // Get the geometry
  54   AliTRDgeometry *TRDgeometry   = TRD->GetGeometry();
  55   // Get the pointer to the hit container
  56   TObjArray      *RecPointArray = TRD->RecPoints();
  57   // Set the branch address
  58   RecTree->GetBranch("TRDrecPoints")->SetAddress(&RecPointArray);
  59
  60   Int_t nEntries = RecTree->GetEntries();
  61   cout << "Number of entries in reconstruction tree = " << nEntries << endl;
  62
  63   // Loop through all entries in the tree
  64   Int_t nbytes;
  65   for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) {
  66
  67     cout << "iEntry = " << iEntry << endl;
  68
  69     // Import the tree
  70     nbytes += RecTree->GetEvent(iEntry);
  71
  72     // Get the number of points in the detector
  73     Int_t nRecPoint = RecPointArray->GetEntriesFast();
  74     cout << " nRecPoint = " << nRecPoint << endl;
  75
  76     // Loop through all TRD digits
  77     for (Int_t iRecPoint = 0; iRecPoint < nRecPoint; iRecPoint++) {
  78
  79       // Get the information for this digit
  80       AliTRDrecPoint *RecPoint = (AliTRDrecPoint *) RecPointArray->UncheckedAt(iRecPoint);
  81       Int_t    detector = RecPoint->GetDetector();
  82       Float_t  row      = RecPoint->GetLocalRow();
  83       Float_t  col      = RecPoint->GetLocalCol();
  84
  85       Int_t    sector   = TRDgeometry->GetSector(detector);
  86       Int_t    plane    = TRDgeometry->GetPlane(detector);
  87       Int_t    chamber  = TRDgeometry->GetChamber(detector);
  88
  89       TVector3 Pos;
  90       TMatrix  Cov;
  91       RecPoint->GetGlobalPosition(Pos,Cov);
  92       HGlobal->Fill(Pos.X(),Pos.Y());
  93       if ((sector == 17) && (plane == 0) && (chamber == 2)) {
  94         HLocal->Fill(row,col);
  95       }
  96
  97     }
  98
  99   }
 100
 101   TCanvas *C = new TCanvas("C","recPoints",10,10,400,600);
 102   C->Divide(1,2);
 103   C->cd(1);
 104   HLocal->Draw("BOX");
 105   C->cd(2);
 106   HGlobal->Draw("BOX");
 107
 108 }