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};
8 AliABSO *pABSO = (AliABSO*) gAlice->GetModule("ABSO");
11 nl[0] = pABSO->NumberOfLayers(0);
12 nl[1] = pABSO->NumberOfLayers(1);
14 Int_t mat[2][15], mat[2][15];
15 Float_t zmin[2][15], zmin[2][15], zmax[2][15], zmax[2][15];
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);
26 Float_t l = zmax[0][nl[0]-1];
32 // 1. Limiting angular resolution in the Branson formalism
35 Float_t a, z, dens, radl, absl;
38 for (j=0; j< 2; j++) {
39 printf("\n A Z dens radL absL");
40 printf("\n==================================================================================");
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];
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);
48 printf("\n %3d %14s %12.3f %12.3f %12.3f %12.3f %12.3f",
49 i+1, name, a, z, dens, radl, absl);
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.);
59 // 2. Limiting Resolution as a function of thickness of extra Fe
61 Float_t zmin0[15], zmax0[15];
62 Float_t x[100], y[100];
64 for (i=0; i< nl[0]; i++) {
65 zmin0[i] = zmin[0][i];
66 zmax0[i] = zmax[0][i];
69 Float_t ds = (zmax[0][3]-zmin[0][2])/100.;
72 for (j=0; j< 100; j++) {
74 Float_t zFe = zmin[0][2]+Float_t(j)*ds;
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];
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);
85 Float_t deltaThetaB = 0.0136 * TMath::Sqrt(f0-f1*f1/f2);
86 x[j] = zmax[0][3]-zFe;
87 y[j] = deltaThetaB*1000.;
89 TGraph *sGraph = new TGraph(100, x, y);
90 TCanvas *c = new TCanvas("c"," ",400,10,600,700);
96 // 3. Limiting Resolution as a function of density of Carbon
98 Float_t dRho = (2.5-1.5)/100.;
100 for (j=0; j< 100; j++) {
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];
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);
112 deltaThetaB = 0.0136 * TMath::Sqrt(f0-f1*f1/f2);
114 y[j] = deltaThetaB*1000.;
116 TGraph *dGraph = new TGraph(100, x, y);
125 TCanvas *c1 = new TCanvas("c1"," ",400,10,600,700);
128 const Float_t kAvo=0.60221367;
132 strncpy(tmp, "LOSS", 4);
136 Float_t ekin[100], dedx[100], de[100], pcut[100], deRad[100];
139 Float_t eps = (emax-emin)/100.;
141 for (j=0; j< 100; j++) {
142 ekin[j] = emin + Float_t(j)*eps;
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.;
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;
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");
168 TGraph *eGraph1 = new TGraph(100, ekin, de);
169 TGraph *eGraphR1 = new TGraph(100, ekin, deRad);
171 eGraph1 ->Draw("LP");
172 eGraphR1->Draw("LP");
174 hr = new TH2F("hr2","Several graphs in the same pad",2,15,200,2,0,5);
177 eGraph1 ->Draw("LP");
178 eGraphR1->Draw("LP");
182 for (j=0; j< 100; j++) {
183 ekin[j] = emin + Float_t(j)*eps;
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.;
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;
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");
209 TGraph *eGraph2 = new TGraph(100, ekin, de);
210 TGraph *eGraphR2 = new TGraph(100, ekin, deRad);
212 eGraph2 ->Draw("LP");
213 eGraphR2->Draw("LP");
215 hr2 = new TH2F("hr4","Several graphs in the same pad",2,15,200,2,0,5);
218 eGraph2 ->Draw("LP");
219 eGraphR2->Draw("LP");