Config.C for MUON studies
[u/mrichter/AliRoot.git] / macros / anal.C
1 void anal (Int_t evNumber=0) 
2 {
3 /////////////////////////////////////////////////////////////////////////
4 //   This macro is a small example of a ROOT macro
5 //   illustrating how to read the output of GALICE
6 //   and fill some histograms.
7 //   
8 //     Root > .L anal.C   //this loads the macro in memory
9 //     Root > anal();     //by default process first event   
10 //     Root > anal(2);    //process third event
11 //Begin_Html
12 /*
13 <img src="picts/anal.gif">
14 */
15 //End_Html
16 /////////////////////////////////////////////////////////////////////////
17
18
19 // Dynamically link some shared libs
20    if (gClassTable->GetID("AliRun") < 0) {
21       gROOT->LoadMacro("loadlibs.C");
22       loadlibs();
23    }
24       
25 // Connect the Root Galice file containing Geometry, Kine and Hits
26    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
27    if (!file) file = new TFile("galice.root");
28
29 // Get AliRun object from file or create it if not on file
30    if (!gAlice) {
31       gAlice = (AliRun*)file->Get("gAlice");
32       if (gAlice) printf("AliRun object found on file\n");
33       if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
34    }
35       
36 // Import the Kine and Hits Trees for the event evNumber in the file
37    Int_t nparticles = gAlice->GetEvent(evNumber);
38    if (nparticles <= 0) return;
39    Float_t x,y,z,mass,e;
40    Int_t nbytes = 0;
41    Int_t j,hit,ipart;
42    Int_t nhits;
43    Int_t sector,plane;
44    TParticle *particle;
45    AliTPChit  *tpcHit;
46    AliTRDhit  *trdHit;
47
48 // Get pointers to Alice detectors and Hits containers
49    AliDetector *TPC  = gAlice->GetDetector("TPC");
50    AliDetector *TRD  = gAlice->GetDetector("TRD");
51    TClonesArray *Particles = gAlice->Particles();
52    if (TPC) TClonesArray *TPChits   = TPC->Hits();
53    if (TRD) TClonesArray *TRDhits   = TRD->Hits();
54
55    TTree *TH = gAlice->TreeH();
56    Int_t ntracks    = TH->GetEntries();
57
58    // Create histograms
59    TH1F *hSectors = new TH1F("hSectors","Number of tracks hits per sector",75,1,76);
60    TH1F *hTRD     = new TH1F("hTRD","Number of planes crossed per track in TRD",6,1,7);
61    TH1F *hTRDprim = new TH1F("hTRDprim","Number of planes crossed per primary track in TRD",6,1,7);
62    
63 // Start loop on tracks in the hits containers
64    for (Int_t track=0; track<ntracks;track++) {
65      gAlice->ResetHits();
66      nbytes += TH->GetEvent(track);
67      // ======>Histogram TPC
68      // histogram number of tracks per sector
69      if (TPC) {
70        nhits = TPChits->GetEntriesFast();
71        for (hit=0;hit<nhits;hit++) {
72          tpcHit   = (AliTPChit*)TPChits->UncheckedAt(hit);
73          sector = tpcHit->fPadRow;
74          hSectors->Fill(sector);
75        }
76      }        
77      // =======>Histogram TRD
78      // histogram number of planes crossed per track
79      // same for primary tracks only
80      if (TRD) {
81        nhits = TRDhits->GetEntriesFast();
82        for (hit=0;hit<nhits;hit++) {
83          trdHit   = (AliTRDhit*)TRDhits->UncheckedAt(hit);
84          ipart    = trdHit->fTrack;
85          particle = (TParticle*)Particles->UncheckedAt(ipart);
86          plane = trdHit->fPlane;
87          hTRD->Fill(plane);
88          if (particle->GetFirstMother() < 0) hTRDprim->Fill(plane);
89        }
90      }        
91    }
92
93 //Create a canvas, set the view range, show histograms
94    TCanvas *c1 = new TCanvas("c1","Alice TPC and ITS hits",400,10,600,700);
95    c1->Divide(1,2);
96    c1->cd(1);
97    gPad->SetFillColor(33);
98    hSectors->SetFillColor(42);
99    hSectors->Fit("pol1");
100    c1->cd(2);
101    gPad->SetFillColor(33);
102    hTRD->SetFillColor(42);
103    hTRD->Draw();
104    hTRDprim->SetFillColor(46);
105    hTRDprim->Draw("same");
106    c1->Print("anal.ps");
107 }