| Geant4 Cross Reference |
1 {
2 // Read reference data in dolan.txt
3 FILE * fg1=fopen("dolan.txt", "r");
4 Int_t n_points_dolan = 8;
5 Float_t x1[n_points_dolan], y1[n_points_dolan];
6 Float_t x, y;
7 Int_t ncols_dolan;
8 Int_t nlines1 = 0;
9
10 while(1)
11 {
12 ncols_dolan = fscanf(fg1, "%f %f", &x, &y);
13 if (ncols_dolan<0) break;
14 x1[nlines1]= x;
15 y1[nlines1] = y;
16 nlines1++;
17 }
18
19 fclose(fg1);
20
21 FILE *fg2=fopen("geant4_dose_Oncura_livermore.txt", "r");
22 Int_t n_points_geant4 = 398;
23 Float_t x2[n_points_geant4], y2[n_points_geant4], ratio_liv[n_points_dolan];
24 Int_t ncols_geant4;
25 Int_t nlines2 = 0;
26
27 while(1)
28 {
29 ncols_geant4 = fscanf(fg2, "%f %f", &x, &y);
30 if (ncols_geant4<0) break;
31 x2[nlines2] = x;
32 y2[nlines2] = y;
33
34 for (int i=0; i<n_points_dolan; i++)
35 {
36 if (x1[i]==x2[nlines2])
37 {
38 ratio_liv[i]= y2[nlines2]/y1[i];
39 // std::cout << "dolan: " << x1[i] << "," << y1[i]
40 // << ", livermore: "<< x2[nlines2] << "," << y2[nlines2]
41 // << " ratio:" << ratio_liv[i] << std::endl;
42 }
43 }
44 nlines2++;
45 }
46
47 fclose(fg2);
48
49 FILE *fg3=fopen("geant4_dose_Oncura_penelope.txt", "r");
50 Float_t x3[n_points_geant4], y3[n_points_geant4], ratio_pen[n_points_dolan];
51 Int_t ncols_geant4_penelope;
52 Int_t nlines3 =0;
53
54 while(1)
55 {
56 ncols_geant4_penelope = fscanf(fg3,"%f %f",&x, &y);
57 if (ncols_geant4_penelope<0) break;
58 // std::cout << "x " << x << std::endl;
59 x3[nlines3]=x;
60 y3[nlines3]=y;
61
62 for (int i=0; i<n_points_dolan; i++)
63 {
64 if (x1[i]==x3[nlines3])
65 {
66 ratio_pen[i]= y3[nlines3]/y1[i];
67 // std::cout << "dolan: " << x1[i] << "," << y1[i]
68 // << ", penelope: "<< x3[nlines3] << "," << y3[nlines3]
69 // << " ratio:" << ratio_pen[i] << std::endl;
70 }
71 }
72
73 nlines3++;
74 }
75
76 fclose(fg3);
77
78 FILE *fg4=fopen("geant4_dose_Oncura_opt0.txt", "r");
79 Float_t x4[n_points_geant4], y4[n_points_geant4], ratio_opt0[n_points_dolan];
80 Int_t ncols_geant4_opt0;
81 Int_t nlines4 =0;
82
83 while(1)
84 {
85 ncols_geant4_opt0 = fscanf(fg4,"%f %f",&x, &y);
86 if (ncols_geant4_opt0<0) break;
87 // std::cout << "x " << x << std::endl;
88 x4[nlines4]=x;
89 y4[nlines4]=y;
90 for (int i=0; i<n_points_dolan; i++)
91 {
92 if (x1[i]==x4[nlines4])
93 {
94 ratio_opt0[i]= y4[nlines4]/y1[i];
95 // std::cout << "dolan: " << x1[i] << "," << y1[i]
96 // << ", opt0: "<< x4[nlines4] << "," << y4[nlines4]
97 // << " ratio:" << ratio_opt0[i] << std::endl;
98 }
99 }
100
101 nlines4++;
102 }
103
104 fclose(fg4);
105
106 FILE *fg5=fopen("geant4_dose_Oncura_opt3.txt", "r");
107 Float_t x5[n_points_geant4], y5[n_points_geant4], ratio_opt3[n_points_dolan];
108 Int_t ncols_geant4_opt3;
109 Int_t nlines5 =0;
110
111 while(1)
112 {
113 ncols_geant4_opt3 = fscanf(fg5,"%f %f",&x, &y);
114 if (ncols_geant4_opt3<0) break;
115 // std::cout << "x " << x << std::endl;
116 x5[nlines5]=x;
117 y5[nlines5]=y;
118
119 for (int i=0; i<n_points_dolan; i++)
120 {
121 if (x1[i]==x5[nlines5])
122 {
123 ratio_opt3[i]= y5[nlines5]/y1[i];
124 std::cout << "dolan: " << x1[i] << "," << y1[i]
125 << ", opt3: "<< x5[nlines5] << "," << y5[nlines5]
126 << " ratio:" << ratio_opt3[i] << std::endl;
127 }
128 }
129 nlines5++;
130 }
131
132 fclose(fg5);
133
134 FILE *fg6=fopen("geant4_dose_Oncura_opt4.txt", "r");
135 Float_t x6[n_points_geant4], y6[n_points_geant4], ratio_opt4[n_points_dolan];
136 Int_t ncols_geant4_opt4;
137 Int_t nlines6 =0;
138
139 while(1)
140 {
141 ncols_geant4_opt4 = fscanf(fg6,"%f %f",&x, &y);
142 if (ncols_geant4_opt4<0) break;
143 // std::cout << "x " << x << std::endl;
144 x6[nlines6]=x;
145 y6[nlines6]=y;
146
147 for (int i=0; i<n_points_dolan; i++)
148 {
149 if (x1[i]==x6[nlines6])
150 {
151 ratio_opt4[i]= y6[nlines6]/y1[i];
152 // std::cout << "dolan: " << x1[i] << "," << y1[i]
153 // << ", opt4: "<< x6[nlines6] << "," << y6[nlines6]
154 // << " ratio:" << ratio_opt4[i] << std::endl;
155 }
156 }
157 nlines6++;
158 }
159
160 fclose(fg6);
161
162 TGraph *gr1 = new TGraph (nlines1, x1, y1);
163 TGraph *gr2 = new TGraph (nlines2, x2, y2);
164 TGraph *gr3 = new TGraph (nlines3, x3, y3);
165 TGraph *gr4 = new TGraph (nlines4, x4, y4);
166 TGraph *gr5 = new TGraph (nlines5, x5, y5);
167 TGraph *gr6 = new TGraph (nlines6, x6, y6);
168 std::cout<< "Livermore" << std::endl;
169
170 for (Int_t j=0; j < nlines1; j++)
171 {
172 std::cout << x1[j] << ", " << ratio_liv[j] << std::endl;
173 }
174 std::cout<< "penelope" << std::endl;
175 for (Int_t j=0; j < nlines1; j++)
176 {
177 std::cout << x1[j] << ", " << ratio_pen[j] << std::endl;
178 }
179
180 std::cout<< "opt0" << std::endl;
181 for (Int_t j=0; j < nlines1; j++)
182 {
183 std::cout << x1[j] << ", " << ratio_opt0[j] << std::endl;
184 }
185
186 std::cout<< "opt3" << std::endl;
187 for (Int_t j=0; j < nlines1; j++)
188 {
189 std::cout << x1[j] << ", " << ratio_opt3[j] << std::endl;
190 }
191 std::cout<< "opt4" << std::endl;
192 for (Int_t j=0; j < nlines1; j++)
193 {
194 std::cout << x1[j] << ", " << ratio_opt4[j] << std::endl;
195 }
196
197 TCanvas *c1 = new TCanvas("c1", "Graph Draw Options", 200, 10, 600, 400);
198
199 gPad->SetLogy();
200
201 // Draw the graph with axis, continuous line, and put a '*' at each point
202
203 gr1->SetTitle("Dose rate distribution");
204 gr1->GetXaxis()->SetTitle("Distance from the centre (cm)");
205 gr1->GetYaxis()->SetTitle("Normalised dose rate distribution");
206 gr1->SetLineWidth(1);
207 gr1->SetMarkerColor(1);
208 gr1->SetMarkerStyle(20);
209 gr1->Draw("AP");
210
211 gr2->SetLineWidth(1);
212 gr2->SetMarkerColor(2);
213 gr2->SetMarkerStyle(21);
214 gr2->SetMarkerSize(0.5);
215 gr2->SetLineColor(2);
216 gr2->Draw("CP");
217
218 gr3->SetLineWidth(0.3);
219 gr3->SetMarkerColor(3);
220 gr3->SetMarkerStyle(21);
221 gr3->SetMarkerSize(0.2);
222 gr3->SetLineColor(3);
223 gr3->Draw("CP");
224
225 gr4->SetLineWidth(0.3);
226 gr4->SetMarkerColor(4);
227 gr4->SetMarkerStyle(21);
228 gr4->SetMarkerSize(0.2);
229 gr4->SetLineColor(4);
230 gr4->Draw("CP");
231
232 gr5->SetLineWidth(0.3);
233 gr5->SetMarkerColor(6);
234 gr5->SetMarkerStyle(21);
235 gr5->SetMarkerSize(0.2);
236 gr5->SetLineColor(6);
237 gr5->Draw("CP");
238
239 gr6->SetLineWidth(0.3);
240 gr6->SetMarkerColor(8);
241 gr6->SetMarkerStyle(21);
242 gr6->SetMarkerSize(0.2);
243 gr6->SetLineColor(8);
244 gr6->Draw("CP");
245
246 TLegend *leg = new TLegend(0.3, 0.5, 0.6, 0.8);
247 leg->SetFillColor(0);
248 leg->AddEntry(gr1, "Reference data", "lp");
249 leg->AddEntry(gr2, "Geant4 - Oncura - Livermore", "lp");
250 leg->AddEntry(gr3, "Geant4 - Penelope", "lp");
251 leg->AddEntry(gr4, "Geant4 - Standard opt0", "lp");
252 leg->AddEntry(gr5, "Geant4 - Standard opt3", "lp");
253 leg->AddEntry(gr6, "Geant4 - Standard opt4", "lp");
254 leg->Draw();
255
256 }
257
258
259