]>
Commit | Line | Data |
---|---|---|
fd2bbc37 | 1 | void EvalAbso() |
2 | { | |
3 | // Material numbers | |
4 | enum {kC=6, kAl=9, kFe=10, kCu=11, kW=12, kPb=13, | |
5 | kNiCuW=21, kVacuum=16, kAir=15, kConcrete=17, | |
6 | kPolyCH2=18, kSteel=10, kInsulation=14, kPolyCc=20}; | |
7 | ||
8 | AliABSO *pABSO = (AliABSO*) gAlice->GetModule("ABSO"); | |
9 | ||
10 | Int_t nl[2]; | |
11 | nl[0] = pABSO->NumberOfLayers(0); | |
12 | nl[1] = pABSO->NumberOfLayers(1); | |
13 | ||
14 | Int_t mat[2][15], mat[2][15]; | |
15 | Float_t zmin[2][15], zmin[2][15], zmax[2][15], zmax[2][15]; | |
16 | Int_t i, j; | |
17 | ||
18 | for (j=0; j< 2; j++) { | |
19 | for (i=0; i< nl[j]; i++) { | |
20 | mat[j][i] = pABSO->MaterialOfLayer(j,i)-1599; | |
21 | zmax[j][i] = pABSO->ZPositionOfLayer(j,i); | |
22 | if (i+1 < nl[j]) zmin[j][i+1] = pABSO->ZPositionOfLayer(j,i); | |
23 | } | |
24 | zmin[j][0]=0.; | |
25 | } | |
26 | Float_t l = zmax[0][nl[0]-1]; | |
27 | ||
28 | ||
29 | // | |
30 | ||
31 | // | |
32 | // 1. Limiting angular resolution in the Branson formalism | |
33 | // | |
34 | Float_t f0,f1,f2; | |
35 | Float_t a, z, dens, radl, absl; | |
36 | char name[21]; | |
37 | ||
38 | for (j=0; j< 2; j++) { | |
92664bc8 | 39 | printf("\n A Z ZPos DZ dens radL absL"); |
fd2bbc37 | 40 | printf("\n=================================================================================="); |
41 | f0=f1=f2=0; | |
42 | for (i=0; i< nl[j]; i++) { | |
43 | pABSO->AliGetMaterial(mat[j][i], name, a, z, dens, radl, absl); | |
44 | Float_t dz = zmax[j][i]-zmin[j][i]; | |
45 | f0 += dz/radl; | |
46 | f1 += dz/(2.*radl)*(2.*zmin[j][i]+dz); | |
47 | f2 += dz/(3.*radl)*(3.*zmin[j][i]*zmin[j][i]+3.*dz*zmin[j][i]+dz*dz); | |
92664bc8 | 48 | printf("\n %3d %14s %12.3f %12.3f %12.3f %12.3f %12.3f %12.3f %12.3f", |
49 | i+1, name, a, z, (zmax[j][i]+zmin[j][i])/2., dz, dens, radl, absl); | |
fd2bbc37 | 50 | } |
51 | // | |
52 | Float_t deltaThetaB = 0.0136 * TMath::Sqrt(f0-f1*f1/f2); | |
53 | Float_t lBranson = f2/f1; | |
54 | printf("\n=================================================================================="); | |
55 | printf("\n Branson plane: %12.3f (cm)", lBranson); | |
56 | printf("\n Limiting resolution (Branson) : %12.3f (mrad)", deltaThetaB*1000.); | |
57 | } | |
58 | // | |
59 | // 2. Limiting Resolution as a function of thickness of extra Fe | |
60 | // | |
61 | Float_t zmin0[15], zmax0[15]; | |
62 | Float_t x[100], y[100]; | |
63 | ||
64 | for (i=0; i< nl[0]; i++) { | |
65 | zmin0[i] = zmin[0][i]; | |
66 | zmax0[i] = zmax[0][i]; | |
67 | } | |
68 | ||
69 | Float_t ds = (zmax[0][3]-zmin[0][2])/100.; | |
70 | Int_t j; | |
71 | ||
72 | for (j=0; j< 100; j++) { | |
73 | f0=f1=f2=0; | |
74 | Float_t zFe = zmin[0][2]+Float_t(j)*ds; | |
75 | zmax0[2]=zFe; | |
76 | zmin0[3]=zFe; | |
77 | for (i=0; i< nl[0]; i++) { | |
78 | Float_t a, z, dens, radl, absl; | |
79 | pABSO->AliGetMaterial(mat[0][i], name, a, z, dens, radl, absl); | |
80 | Float_t dz = zmax0[i]-zmin0[i]; | |
81 | f0 += dz/radl; | |
82 | f1 += dz/(2.*radl)*(2.*zmin0[i]+dz); | |
83 | f2 += dz/(3.*radl)*(3.*zmin0[i]*zmin0[i]+3.*dz*zmin0[i]+dz*dz); | |
84 | } | |
85 | Float_t deltaThetaB = 0.0136 * TMath::Sqrt(f0-f1*f1/f2); | |
86 | x[j] = zmax[0][3]-zFe; | |
87 | y[j] = deltaThetaB*1000.; | |
88 | } | |
89 | TGraph *sGraph = new TGraph(100, x, y); | |
90 | TCanvas *c = new TCanvas("c"," ",400,10,600,700); | |
91 | c->Divide(2,2); | |
92 | c->cd(1); | |
93 | sGraph->Draw("AC"); | |
94 | c->Update(); | |
95 | // | |
96 | // 3. Limiting Resolution as a function of density of Carbon | |
97 | // | |
98 | Float_t dRho = (2.5-1.5)/100.; | |
99 | ||
100 | for (j=0; j< 100; j++) { | |
101 | f0=f1=f2=0; | |
102 | Float_t rho = 1.5 + Float_t(j)*dRho; | |
103 | for (i=0; i< nl[0]; i++) { | |
104 | Float_t a, z, dens, radl, absl; | |
105 | pABSO->AliGetMaterial(mat[0][i], name, a, z, dens, radl, absl); | |
106 | if (i==1) radl*=1.75/rho; | |
107 | Float_t dz = zmax[0][i]-zmin[0][i]; | |
108 | f0 += dz/radl; | |
109 | f1 += dz/(2.*radl)*(2.*zmin[0][i]+dz); | |
110 | f2 += dz/(3.*radl)*(3.*zmin[0][i]*zmin[0][i]+3.*dz*zmin[0][i]+dz*dz); | |
111 | } | |
112 | deltaThetaB = 0.0136 * TMath::Sqrt(f0-f1*f1/f2); | |
113 | x[j] = rho; | |
114 | y[j] = deltaThetaB*1000.; | |
115 | } | |
116 | TGraph *dGraph = new TGraph(100, x, y); | |
117 | c->cd(2); | |
118 | dGraph->Draw("AC"); | |
119 | c->Update(); | |
120 | ||
121 | // | |
122 | // 4. Energy Loss | |
123 | // | |
124 | // | |
125 | TCanvas *c1 = new TCanvas("c1"," ",400,10,600,700); | |
126 | c1->Divide(2,2); | |
127 | ||
128 | const Float_t kAvo=0.60221367; | |
129 | ||
130 | char* tmp; | |
131 | tmp = new char[5]; | |
132 | strncpy(tmp, "LOSS", 4); | |
133 | tmp[4]='\0'; | |
134 | ||
135 | Int_t ixst; | |
136 | Float_t ekin[100], dedx[100], de[100], pcut[100], deRad[100]; | |
137 | Float_t emin = 2; | |
138 | Float_t emax = 200; | |
139 | Float_t eps = (emax-emin)/100.; | |
140 | // 3 < theta < 9 | |
141 | for (j=0; j< 100; j++) { | |
142 | ekin[j] = emin + Float_t(j)*eps; | |
143 | de[j]=0; | |
144 | deRad[j]=0; | |
145 | } | |
146 | // all losses | |
147 | for (i=0; i< nl[0]; i++) { | |
148 | ((TGeant3*) gMC)->Gftmat(pABSO->GetMatId(mat[0][i]), 5, tmp, 100, ekin, dedx, pcut, ixst); | |
149 | for (j=0; j< 100; j++) de[j]+=dedx[j]*(zmax[0][i]-zmin[0][i])/1000.; | |
150 | } | |
151 | // radative losses | |
152 | for (i=0; i< nl[0]; i++) { | |
153 | Float_t a, z, dens, radl, absl; | |
154 | pABSO->AliGetMaterial(mat[0][i], name, a, z, dens, radl, absl); | |
155 | for (Int_t j=0; j<100; j++) { | |
156 | Float_t nor= kAvo*dens/a*(zmax[0][i]-zmin[0][i]); | |
157 | deRad[j] += (((TGeant3*) gMC)->Gbrelm(z,ekin[j],1.e10))*nor; | |
158 | deRad[j] += (((TGeant3*) gMC)->Gprelm(z,ekin[j],1.e10))*nor; | |
159 | } | |
160 | } | |
161 | ||
162 | TH2F *hr = new TH2F("hr1","Several graphs in the same pad",2,0,15,2,0,5); | |
163 | hr->SetXTitle("X title"); | |
164 | hr->SetYTitle("Y title"); | |
165 | c1->cd(1); | |
166 | hr->Draw(); | |
167 | ||
168 | TGraph *eGraph1 = new TGraph(100, ekin, de); | |
169 | TGraph *eGraphR1 = new TGraph(100, ekin, deRad); | |
170 | ||
171 | eGraph1 ->Draw("LP"); | |
172 | eGraphR1->Draw("LP"); | |
173 | ||
174 | hr = new TH2F("hr2","Several graphs in the same pad",2,15,200,2,0,5); | |
175 | c1->cd(2); | |
176 | hr->Draw(); | |
177 | eGraph1 ->Draw("LP"); | |
178 | eGraphR1->Draw("LP"); | |
179 | ||
180 | ||
181 | // 2 < theta < 3 | |
182 | for (j=0; j< 100; j++) { | |
183 | ekin[j] = emin + Float_t(j)*eps; | |
184 | de[j]=0; | |
185 | deRad[j]=0; | |
186 | } | |
187 | // all losses | |
188 | for (i=0; i< nl[1]; i++) { | |
189 | ((TGeant3*) gMC)->Gftmat(pABSO->GetMatId(mat[1][i]), 5, tmp, 100, ekin, dedx, pcut, ixst); | |
190 | for (j=0; j< 100; j++) de[j]+=dedx[j]*(zmax[1][i]-zmin[1][i])/1000.; | |
191 | } | |
192 | // radative losses | |
193 | for (i=0; i< nl[1]; i++) { | |
194 | Float_t a, z, dens, radl, absl; | |
195 | pABSO->AliGetMaterial(mat[1][i], name, a, z, dens, radl, absl); | |
196 | for (Int_t j=0; j<100; j++) { | |
197 | Float_t nor= kAvo*dens/a*(zmax[1][i]-zmin[1][i]); | |
198 | deRad[j] += (((TGeant3*) gMC)->Gbrelm(z,ekin[j],1.e10))*nor; | |
199 | deRad[j] += (((TGeant3*) gMC)->Gprelm(z,ekin[j],1.e10))*nor; | |
200 | } | |
201 | } | |
202 | ||
203 | TH2F *hr2 = new TH2F("hr3","Several graphs in the same pad",2,0,15,2,0,5); | |
204 | hr2->SetXTitle("X title"); | |
205 | hr2->SetYTitle("Y title"); | |
206 | c1->cd(3); | |
207 | hr2->Draw(); | |
208 | ||
209 | TGraph *eGraph2 = new TGraph(100, ekin, de); | |
210 | TGraph *eGraphR2 = new TGraph(100, ekin, deRad); | |
211 | ||
212 | eGraph2 ->Draw("LP"); | |
213 | eGraphR2->Draw("LP"); | |
214 | ||
215 | hr2 = new TH2F("hr4","Several graphs in the same pad",2,15,200,2,0,5); | |
216 | c1->cd(4); | |
217 | hr2->Draw(); | |
218 | eGraph2 ->Draw("LP"); | |
219 | eGraphR2->Draw("LP"); | |
220 | ||
221 | ||
222 | c->Update(); | |
223 | ||
224 | } | |
225 | ||
226 | ||
227 | ||
228 |