| Geant4 Cross Reference |
1 // ------------------------------------------- 1 // -------------------------------------------------------------------
>> 2 // $Id: plot.C 69157 2013-04-19 14:05:11Z gcosmo $
2 // ------------------------------------------- 3 // -------------------------------------------------------------------
3 // 4 //
4 // ******************************************* 5 // *********************************************************************
5 // To execute this macro under ROOT, 6 // To execute this macro under ROOT,
6 // 1 - launch ROOT (usually type 'root' at y 7 // 1 - launch ROOT (usually type 'root' at your machine's prompt)
7 // 2 - type '.X plot.C' at the ROOT session 8 // 2 - type '.X plot.C' at the ROOT session prompt
8 // This macro needs five files : dose.txt, sto 9 // This macro needs five files : dose.txt, stoppingPower.txt, range.txt,
9 // 3DDose.txt and beamPosition.txt 10 // 3DDose.txt and beamPosition.txt
10 // written by S. Incerti and O. Boissonnade, 1 11 // written by S. Incerti and O. Boissonnade, 10/04/2006
11 // ******************************************* 12 // *********************************************************************
12 { 13 {
13 14
14 gROOT->Reset(); 15 gROOT->Reset();
15 16
16 //**************** 17 //****************
17 // DOSE IN NUCLEUS 18 // DOSE IN NUCLEUS
18 //**************** 19 //****************
19 20
20 gStyle->SetOptStat(0000); 21 gStyle->SetOptStat(0000);
21 gStyle->SetOptFit(); 22 gStyle->SetOptFit();
22 gStyle->SetPalette(1); 23 gStyle->SetPalette(1);
23 gROOT->SetStyle("Plain"); 24 gROOT->SetStyle("Plain");
24 Double_t scale; 25 Double_t scale;
25 26
26 27
27 c1 = new TCanvas ("c1","",20,20,1200,900); 28 c1 = new TCanvas ("c1","",20,20,1200,900);
28 c1->Divide(4,3); << 29 c1.Divide(4,3);
29 30
30 //********************* 31 //*********************
31 // INTENSITY HISTOGRAMS 32 // INTENSITY HISTOGRAMS
32 //********************* 33 //*********************
>> 34 jump:
33 35
34 FILE * fp = fopen("phantom.dat","r"); 36 FILE * fp = fopen("phantom.dat","r");
35 Float_t xVox, yVox, zVox, tmp, den, dose; 37 Float_t xVox, yVox, zVox, tmp, den, dose;
36 38
37 Float_t X,Y,Z; 39 Float_t X,Y,Z;
38 Float_t vox = 0, mat = 0; 40 Float_t vox = 0, mat = 0;
39 Float_t voxelSizeX, voxelSizeY, voxelSizeZ; 41 Float_t voxelSizeX, voxelSizeY, voxelSizeZ;
40 42
41 TH1F *h1 = new TH1F("h1","Nucleus marker inte 43 TH1F *h1 = new TH1F("h1","Nucleus marker intensity",100,1,300);
42 TH1F *h11 = new TH1F("h11 ","",100,1,300); 44 TH1F *h11 = new TH1F("h11 ","",100,1,300);
43 45
44 TH1F *h2 = new TH1F("h2","Cytoplasm marker in 46 TH1F *h2 = new TH1F("h2","Cytoplasm marker intensity",100,1,300);
45 TH1F *h20 = new TH1F("h20 ","",100,1,300); 47 TH1F *h20 = new TH1F("h20 ","",100,1,300);
46 48
47 TNtuple *ntupleYXN = new TNtuple("NUCLEUS","nt 49 TNtuple *ntupleYXN = new TNtuple("NUCLEUS","ntuple","Y:X:vox");
48 TNtuple *ntupleZX = new TNtuple("CYTOPLASM","n 50 TNtuple *ntupleZX = new TNtuple("CYTOPLASM","ntuple","Z:X:vox");
49 TNtuple *ntupleYX = new TNtuple("CYTOPLASM","n 51 TNtuple *ntupleYX = new TNtuple("CYTOPLASM","ntuple","Y:X:vox");
50 52
51 Int_t nlines=0; << 53 nlines=0;
52 Int_t ncols=0; << 54 ncols=0;
53 55
54 while (1) 56 while (1)
55 { 57 {
56 if ( nlines == 0 ) ncols = fscanf(fp,"%f 58 if ( nlines == 0 ) ncols = fscanf(fp,"%f %f %f",&tmp,&tmp,&tmp);
57 if ( nlines == 1 ) ncols = fscanf(fp,"%f 59 if ( nlines == 1 ) ncols = fscanf(fp,"%f %f %f",&voxelSizeX,&voxelSizeY,&voxelSizeZ);
58 if ( nlines == 2 ) ncols = fscanf(fp,"%f 60 if ( nlines == 2 ) ncols = fscanf(fp,"%f %f %f",&tmp,&tmp,&tmp);
59 if ( nlines >= 3 ) ncols = fscanf(fp,"%f 61 if ( nlines >= 3 ) ncols = fscanf(fp,"%f %f %f %f %f %f", &X, &Y, &Z, &mat, &den, &vox);
60 if (ncols < 0) break; 62 if (ncols < 0) break;
61 63
62 X= X*voxelSizeX; 64 X= X*voxelSizeX;
63 Y= Y*voxelSizeY; 65 Y= Y*voxelSizeY;
64 Z= Z*voxelSizeZ; 66 Z= Z*voxelSizeZ;
65 67
66 if ( mat == 2 ) // noyau 68 if ( mat == 2 ) // noyau
67 { 69 {
68 if (den==1) h1->Fill( vox ); 70 if (den==1) h1->Fill( vox );
69 if (den==2) h11->Fill( vox ); 71 if (den==2) h11->Fill( vox );
70 ntupleYXN->Fill(Y,X,vox); 72 ntupleYXN->Fill(Y,X,vox);
71 } 73 }
72 if ( mat == 1 ) // cytoplasm 74 if ( mat == 1 ) // cytoplasm
73 { 75 {
74 if (den==1) h2->Fill( vox ); 76 if (den==1) h2->Fill( vox );
75 if (den==2) h20->Fill( vox ); 77 if (den==2) h20->Fill( vox );
76 ntupleZX->Fill(Z,X,vox); 78 ntupleZX->Fill(Z,X,vox);
77 ntupleYX->Fill(Y,X,vox); 79 ntupleYX->Fill(Y,X,vox);
78 } 80 }
79 nlines++; 81 nlines++;
80 82
81 } 83 }
82 fclose(fp); 84 fclose(fp);
83 85
84 // HISTO NUCLEUS 86 // HISTO NUCLEUS
85 87
86 c1->cd(1); << 88 c1.cd(1);
87 h1->Draw(); 89 h1->Draw();
88 h1->GetXaxis()->SetLabelSize(0.025); 90 h1->GetXaxis()->SetLabelSize(0.025);
89 h1->GetYaxis()->SetLabelSize(0.025); 91 h1->GetYaxis()->SetLabelSize(0.025);
90 h1->GetXaxis()->SetTitleSize(0.035); 92 h1->GetXaxis()->SetTitleSize(0.035);
91 h1->GetYaxis()->SetTitleSize(0.035); 93 h1->GetYaxis()->SetTitleSize(0.035);
92 h1->GetXaxis()->SetTitleOffset(1.4); 94 h1->GetXaxis()->SetTitleOffset(1.4);
93 h1->GetYaxis()->SetTitleOffset(1.4); 95 h1->GetYaxis()->SetTitleOffset(1.4);
94 h1->GetXaxis()->SetTitle("Voxel intensity (0 96 h1->GetXaxis()->SetTitle("Voxel intensity (0-255)");
95 h1->GetYaxis()->SetTitle("Number of events") 97 h1->GetYaxis()->SetTitle("Number of events");
96 h1->SetLineColor(3); 98 h1->SetLineColor(3);
97 h1->SetFillColor(3); // green 99 h1->SetFillColor(3); // green
98 100
99 h11->SetLineColor(8); 101 h11->SetLineColor(8);
100 h11->SetFillColor(8); // dark green 102 h11->SetFillColor(8); // dark green
101 h11->Draw("same"); 103 h11->Draw("same");
102 104
103 // HISTO CYTOPLASM 105 // HISTO CYTOPLASM
104 106
105 c1->cd(5); << 107 c1.cd(5);
106 h2->Draw(); 108 h2->Draw();
107 h2->GetXaxis()->SetLabelSize(0.025); 109 h2->GetXaxis()->SetLabelSize(0.025);
108 h2->GetYaxis()->SetLabelSize(0.025); 110 h2->GetYaxis()->SetLabelSize(0.025);
109 h2->GetXaxis()->SetTitleSize(0.035); 111 h2->GetXaxis()->SetTitleSize(0.035);
110 h2->GetYaxis()->SetTitleSize(0.035); 112 h2->GetYaxis()->SetTitleSize(0.035);
111 h2->GetXaxis()->SetTitleOffset(1.4); 113 h2->GetXaxis()->SetTitleOffset(1.4);
112 h2->GetYaxis()->SetTitleOffset(1.4); 114 h2->GetYaxis()->SetTitleOffset(1.4);
113 h2->GetXaxis()->SetTitle("Voxel intensity (0 115 h2->GetXaxis()->SetTitle("Voxel intensity (0-255)");
114 h2->GetYaxis()->SetTitle("Number of events") 116 h2->GetYaxis()->SetTitle("Number of events");
115 h2->SetLineColor(2); 117 h2->SetLineColor(2);
116 h2->SetFillColor(2); // red 118 h2->SetFillColor(2); // red
117 119
118 h20->SetLineColor(5); 120 h20->SetLineColor(5);
119 h20->SetFillColor(5); // yellow (nucleoli) 121 h20->SetFillColor(5); // yellow (nucleoli)
120 h20->Draw("same"); 122 h20->Draw("same");
121 123
122 //************************* 124 //*************************
123 // CUMULATED CELL INTENSITY 125 // CUMULATED CELL INTENSITY
124 //************************* 126 //*************************
125 127
126 gStyle->SetOptStat(0000); 128 gStyle->SetOptStat(0000);
127 gStyle->SetOptFit(); 129 gStyle->SetOptFit();
128 gStyle->SetPalette(1); 130 gStyle->SetPalette(1);
129 gROOT->SetStyle("Plain"); 131 gROOT->SetStyle("Plain");
130 132
131 //CYTOPLASM 133 //CYTOPLASM
132 134
133 c1->cd(7); // axe YX << 135 c1.cd(7); // axe YX
134 TH2F *hist = new TH2F("hist","hist",50,-20,2 136 TH2F *hist = new TH2F("hist","hist",50,-20,20,50,-20,20);
135 ntupleYX->Draw("Y:X>>hist","vox","colz"); 137 ntupleYX->Draw("Y:X>>hist","vox","colz");
136 gPad->SetLogz(); 138 gPad->SetLogz();
137 hist->Draw("colz"); 139 hist->Draw("colz");
138 hist->GetXaxis()->SetLabelSize(0.025); 140 hist->GetXaxis()->SetLabelSize(0.025);
139 hist->GetYaxis()->SetLabelSize(0.025); 141 hist->GetYaxis()->SetLabelSize(0.025);
140 hist->GetZaxis()->SetLabelSize(0.025); 142 hist->GetZaxis()->SetLabelSize(0.025);
141 hist->GetXaxis()->SetTitleSize(0.035); 143 hist->GetXaxis()->SetTitleSize(0.035);
142 hist->GetYaxis()->SetTitleSize(0.035); 144 hist->GetYaxis()->SetTitleSize(0.035);
143 hist->GetXaxis()->SetTitleOffset(1.4); 145 hist->GetXaxis()->SetTitleOffset(1.4);
144 hist->GetYaxis()->SetTitleOffset(1.4); 146 hist->GetYaxis()->SetTitleOffset(1.4);
145 hist->GetXaxis()->SetTitle("Y (um)"); << 147 hist->GetXaxis()->SetTitle("Y (µm)");
146 hist->GetYaxis()->SetTitle("X (um)"); << 148 hist->GetYaxis()->SetTitle("X (µm)");
147 hist->SetTitle("Cytoplasm intensity on trans 149 hist->SetTitle("Cytoplasm intensity on transverse section");
148 150
149 //NUCLEUS 151 //NUCLEUS
150 152
151 c1->cd(3); // axe YX << 153 c1.cd(3); // axe YX
152 TH2F *hist2 = new TH2F("hist2","hist2",50,-2 << 154 TH2F *hist = new TH2F("hist","hist",50,-20,20,50,-20,20);
153 ntupleYXN->Draw("Y:X>>hist2","vox","colz"); << 155 ntupleYXN->Draw("Y:X>>hist","vox","colz");
154 gPad->SetLogz(); 156 gPad->SetLogz();
155 hist2->Draw("colz"); << 157 hist->Draw("colz");
156 hist2->GetXaxis()->SetLabelSize(0.025); << 158 hist->GetXaxis()->SetLabelSize(0.025);
157 hist2->GetYaxis()->SetLabelSize(0.025); << 159 hist->GetYaxis()->SetLabelSize(0.025);
158 hist2->GetZaxis()->SetLabelSize(0.025); << 160 hist->GetZaxis()->SetLabelSize(0.025);
159 hist2->GetXaxis()->SetTitleSize(0.035); << 161 hist->GetXaxis()->SetTitleSize(0.035);
160 hist2->GetYaxis()->SetTitleSize(0.035); << 162 hist->GetYaxis()->SetTitleSize(0.035);
161 hist2->GetXaxis()->SetTitleOffset(1.4); << 163 hist->GetXaxis()->SetTitleOffset(1.4);
162 hist2->GetYaxis()->SetTitleOffset(1.4); << 164 hist->GetYaxis()->SetTitleOffset(1.4);
163 hist2->GetXaxis()->SetTitle("Y (um)"); << 165 hist->GetXaxis()->SetTitle("Y (µm)");
164 hist2->GetYaxis()->SetTitle("X (um)"); << 166 hist->GetYaxis()->SetTitle("X (µm)");
165 hist2->SetTitle("Nucleus intensity on transv << 167 hist->SetTitle("Nucleus intensity on transverse section");
166 168
167 // 169 //
168 170
169 system ("rm -rf microbeam.root"); 171 system ("rm -rf microbeam.root");
170 system ("hadd -O microbeam.root microbeam_*.ro << 172 system ("hadd microbeam.root microbeam_*.root");
171 173
172 TFile f("microbeam.root"); 174 TFile f("microbeam.root");
173 175
174 TNtuple* ntuple0; 176 TNtuple* ntuple0;
175 TNtuple* ntuple1; 177 TNtuple* ntuple1;
176 TNtuple* ntuple2; 178 TNtuple* ntuple2;
177 TNtuple* ntuple3; 179 TNtuple* ntuple3;
178 TNtuple* ntuple4; 180 TNtuple* ntuple4;
179 181
180 ntuple0 = (TNtuple*)f.Get("ntuple0"); << 182 ntuple0 = (TNtuple*)f->Get("ntuple0");
181 ntuple1 = (TNtuple*)f.Get("ntuple1"); << 183 ntuple1 = (TNtuple*)f->Get("ntuple1");
182 ntuple2 = (TNtuple*)f.Get("ntuple2"); << 184 ntuple2 = (TNtuple*)f->Get("ntuple2");
183 ntuple3 = (TNtuple*)f.Get("ntuple3"); << 185 ntuple3 = (TNtuple*)f->Get("ntuple3");
184 ntuple4 = (TNtuple*)f.Get("ntuple4"); << 186 ntuple4 = (TNtuple*)f->Get("ntuple4");
185 << 187
186 TH1F *h1bis = new TH1F("h1bis","Dose distribu << 188 TH1F *h1 = new TH1F("h1","Dose distribution in Nucleus",100,0.001,1.);
187 TH1F *h10 = new TH1F("h10bis","Dose distributi << 189 TH1F *h10 = new TH1F("h10","Dose distribution in Cytoplasm",100,0.001,.2);
188 << 190
189 c1->cd(2); << 191 c1.cd(2);
190 << 192
191 ntuple3->Project("h1bis","doseN"); << 193 ntuple3->Project("h1","doseN");
192 scale = 1/h1bis->Integral(); << 194 scale = 1/h1->Integral();
193 h1bis->Scale(scale); << 195 h1->Scale(scale);
194 h1bis->Draw(); << 196 h1->Draw();
195 h1bis->GetXaxis()->SetLabelSize(0.025); << 197 h1->GetXaxis()->SetLabelSize(0.025);
196 h1bis->GetYaxis()->SetLabelSize(0.025); << 198 h1->GetYaxis()->SetLabelSize(0.025);
197 h1bis->GetXaxis()->SetTitleSize(0.035); << 199 h1->GetXaxis()->SetTitleSize(0.035);
198 h1bis->GetYaxis()->SetTitleSize(0.035); << 200 h1->GetYaxis()->SetTitleSize(0.035);
199 h1bis->GetXaxis()->SetTitleOffset(1.4); << 201 h1->GetXaxis()->SetTitleOffset(1.4);
200 h1bis->GetYaxis()->SetTitleOffset(1.4); << 202 h1->GetYaxis()->SetTitleOffset(1.4);
201 h1bis->GetXaxis()->SetTitle("Absorbed dose ( << 203 h1->GetXaxis()->SetTitle("Absorbed dose (Gy)");
202 h1bis->GetYaxis()->SetTitle("Fraction of eve << 204 h1->GetYaxis()->SetTitle("Fraction of events");
203 h1bis->SetLineColor(3); << 205 h1->SetLineColor(3);
204 h1bis->SetFillColor(3); << 206 h1->SetFillColor(3);
205 207
206 //***************** 208 //*****************
207 // DOSE IN CYTOPLASM 209 // DOSE IN CYTOPLASM
208 //***************** 210 //*****************
209 211
210 c1->cd(6); << 212 c1.cd(6);
211 ntuple3->Project("h10bis","doseC"); << 213 ntuple3->Project("h10","doseC");
212 scale = 1/h10->Integral(); 214 scale = 1/h10->Integral();
213 h10->Scale(scale); 215 h10->Scale(scale);
214 h10->Draw(); 216 h10->Draw();
215 h10->GetXaxis()->SetLabelSize(0.025); 217 h10->GetXaxis()->SetLabelSize(0.025);
216 h10->GetYaxis()->SetLabelSize(0.025); 218 h10->GetYaxis()->SetLabelSize(0.025);
217 h10->GetXaxis()->SetTitleSize(0.035); 219 h10->GetXaxis()->SetTitleSize(0.035);
218 h10->GetYaxis()->SetTitleSize(0.035); 220 h10->GetYaxis()->SetTitleSize(0.035);
219 h10->GetXaxis()->SetTitleOffset(1.4); 221 h10->GetXaxis()->SetTitleOffset(1.4);
220 h10->GetYaxis()->SetTitleOffset(1.4); 222 h10->GetYaxis()->SetTitleOffset(1.4);
221 h10->GetXaxis()->SetTitle("Absorbed dose (Gy 223 h10->GetXaxis()->SetTitle("Absorbed dose (Gy)");
222 h10->GetYaxis()->SetTitle("Fraction of event 224 h10->GetYaxis()->SetTitle("Fraction of events");
223 h10->SetLineColor(2); 225 h10->SetLineColor(2);
224 h10->SetFillColor(2); 226 h10->SetFillColor(2);
225 227
226 //******************************** 228 //********************************
227 // STOPPING POWER AT CELL ENTRANCE 229 // STOPPING POWER AT CELL ENTRANCE
228 //******************************** 230 //********************************
229 231
230 gStyle->SetOptStat(0000); 232 gStyle->SetOptStat(0000);
231 gStyle->SetOptFit(); 233 gStyle->SetOptFit();
232 gStyle->SetPalette(1); 234 gStyle->SetPalette(1);
233 gROOT->SetStyle("Plain"); 235 gROOT->SetStyle("Plain");
234 236
235 Float_t d; 237 Float_t d;
236 238
237 TH1F *h2bis = new TH1F("h2bis","Beam stopping << 239 TH1F *h2 = new TH1F("h2","Beam stopping power at cell entrance",200,0,300);
238 240
239 c1->cd(9); << 241 c1.cd(9);
240 ntuple0->Project("h2bis","sp"); << 242 ntuple0->Project("h2","sp");
241 scale = 1/h2bis->Integral(); << 243 scale = 1/h2->Integral();
242 h2bis->Scale(scale); << 244 h2->Scale(scale);
243 h2bis->Draw(); << 245 h2->Draw();
244 h2bis->GetXaxis()->SetLabelSize(0.025); << 246 h2->GetXaxis()->SetLabelSize(0.025);
245 h2bis->GetYaxis()->SetLabelSize(0.025); << 247 h2->GetYaxis()->SetLabelSize(0.025);
246 h2bis->GetXaxis()->SetTitleSize(0.035); << 248 h2->GetXaxis()->SetTitleSize(0.035);
247 h2bis->GetYaxis()->SetTitleSize(0.035); << 249 h2->GetYaxis()->SetTitleSize(0.035);
248 h2bis->GetXaxis()->SetTitleOffset(1.4); << 250 h2->GetXaxis()->SetTitleOffset(1.4);
249 h2bis->GetYaxis()->SetTitleOffset(1.4); << 251 h2->GetYaxis()->SetTitleOffset(1.4);
250 h2bis->GetXaxis()->SetTitle("dE/dx (keV/um << 252 h2->GetXaxis()->SetTitle("dE/dx (keV/µm)");
251 h2bis->GetYaxis()->SetTitle("Fraction of eve << 253 h2->GetYaxis()->SetTitle("Fraction of events");
252 h2bis->SetTitle("dE/dx at cell entrance"); << 254 h2->SetTitle("dE/dx at cell entrance");
253 h2bis->SetFillColor(4); << 255 h2->SetFillColor(4);
254 h2bis->SetLineColor(4); << 256 h2->SetLineColor(4);
255 h2bis->Fit("gaus"); << 257 h2->Fit("gaus");
256 gaus->SetLineColor(6); 258 gaus->SetLineColor(6);
257 h2bis->Fit("gaus"); << 259 h2->Fit("gaus");
>> 260
258 261
259 //************** 262 //**************
260 // RANGE IN CELL 263 // RANGE IN CELL
261 //************** 264 //**************
262 265
263 Double_t Xc,Zc,X1,Y1,Z1,X2,Y2,Z2; 266 Double_t Xc,Zc,X1,Y1,Z1,X2,Y2,Z2;
264 267
265 // X position of target in World 268 // X position of target in World
266 Xc = -1295.59e3 - 955e3*sin(10*TMath::Pi()/180 269 Xc = -1295.59e3 - 955e3*sin(10*TMath::Pi()/180);
267 270
268 // Z position of target in World 271 // Z position of target in World
269 Zc = -1327e3 + 955e3*cos(10*TMath::Pi()/180); 272 Zc = -1327e3 + 955e3*cos(10*TMath::Pi()/180);
270 273
271 // Line alignment (cf MicrobeamEMField.cc) 274 // Line alignment (cf MicrobeamEMField.cc)
272 Xc = Xc + 5.24*cos(10*TMath::Pi()/180); 275 Xc = Xc + 5.24*cos(10*TMath::Pi()/180);
273 Zc = Zc + 5.24*sin(10*TMath::Pi()/180); 276 Zc = Zc + 5.24*sin(10*TMath::Pi()/180);
274 277
275 TNtuple *ntupleR = new TNtuple("Rmax","ntuple" 278 TNtuple *ntupleR = new TNtuple("Rmax","ntuple","Z2:Y2:X2");
276 Double_t x,y,z,xx,zz; 279 Double_t x,y,z,xx,zz;
277 ntuple2->SetBranchAddress("x",&x); 280 ntuple2->SetBranchAddress("x",&x);
278 ntuple2->SetBranchAddress("y",&y); 281 ntuple2->SetBranchAddress("y",&y);
279 ntuple2->SetBranchAddress("z",&z); 282 ntuple2->SetBranchAddress("z",&z);
280 Int_t nentries = (Int_t)ntuple2->GetEntries(); 283 Int_t nentries = (Int_t)ntuple2->GetEntries();
281 for (Int_t i=0;i<nentries;i++) 284 for (Int_t i=0;i<nentries;i++)
282 { 285 {
283 ntuple2->GetEntry(i); 286 ntuple2->GetEntry(i);
284 X1=x; 287 X1=x;
285 Y1=y; 288 Y1=y;
286 Z1=z; 289 Z1=z;
287 xx = X1-Xc; 290 xx = X1-Xc;
288 zz = Z1-Zc; 291 zz = Z1-Zc;
289 Z2 = zz*cos(10*TMath::Pi()/180)-xx*sin(1 292 Z2 = zz*cos(10*TMath::Pi()/180)-xx*sin(10*TMath::Pi()/180);
290 X2 = zz*sin(10*TMath::Pi()/180)+xx*cos(1 293 X2 = zz*sin(10*TMath::Pi()/180)+xx*cos(10*TMath::Pi()/180);
291 Y2 = Y1; 294 Y2 = Y1;
292 ntupleR->Fill(Z2,Y2,X2); 295 ntupleR->Fill(Z2,Y2,X2);
293 } 296 }
294 297
295 c1->cd(10); << 298 c1.cd(10);
296 ntupleR->Draw("X2:Z2","abs(X2)<50","surf3"); 299 ntupleR->Draw("X2:Z2","abs(X2)<50","surf3");
297 gPad->SetLogz(); 300 gPad->SetLogz();
>> 301 /*
>> 302 htemp->GetXaxis()->SetLabelSize(0.025);
>> 303 htemp->GetYaxis()->SetLabelSize(0.025);
>> 304 htemp->GetZaxis()->SetLabelSize(0.025);
>> 305 htemp->GetXaxis()->SetTitleSize(0.035);
>> 306 htemp->GetYaxis()->SetTitleSize(0.035);
>> 307 htemp->GetXaxis()->SetTitleOffset(1.4);
>> 308 htemp->GetYaxis()->SetTitleOffset(1.4);
>> 309 htemp->GetXaxis()->SetTitle("Z (µm)");
>> 310 htemp->GetYaxis()->SetTitle("X (µm)");
>> 311 htemp->SetTitle("Range in cell");
>> 312 */
298 313
299 //**************** 314 //****************
300 // ENERGY DEPOSITS 315 // ENERGY DEPOSITS
301 //**************** 316 //****************
302 317
>> 318
303 gStyle->SetOptStat(0000); 319 gStyle->SetOptStat(0000);
304 gStyle->SetOptFit(); 320 gStyle->SetOptFit();
305 gStyle->SetPalette(1); 321 gStyle->SetPalette(1);
306 gROOT->SetStyle("Plain"); 322 gROOT->SetStyle("Plain");
307 323
308 c1->cd(11); << 324 c1.cd(11);
309 TH2F *histbis = new TH2F("histbis","histbis" << 325 TH2F *hist = new TH2F("hist","hist",50,-20,20,50,-20,20);
310 ntuple4->Draw("y*0.359060:x*0.359060>>histbi << 326 ntuple4->Draw("y*0.359060:x*0.359060>>hist","doseV","contz");
311 gPad->SetLogz(); 327 gPad->SetLogz();
312 histbis->Draw("contz"); << 328 hist->Draw("contz");
313 histbis->GetXaxis()->SetLabelSize(0.025); << 329 hist->GetXaxis()->SetLabelSize(0.025);
314 histbis->GetYaxis()->SetLabelSize(0.025); << 330 hist->GetYaxis()->SetLabelSize(0.025);
315 histbis->GetZaxis()->SetLabelSize(0.025); << 331 hist->GetZaxis()->SetLabelSize(0.025);
316 histbis->GetXaxis()->SetTitleSize(0.035); << 332 hist->GetXaxis()->SetTitleSize(0.035);
317 histbis->GetYaxis()->SetTitleSize(0.035); << 333 hist->GetYaxis()->SetTitleSize(0.035);
318 histbis->GetXaxis()->SetTitleOffset(1.4); << 334 hist->GetXaxis()->SetTitleOffset(1.4);
319 histbis->GetYaxis()->SetTitleOffset(1.4); << 335 hist->GetYaxis()->SetTitleOffset(1.4);
320 histbis->GetXaxis()->SetTitle("Y (um)"); << 336 hist->GetXaxis()->SetTitle("Y (µm)");
321 histbis->GetYaxis()->SetTitle("X (um)"); << 337 hist->GetYaxis()->SetTitle("X (µm)");
322 histbis->SetTitle("Energy deposit -transvers << 338 hist->SetTitle("Mean energy deposit -transverse- (z axis in eV)");
323 << 339
324 c1->cd(12); << 340
325 TH2F *histter = new TH2F("histter","histter" << 341 c1.cd(12);
326 ntuple4->Draw("x*0.359060:(z+1500/0.162810+2 << 342 TH2F *hist = new TH2F("hist","hist",50,-20,20,50,-20,20);
>> 343 ntuple4->Draw("x*0.359060:z*0.162810>>hist","doseV","contz");
327 gPad->SetLogz(); 344 gPad->SetLogz();
328 histter->Draw("contz"); << 345 hist->Draw("contz");
329 histter->GetXaxis()->SetLabelSize(0.025); << 346 hist->GetXaxis()->SetLabelSize(0.025);
330 histter->GetYaxis()->SetLabelSize(0.025); << 347 hist->GetYaxis()->SetLabelSize(0.025);
331 histter->GetZaxis()->SetLabelSize(0.025); << 348 hist->GetZaxis()->SetLabelSize(0.025);
332 histter->GetXaxis()->SetTitleSize(0.035); << 349 hist->GetXaxis()->SetTitleSize(0.035);
333 histter->GetYaxis()->SetTitleSize(0.035); << 350 hist->GetYaxis()->SetTitleSize(0.035);
334 histter->GetXaxis()->SetTitleOffset(1.4); << 351 hist->GetXaxis()->SetTitleOffset(1.4);
335 histter->GetYaxis()->SetTitleOffset(1.4); << 352 hist->GetYaxis()->SetTitleOffset(1.4);
336 histter->GetXaxis()->SetTitle("Z (um)"); << 353 hist->GetXaxis()->SetTitle("Z (µm)");
337 histter->GetYaxis()->SetTitle("X (um)"); << 354 hist->GetYaxis()->SetTitle("X (µm)");
338 histter->SetTitle("Energy deposit -longitudi << 355 hist->SetTitle("Mean energy deposit -longitudinal- (z axis in eV)");
339 356
340 //******************************* 357 //*******************************
341 // BEAM POSITION AT CELL ENTRANCE 358 // BEAM POSITION AT CELL ENTRANCE
342 //******************************* 359 //*******************************
343 360
344 gStyle->SetOptStat(0000); 361 gStyle->SetOptStat(0000);
345 gStyle->SetOptFit(); 362 gStyle->SetOptFit();
346 gStyle->SetPalette(1); 363 gStyle->SetPalette(1);
347 gROOT->SetStyle("Plain"); 364 gROOT->SetStyle("Plain");
348 365
349 TH1F *h77 = new TH1F("hx","h1",200,-10,10); 366 TH1F *h77 = new TH1F("hx","h1",200,-10,10);
350 TH1F *h88 = new TH1F("hy","h1",200,-10,10); 367 TH1F *h88 = new TH1F("hy","h1",200,-10,10);
351 368
352 c1->cd(4); << 369 c1.cd(4);
353 ntuple1->Project("hx","x"); 370 ntuple1->Project("hx","x");
354 scale = 1/h77->Integral(); 371 scale = 1/h77->Integral();
355 h77->Scale(scale); 372 h77->Scale(scale);
356 h77->Draw(); 373 h77->Draw();
357 h77->GetXaxis()->SetLabelSize(0.025); 374 h77->GetXaxis()->SetLabelSize(0.025);
358 h77->GetYaxis()->SetLabelSize(0.025); 375 h77->GetYaxis()->SetLabelSize(0.025);
359 h77->GetXaxis()->SetTitleSize(0.035); 376 h77->GetXaxis()->SetTitleSize(0.035);
360 h77->GetYaxis()->SetTitleSize(0.035); 377 h77->GetYaxis()->SetTitleSize(0.035);
361 h77->GetXaxis()->SetTitleOffset(1.4); 378 h77->GetXaxis()->SetTitleOffset(1.4);
362 h77->GetYaxis()->SetTitleOffset(1.4); 379 h77->GetYaxis()->SetTitleOffset(1.4);
363 h77->GetXaxis()->SetTitle("Position (um)") << 380 h77->GetXaxis()->SetTitle("Position (µm)");
364 h77->GetYaxis()->SetTitle("Fraction of event 381 h77->GetYaxis()->SetTitle("Fraction of events");
365 h77->SetTitle("Beam X position on cell"); 382 h77->SetTitle("Beam X position on cell");
366 h77->SetFillColor(4); 383 h77->SetFillColor(4);
367 h77->SetLineColor(4); 384 h77->SetLineColor(4);
>> 385 gaus->SetLineColor(6);
368 h77->Fit("gaus"); 386 h77->Fit("gaus");
369 387
370 c1->cd(8); << 388 c1.cd(8);
371 ntuple1->Project("hy","y"); 389 ntuple1->Project("hy","y");
372 scale = 1/h88->Integral(); 390 scale = 1/h88->Integral();
373 h88->Scale(scale); 391 h88->Scale(scale);
374 h88->Draw(); 392 h88->Draw();
375 h88->GetXaxis()->SetLabelSize(0.025); 393 h88->GetXaxis()->SetLabelSize(0.025);
376 h88->GetYaxis()->SetLabelSize(0.025); 394 h88->GetYaxis()->SetLabelSize(0.025);
377 h88->GetXaxis()->SetTitleSize(0.035); 395 h88->GetXaxis()->SetTitleSize(0.035);
378 h88->GetYaxis()->SetTitleSize(0.035); 396 h88->GetYaxis()->SetTitleSize(0.035);
379 h88->GetXaxis()->SetTitleOffset(1.4); 397 h88->GetXaxis()->SetTitleOffset(1.4);
380 h88->GetYaxis()->SetTitleOffset(1.4); 398 h88->GetYaxis()->SetTitleOffset(1.4);
381 h88->GetXaxis()->SetTitle("Position (um)") << 399 h88->GetXaxis()->SetTitle("Position (µm)");
382 h88->GetYaxis()->SetTitle("Fraction of event 400 h88->GetYaxis()->SetTitle("Fraction of events");
383 h88->SetTitle("Beam Y position on cell"); 401 h88->SetTitle("Beam Y position on cell");
384 h88->SetFillColor(4); 402 h88->SetFillColor(4);
385 h88->SetLineColor(4); 403 h88->SetLineColor(4);
>> 404 gaus->SetLineColor(6);
386 h88->Fit("gaus"); 405 h88->Fit("gaus");
>> 406
387 } 407 }
388 408