ITS/ITSOccupancy.C

   1 /******************************************************************************
   2
   3   "ITSOccupancy.C"
   4
   5   this macro calculates the mean occupancy of each ITS layer, counting the
   6   "fired" digits of each module, and making two overall calculi:
   7      1) the calculus of the mean overall layer occupancy (as the ratio
   8              between the total number of active digits and the total number of
   9                   channels in the module;
  10           2) a distribution of the occupancies as a funcion of z, to obtain some
  11                   kind of most significand data for this value, along the z axis
  12
  13   The macro creates also the GIF of the TCanvas where the histograms are plotted
  14
  15   The arguments that must be inserted are the following:
  16      1) the answer to the following question (as a Bool_t)
  17                   "Can I get the total channels distribution from the file totals.root?"
  18                   if the answer is NO (false) these histograms are re-created and re-stored
  19                   in that file.
  20
  21           2) the name of the file to analyze;
  22              NOTE: in this argument you MUST omit the ".root" extension, because the
  23                         macro uses this argument to recall the file, but also to name the
  24                                   GIF & PS of the output canvas, in order to save the plot of the
  25                                   distributions to the disk. If you put, as argument, the string
  26                                   "file.root", the macro will look for an unexisting file named
  27                                   "file.root.root"...
  28
  29   It is also possible (and recommended) to compile this macro,
  30   to execute it faster. To do this, it is necessary to:
  31   1) make sure that de "#define __COMPILED__" instruction below is UNcommented
  32      and the "#define __NOCOMPILED__" instruction is commented (remember that
  33           if you comment or uncomment both these instructions, the macro can't work)
  34   2) type, at the root prompt, the instruction
  35      gSystem->SetIncludePath("-I- -I$ALICE_ROOT/include -I$ALICE_ROOT/ITS -g")
  36           to adjoin the ALICE include libraries to the main include directory
  37   3) load the function via the ".L ITSOccupancy.C++" statement (only the first
  38      time, because the first time the macro must be compiled). Frome the
  39           second time you use the macro, you must only load it via the
  40           ".L ITSOccupancy.C+" instruction (note the number of '+' signs in each case
  41   4) call the function with only its name, es: ITSOccupancy("HP", "galice.root")
  42
  43  Author: Alberto Pulvirenti
  44 ******************************************************************************/
  45
  46 #define __COMPILED__
  47 //#define __NOCOMPILED__
  48
  49 #ifdef __COMPILED__
  50         #include <iostream.h>
  51         #include <TROOT.h>
  52         #include <TArrayI.h>
  53         #include <TCanvas.h>
  54         #include <TClassTable.h>
  55         #include <TClonesArray.h>
  56         #include <TFile.h>
  57         #include <TObject.h>
  58         #include <TObjArray.h>
  59         #include <TTree.h>
  60         #include <TMath.h>
  61         #include <TString.h>
  62         #include <TH1.h>
  63         #include <AliRun.h>
  64         #include <AliITS.h>
  65         #include <AliITSgeom.h>
  66         #include <AliITSDetType.h>
  67         #include <AliITSRecPoint.h>
  68         #include <AliITSdigit.h>
  69         #include <AliITShit.h>
  70         #include <AliITSmodule.h>
  71         #include <AliITSsegmentation.h>
  72         #include <AliITSsegmentationSPD.h>
  73         #include <AliITSsegmentationSDD.h>
  74         #include <AliITSsegmentationSSD.h>
  75 #endif
  76
  77 Int_t ITSOccupancy(char *name = "galice", Int_t evNum = 0) {
  78
  79         // Evaluate the name of the file to open...
  80         Int_t filelen = strlen(name);
  81         Text_t *fileroot = new Text_t[filelen + 6];
  82         // ...and copy this name to the data ROOT file
  83         strcpy (fileroot, name);
  84         strcat (fileroot, ".root");
  85
  86         // Open the Root input file containing the AliRun data and
  87         // the Root output data file
  88         TFile *froot = (TFile*)gROOT->GetListOfFiles()->FindObject(fileroot);
  89         if (!froot) froot = new TFile(fileroot, "READ");
  90         // Load the gAlice shared libs if not already in memory and
  91         // clears an eventually existing AliRun object
  92 #ifdef __NOCOMPILED__
  93         if (gClassTable->GetID("AliRun") < 0) {
  94         gROOT->LoadMacro("loadlibs.C");
  95         loadlibs();
  96         }
  97 #endif
  98    if (gAlice) {
  99       delete gAlice;
 100            gAlice = 0;
 101    }
 102         // Get the AliRun object from file (if there's none, the macro is stopped)
 103         Int_t nparticles = 0;
 104         gAlice = (AliRun*)froot->Get("gAlice");
 105         if (gAlice) {
 106                 cout << "Found an AliRun object in `" << fileroot << "'" << endl;
 107                 // In this case, select the event number specified (default is 0)
 108                 nparticles = gAlice->GetEvent(evNum);
 109                 if (nparticles)
 110                         cout << "\nNumber of particles   = " << nparticles << endl;
 111                 else {
 112                         cout << "NO PARTICLES!!!!" << endl;
 113                         return 0;
 114                 }
 115         }
 116         else {
 117                 cout<<"Not found any AliRun object in `" << fileroot << "'" << endl;
 118                 return 0;
 119         }
 120         // Initialize the ITS and its geometry
 121         AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
 122         AliITSgeom *gm = ITS->GetITSgeom();
 123
 124         // Variables and objects definition and setting
 125         Float_t zmax[] = {20.0,20.0,25.0,35.0,49.5,54.0}; // z edges for TH1F (cm)
 126         Int_t nbins[]  = {20,20,20,28,22,24};             // bins number for TH1Fs
 127
 128         // Histos for plotting occupancy distributions
 129         TH1F *above[6], *below[6];
 130
 131         for (Int_t lay = 0; lay < 6; lay++) {
 132                 Int_t nlads = gm->GetNladders(lay+1);
 133                 Int_t ndets = gm->GetNdetectors(lay+1);
 134                 Int_t dtype = lay / 2;
 135                 Int_t minID = gm->GetModuleIndex(lay+1, 1, 1);
 136                 Int_t maxID = gm->GetModuleIndex(lay+1, nlads, ndets);
 137                 Text_t ffname[20];
 138                 sprintf(ffname, "h_%d", lay+1);
 139                 below[lay] = new TH1F(ffname, "Total z distribution of digits", nbins[lay], -zmax[lay], zmax[lay]);
 140                 cout << "Evaluating total channels number of layer " << lay+1 << endl;
 141                 for (Int_t I = minID; I <= maxID; I++) {
 142                         AliITSDetType *det = ITS->DetType(dtype);
 143                         AliITSsegmentation *seg = det->GetSegmentationModel();
 144                         Int_t NX = seg->Npx();
 145                         Int_t NZ = seg->Npz();
 146                         for (Int_t ix = 0; ix <= NX; ix++) {
 147                                 for (Int_t iz = 0; iz <= NZ; iz++) {
 148                                         Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
 149                                         seg->DetToLocal(ix, iz, lx[0], lx[2]);
 150                                         gm->LtoG(I, lx, gx);
 151                                         below[lay]->Fill(gx[2]);
 152                                 }
 153                         }
 154                 }
 155                 sprintf(ffname, "H_%d", lay+1);
 156                 above[lay] = new TH1F(ffname, "histo", nbins[lay], -zmax[lay], zmax[lay]);
 157         }
 158
 159         // Counting the hits, digits and recpoints contained in the ITS
 160         TTree *TD = gAlice->TreeD();
 161         ITS->ResetDigits();
 162         Float_t mean[6];
 163
 164         for (Int_t L = 0; L < 6; L++) {
 165
 166                 cout << "Layer " << L + 1 << ": " << flush;
 167
 168                 // To avoid two nested loops, are calculated
 169                 // the ID of the first and last module of the L
 170                 // (remember the L goest from 0 to 5 (not from 1 to 6)
 171                 Int_t first, last;
 172                 first = gm->GetModuleIndex(L + 1, 1, 1);
 173                 last = gm->GetModuleIndex(L + 1, gm->GetNladders(L + 1), gm->GetNdetectors(L + 1));
 174
 175                 // Start loop on modules
 176                 for (Int_t ID = first; ID <= last; ID++) {
 177                         Int_t dtype = L / 2;
 178                         AliITSDetType *det = ITS->DetType(dtype);
 179                         AliITSsegmentation *seg = det->GetSegmentationModel();
 180                         if (dtype == 2) seg->SetLayer(L+1);
 181
 182                         TD->GetEvent(ID);
 183                         TClonesArray *digits_array = ITS->DigitsAddress(dtype);
 184                         Int_t digits_num = digits_array->GetEntries();
 185                         // Get the coordinates of the module
 186                         for (Int_t j = 0; j < digits_num; j++) {
 187                                 Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
 188                                 AliITSdigit *digit = (AliITSdigit*)digits_array->UncheckedAt(j);
 189                                 Int_t iz=digit->fCoord1;  // cell number z
 190                                 Int_t ix=digit->fCoord2;  // cell number x
 191                         // Get local coordinates of the element (microns)
 192                         seg->DetToLocal(ix, iz, lx[0], lx[2]);
 193                                 gm->LtoG(ID, lx, gx);
 194                                 above[L]->Fill(gx[2]);
 195                         }
 196                 }
 197
 198                 above[L]->Divide(above[L], below[L], 100.0, 1.0);
 199                 mean[L] = above[L]->GetSumOfWeights() / above[L]->GetNbinsX();
 200                 cout.setf(ios::fixed);
 201                 cout.precision(2);
 202                 cout << " digits occupancy = " << mean[L] << "%" << endl;
 203         }
 204
 205         TCanvas *view = new TCanvas("view", "Digits occupancy distributions", 600, 900);
 206         view->Divide(2, 3);
 207
 208         for (Int_t I = 1; I < 7; I++) {
 209                 view->cd(I);
 210                 Text_t title[50];
 211                 sprintf(title, "Layer %d: %4.2f%c", I, mean[I-1], '%');
 212                 title[6] = (Text_t)I + '0';
 213                 above[I-1]->SetTitle(title);
 214                 above[I-1]->SetStats(kFALSE);
 215                 above[I-1]->SetXTitle("z (cm)");
 216                 above[I-1]->SetYTitle("%");
 217                 above[I-1]->Draw();
 218                 view->Update();
 219         }
 220         Text_t *filegif  = new Text_t[filelen + 23];
 221         Text_t *fileps  = new Text_t[filelen + 23];
 222         strcpy (filegif, name);
 223         strcpy (fileps, name);
 224         strcat (filegif, "_digit_occupancy.gif");
 225         strcat (fileps, "_digit_occupancy.eps");
 226         view->SaveAs(filegif);
 227         view->SaveAs(fileps);
 228
 229
 230         return 1;
 231 }