]>
Commit | Line | Data |
---|---|---|
a9e2aefa | 1 | void MUONhits (Int_t evNumber1=0,Int_t evNumber2=0, Int_t ic=1) |
2 | { | |
3 | ///////////////////////////////////////////////////////////////////////// | |
4 | // This macro is a small example of a ROOT macro | |
5 | // illustrating how to read the output of GALICE | |
6 | // and do some analysis. | |
7 | // | |
8 | ///////////////////////////////////////////////////////////////////////// | |
9 | ||
10 | // Dynamically link some shared libs | |
11 | Float_t rmin[14] = {17.5, 17.5, 23.5, 23.5, 33.5, 33.5, 43., 43., 50., 50, | |
12 | 56.1, 56.1, 59.6, 59.6}; | |
13 | Float_t rmax[14] = {81.6, 81.6, 109.3, 109.3, 154.4, 154.4, 197.8, | |
14 | 197.8, 229.5, 229.5, 254., 254, 270., 270.}; | |
15 | ||
16 | ||
17 | if (gClassTable->GetID("AliRun") < 0) { | |
18 | gROOT->LoadMacro("loadlibs.C"); | |
19 | loadlibs(); | |
20 | } | |
21 | // Connect the Root Galice file containing Geometry, Kine and Hits | |
22 | ||
23 | TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root"); | |
24 | ||
25 | ||
26 | if (!file) { | |
27 | printf("\n Creating galice.root \n"); | |
28 | file = new TFile("galice.root"); | |
29 | } else { | |
30 | printf("\n galice.root found in file list"); | |
31 | } | |
32 | file->ls(); | |
33 | // | |
34 | TDatabasePDG* DataBase = new TDatabasePDG(); | |
35 | ||
36 | // Get AliRun object from file or create it if not on file | |
37 | if (!gAlice) { | |
38 | gAlice = (AliRun*)(file->Get("gAlice")); | |
39 | if (gAlice) printf("AliRun object found on file\n"); | |
40 | if (!gAlice) { | |
41 | printf("\n create new gAlice object"); | |
42 | gAlice = new AliRun("gAlice","Alice test program"); | |
43 | } | |
44 | } | |
45 | ||
46 | // Create some histograms | |
47 | ||
48 | TH1F *Hits[14], *HitDensity[14]; | |
49 | TH1F *Hits_g[14], *HitDensity_g[14]; | |
50 | TH1F *Hits_n[14], *HitDensity_n[14]; | |
51 | TH1F *Mom[14], *Mom_g[14], *Mom_n[14]; | |
52 | ||
53 | for (Int_t i=0; i<14; i++) { | |
54 | Hits[i] = new TH1F("hits1","Hit Distribution",30,0,300); | |
55 | Hits[i]->SetFillColor(0); | |
56 | Hits[i]->SetXTitle("R (cm)"); | |
57 | ||
58 | HitDensity[i] = new TH1F("dhits1","Hit Density",30,0,300); | |
59 | HitDensity[i]->SetFillColor(0); | |
60 | HitDensity[i]->SetXTitle("R (cm)"); | |
61 | ||
62 | Hits_g[i] = new TH1F("hits1","Hit Distribution",30,0,300); | |
63 | Hits_g[i]->SetFillColor(0); | |
64 | Hits_g[i]->SetXTitle("R (cm)"); | |
65 | ||
66 | HitDensity_g[i] = new TH1F("dhits1","Hit Density",30,0,300); | |
67 | HitDensity_g[i]->SetFillColor(0); | |
68 | HitDensity_g[i]->SetXTitle("R (cm)"); | |
69 | ||
70 | Mom[i] = new TH1F("mom","Energy",70,-6,1); | |
71 | Mom[i]->SetFillColor(0); | |
72 | Mom[i]->SetXTitle("E (GeV)"); | |
73 | ||
74 | Mom_g[i] = new TH1F("mom","Energy",70,-6,1); | |
75 | Mom_g[i]->SetFillColor(0); | |
76 | Mom_g[i]->SetXTitle("E (GeV)"); | |
77 | ||
78 | Mom_n[i] = new TH1F("mom","Energy",70,-6,1); | |
79 | Mom_n[i]->SetFillColor(0); | |
80 | Mom_n[i]->SetXTitle("E (GeV)"); | |
81 | ||
82 | Hits_n[i] = new TH1F("hits1","Hit Distribution",30,0,300); | |
83 | Hits_n[i]->SetFillColor(0); | |
84 | Hits_n[i]->SetXTitle("R (cm)"); | |
85 | ||
86 | HitDensity_n[i] = new TH1F("dhits1","Hit Density",30,0,300); | |
87 | HitDensity_n[i]->SetFillColor(0); | |
88 | HitDensity_n[i]->SetXTitle("R (cm)"); | |
89 | } | |
90 | ||
91 | TH1F *theta = new TH1F("theta","Theta distribution",180,0,180); | |
92 | ||
93 | TH1F *emult = new TH1F("emult","Event Multiplicity",100,0,1000); | |
94 | AliMUONChamber* iChamber; | |
d6c607d0 | 95 | AliSegmentation* seg; |
a9e2aefa | 96 | |
97 | // | |
98 | // Loop over events | |
99 | // | |
100 | Float_t dpx[2], dpy[2]; | |
101 | Int_t mdig =0; | |
102 | Int_t nhit=0; | |
103 | ||
104 | for (Int_t nev=0; nev<= evNumber2; nev++) { | |
105 | Int_t nparticles = gAlice->GetEvent(nev); | |
106 | if (nev < evNumber1) continue; | |
107 | if (nparticles <= 0) return; | |
108 | ||
109 | AliMUON *MUON = (AliMUON*) gAlice->GetModule("MUON"); | |
110 | printf("\n track %d %d \n ", nev, MUON); | |
111 | ||
112 | TTree *TH = gAlice->TreeH(); | |
113 | Int_t ntracks = TH->GetEntries(); | |
114 | // | |
115 | // Loop over events | |
116 | // | |
117 | Int_t EvMult=0; | |
118 | ||
119 | for (Int_t track=0; track<ntracks;track++) { | |
120 | gAlice->ResetHits(); | |
121 | Int_t nbytes += TH->GetEvent(track); | |
122 | if (MUON) { | |
123 | // | |
124 | // Loop over hits | |
125 | // | |
126 | Int_t NperCh=0; | |
127 | ||
128 | for(AliMUONHit* mHit=(AliMUONHit*)MUON->FirstHit(-1); | |
129 | mHit; | |
130 | mHit=(AliMUONHit*)MUON->NextHit()) | |
131 | { | |
bd21d472 | 132 | Int_t nch = mHit->Chamber(); // chamber number |
133 | Float_t x = mHit->X(); // x-pos of hit | |
134 | Float_t y = mHit->Y(); // y-pos | |
135 | Float_t Eloss = mHit->Eloss(); | |
136 | Float_t Theta = mHit->Theta(); | |
137 | Float_t Particle = mHit->Particle(); | |
1e63a707 | 138 | Float_t P = mHit->Momentum(); |
a9e2aefa | 139 | TParticlePDG* Part = DataBase->GetParticle(Particle); |
140 | Double_t mass = Part->Mass(); | |
141 | ||
142 | if (nch >13) continue; | |
143 | if (nch ==1) EvMult++; | |
bd21d472 | 144 | if (mHit->Age() > 5.e-6) continue; |
a9e2aefa | 145 | |
146 | Float_t r=TMath::Sqrt(x*x+y*y); | |
147 | if (nch ==0) continue; | |
148 | ||
149 | if (r < rmin[nch-1]) continue; | |
150 | if (r > rmax[nch-1]) continue; | |
151 | ||
152 | Float_t wgt=1/(2*10*TMath::Pi()*r)/(evNumber2+1); | |
153 | ||
154 | Float_t Ekin=TMath::Sqrt(P*P+mass*mass)-mass; | |
155 | if (Particle == 2112) { | |
156 | Hits_n[nch-1]->Fill(r,(float) 1); | |
157 | HitDensity_n[nch-1]->Fill(r,wgt); | |
158 | Mom_n[nch-1]->Fill(TMath::Log10(Ekin),1); | |
159 | } else if (Particle == 22) { | |
160 | Hits_g[nch-1]->Fill(r,(float) 1); | |
161 | HitDensity_g[nch-1]->Fill(r,wgt); | |
162 | Mom_g[nch-1]->Fill(TMath::Log10(Ekin),1); | |
163 | } else { | |
164 | Hits[nch-1]->Fill(r,(float) 1); | |
165 | HitDensity[nch-1]->Fill(r,wgt); | |
166 | Mom[nch-1]->Fill(TMath::Log10(Ekin),1); | |
167 | } | |
168 | } // hit loop | |
169 | } // if MUON | |
170 | } // track loop | |
171 | } | |
172 | ||
173 | //Create a canvas, set the view range, show histograms | |
174 | Int_t k; | |
175 | TCanvas *c1 = new TCanvas("c1","Hit Densities",400,10,600,700); | |
176 | c1->Divide(2,4); | |
177 | for (k=0; k<7; k++) { | |
178 | c1->cd(k+1); | |
179 | HitDensity[2*k]->Draw(); | |
180 | } | |
181 | ||
182 | TCanvas *c2 = new TCanvas("c2","Hit Densities (gamma)",400,10,600,700); | |
183 | c2->Divide(2,4); | |
184 | for (k=0; k<7; k++) { | |
185 | c2->cd(k+1); | |
186 | HitDensity_g[2*k]->Draw(); | |
187 | } | |
188 | ||
189 | TCanvas *c3 = new TCanvas("c3","Hit Densities (neutron)",400,10,600,700); | |
190 | c3->Divide(2,4); | |
191 | for (k=0; k<7; k++) { | |
192 | c3->cd(k+1); | |
193 | HitDensity_n[2*k]->Draw(); | |
194 | } | |
195 | ||
196 | TCanvas *c4 = new TCanvas("c4","Energy (charged)",400,10,600,700); | |
197 | c4->Divide(2,4); | |
198 | for (k=0; k<7; k++) { | |
199 | c4->cd(k+1); | |
200 | Mom[2*k]->Draw(); | |
201 | } | |
202 | ||
203 | TCanvas *c5 = new TCanvas("c5","Energy (gamma)",400,10,600,700); | |
204 | c5->Divide(2,4); | |
205 | for (k=0; k<7; k++) { | |
206 | c5->cd(k+1); | |
207 | Mom_g[2*k]->Draw(); | |
208 | } | |
209 | ||
210 | TCanvas *c6 = new TCanvas("c6","Energy (gamma)",400,10,600,700); | |
211 | c6->Divide(2,4); | |
212 | for (k=0; k<7; k++) { | |
213 | c6->cd(k+1); | |
214 | Mom_n[2*k]->Draw(); | |
215 | } | |
216 | } | |
217 | ||
218 | ||
219 | ||
220 | ||
221 | ||
222 | ||
223 |