Rev 38 | Rev 47 | Go to most recent revision | Only display areas with differences | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed
Rev 38 | Rev 40 | ||
---|---|---|---|
1 | #include "TROOT.h" |
1 | #include "TROOT.h" |
2 | #include "TFile.h" |
2 | #include "TFile.h" |
3 | #include "TBenchmark.h" |
3 | #include "TBenchmark.h" |
4 | #include "TH1F.h" |
4 | #include "TH1F.h" |
5 | #include "TH2F.h" |
5 | #include "TH2F.h" |
6 | #include "TCanvas.h" |
6 | #include "TCanvas.h" |
7 | #include "TStyle.h" |
7 | #include "TStyle.h" |
8 | #include "TPad.h" |
8 | #include "TPad.h" |
9 | #include "TF1.h" |
9 | #include "TF1.h" |
10 | #include "TGraph.h" |
10 | #include "TGraph.h" |
11 | #include "TSpectrum.h" |
11 | #include "TSpectrum.h" |
12 | #include "stdio.h" |
12 | #include "stdio.h" |
13 | 13 | ||
14 | #include "include/RTUtil.h" |
14 | #include "include/RTUtil.h" |
15 | 15 | ||
16 | double getNoise(TH2F*, int, int); |
16 | double getNoise(TH2F*, int, int); |
17 | 17 | ||
18 | int sipm(char filename[256] = "test", char plopt[256]="all", int chXstart=0, int chXend=7, int chYstart=0, int chYend=7, bool debug = false) |
18 | int sipm(char filename[256] = "test", char plopt[256]="all", int chXstart=0, int chXend=7, int chYstart=0, int chYend=7, bool debug = false) |
19 | { |
19 | { |
20 | const int c_nChannels = 64; |
20 | const int c_nChannels = 64; |
21 | const double c_xOffset = 0; // mm |
21 | const double c_xOffset = 0; // mm |
22 | const double c_yOffset = 0; |
22 | const double c_yOffset = 0; |
23 | 23 | ||
24 | int map[8][8]={{32,34,53,55,40,42,61,63}, |
24 | int map[8][8]={{32,34,53,55,40,42,61,63}, |
25 | {48,50,37,39,56,58,45,47}, |
25 | {48,50,37,39,56,58,45,47}, |
26 | {33,35,52,54,41,43,60,62}, |
26 | {33,35,52,54,41,43,60,62}, |
27 | {49,51,36,38,57,59,44,46}, |
27 | {49,51,36,38,57,59,44,46}, |
28 | {17,19,4,6,25,27,12,14}, |
28 | {17,19,4,6,25,27,12,14}, |
29 | {1,3,20,22,9,11,28,30}, |
29 | {1,3,20,22,9,11,28,30}, |
30 | {16,18,5,7,24,26,13,15}, |
30 | {16,18,5,7,24,26,13,15}, |
31 | {0,2,21,23,8,10,29,31} |
31 | {0,2,21,23,8,10,29,31} |
32 | }; |
32 | }; |
33 | 33 | ||
34 | char fnameroot[256]; |
34 | char fnameroot[256]; |
35 | TFile* rootfile; |
35 | TFile* rootfile; |
36 | sprintf(fnameroot, "root/%s.root", filename); |
36 | sprintf(fnameroot, "root/%s.root", filename); |
37 | rootfile = (TFile *) gROOT->FindObject(filename); |
37 | rootfile = (TFile *) gROOT->FindObject(filename); |
38 | if(rootfile==NULL) rootfile = new TFile(fnameroot); |
38 | if(rootfile==NULL) rootfile = new TFile(fnameroot); |
39 | if(rootfile==NULL) { |
39 | if(rootfile==NULL) { |
40 | printf("Cannot open root file %s!!!\n",fnameroot); |
40 | printf("Cannot open root file %s!!!\n",fnameroot); |
41 | return(0); |
41 | return(0); |
42 | } |
42 | } |
43 | 43 | ||
44 | // set draw style |
44 | // set draw style |
45 | RTSetStyle(gStyle); |
45 | RTSetStyle(gStyle); |
46 | 46 | ||
47 | if( strstr(plopt, "all") != NULL ) { |
47 | if( strstr(plopt, "all") != NULL ) { |
48 | TCanvas *canvas2 = new TCanvas("canvas2","Hits x;;",2000,2000); |
48 | TCanvas *canvas2 = new TCanvas("canvas2","Hits x;;",2000,2000); |
49 | TCanvas *canvas3 = new TCanvas("canvas3","Hits y;;",2000,2000); |
49 | TCanvas *canvas3 = new TCanvas("canvas3","Hits y;;",2000,2000); |
50 | canvas2->Divide(8,8); |
50 | canvas2->Divide(8,8); |
51 | canvas3->Divide(8,8); |
51 | canvas3->Divide(8,8); |
52 | TH1F* h_hitsx; |
52 | TH1F* h_hitsx; |
53 | TH1F* h_hitsy; |
53 | TH1F* h_hitsy; |
54 | for(int i=0; i<c_nChannels; i++) { |
54 | for(int i=0; i<c_nChannels; i++) { |
55 | canvas2->cd(i+1); |
55 | canvas2->cd(i+1); |
56 | char hname[128]; |
56 | char hname[128]; |
57 | sprintf(hname, "hnhitsx%d", i); |
57 | sprintf(hname, "hnhitsx%d", i); |
58 | h_hitsx = (TH1F*)rootfile->Get(hname); |
58 | h_hitsx = (TH1F*)rootfile->Get(hname); |
59 | h_hitsx->Draw(); |
59 | h_hitsx->Draw(); |
60 | canvas3->cd(i+1); |
60 | canvas3->cd(i+1); |
61 | sprintf(hname, "hnhitsy%d", i); |
61 | sprintf(hname, "hnhitsy%d", i); |
62 | h_hitsy = (TH1F*)rootfile->Get(hname); |
62 | h_hitsy = (TH1F*)rootfile->Get(hname); |
63 | h_hitsy->Draw(); |
63 | h_hitsy->Draw(); |
64 | } |
64 | } |
65 | } |
65 | } |
66 | 66 | ||
67 | if( strstr(plopt, "x") != NULL ) { |
67 | if( strstr(plopt, "x") != NULL ) { |
68 | TCanvas *canvas10 = new TCanvas("canvas10","Ch x;;",500,500); |
68 | TCanvas *canvas10 = new TCanvas("canvas10","Ch x;;",500,500); |
69 | TH1F* h_hitsx; |
69 | TH1F* h_hitsx; |
70 | canvas10->cd(); |
70 | canvas10->cd(); |
71 | char hname[128]; |
71 | char hname[128]; |
72 | sprintf(hname, "hnhitsx%d", chXstart); |
72 | sprintf(hname, "hnhitsx%d", chXstart); |
73 | h_hitsx = (TH1F*)rootfile->Get(hname); |
73 | h_hitsx = (TH1F*)rootfile->Get(hname); |
74 | h_hitsx->Draw(); |
74 | h_hitsx->Draw(); |
75 | } |
75 | } |
76 | 76 | ||
77 | if( strstr(plopt, "y") != NULL ) { |
77 | if( strstr(plopt, "y") != NULL ) { |
78 | TCanvas *canvas11 = new TCanvas("canvas11","Ch x;;",500,500); |
78 | TCanvas *canvas11 = new TCanvas("canvas11","Ch x;;",500,500); |
79 | TH1F* h_hitsy; |
79 | TH1F* h_hitsy; |
80 | canvas11->cd(); |
80 | canvas11->cd(); |
81 | char hname[128]; |
81 | char hname[128]; |
82 | sprintf(hname, "hnhitsy%d", chXstart); |
82 | sprintf(hname, "hnhitsy%d", chXstart); |
83 | h_hitsy = (TH1F*)rootfile->Get(hname); |
83 | h_hitsy = (TH1F*)rootfile->Get(hname); |
84 | h_hitsy->Draw(); |
84 | h_hitsy->Draw(); |
85 | } |
85 | } |
86 | 86 | ||
87 | if( strstr(plopt, "share") != NULL ) { |
87 | if( strstr(plopt, "share") != NULL ) { |
88 | /*TCanvas *canvas4 = new TCanvas("canvas1","canvas1",1000,1000); |
88 | /*TCanvas *canvas4 = new TCanvas("canvas1","canvas1",1000,1000); |
89 | int nChannels = chYend-chYstart+1; |
89 | int nChannels = chYend-chYstart+1; |
90 | int ncols = nChannels/2; |
90 | int ncols = nChannels/2; |
91 | printf("nch %d nch\\2 %d\n", nChannels, ncols); |
91 | printf("nch %d nch\\2 %d\n", nChannels, ncols); |
92 | canvas4->Divide(2,ncols); |
92 | canvas4->Divide(2,ncols); |
93 | TH1F* h_hitsy; |
93 | TH1F* h_hitsy; |
94 | for(int i=chYstart; i<=chYend; i++){ |
94 | for(int i=chYstart; i<=chYend; i++){ |
95 | canvas4->cd(i-chYstart+1); |
95 | canvas4->cd(i-chYstart+1); |
96 | char hname[128]; |
96 | char hname[128]; |
97 | int chPosition = map[0][i]; |
97 | int chPosition = map[0][i]; |
98 | sprintf(hname, "hnhitsy%d", chPosition); |
98 | sprintf(hname, "hnhitsy%d", chPosition); |
99 | h_hitsy = (TH1F*)rootfile->Get(hname); |
99 | h_hitsy = (TH1F*)rootfile->Get(hname); |
100 | h_hitsy->Draw(); |
100 | h_hitsy->Draw(); |
101 | }*/ |
101 | }*/ |
102 | 102 | ||
103 | TCanvas *canvas4 = new TCanvas("canvas4","canvas4",500,500); |
103 | TCanvas *canvas4 = new TCanvas("canvas4","canvas4",500,500); |
104 | canvas4->cd(); |
104 | canvas4->cd(); |
105 | for(int i=chXstart; i<=chXend; i++) { |
105 | for(int i=chXstart; i<=chXend; i++) { |
106 | TH1F* h_hitsx; |
106 | TH1F* h_hitsx; |
107 | char hname[128]; |
107 | char hname[128]; |
108 | int chPosition = map[i][chYstart]; |
108 | int chPosition = map[i][chYstart]; |
109 | sprintf(hname, "hnhitsx%d", chPosition); |
109 | sprintf(hname, "hnhitsx%d", chPosition); |
110 | h_hitsx = (TH1F*)rootfile->Get(hname); |
110 | h_hitsx = (TH1F*)rootfile->Get(hname); |
111 | h_hitsx->SetTitle("Scan X;x [mm]; Entries"); |
111 | h_hitsx->SetTitle("Scan X;x [mm]; Entries"); |
112 | h_hitsx->GetYaxis()->SetTitleOffset(1.3); |
112 | h_hitsx->GetYaxis()->SetTitleOffset(1.3); |
113 | h_hitsx->SetStats(0); |
113 | h_hitsx->SetStats(0); |
114 | if (i == chXstart) |
114 | if (i == chXstart) |
115 | h_hitsx->Draw(); |
115 | h_hitsx->Draw(); |
116 | else { |
116 | else { |
117 | h_hitsx->SetLineColor(i+1); |
117 | h_hitsx->SetLineColor(i+1); |
118 | h_hitsx->Draw("same"); |
118 | h_hitsx->Draw("same"); |
119 | } |
119 | } |
120 | } |
120 | } |
121 | //sprintf(fullname, "ps/%s_Yshare.eps", filename); |
121 | //sprintf(fullname, "ps/%s_Yshare.eps", filename); |
122 | //canvas4->SaveAs(fullname); |
122 | //canvas4->SaveAs(fullname); |
123 | 123 | ||
124 | TCanvas *canvas5 = new TCanvas("canvas5","canvas5",500,500); |
124 | TCanvas *canvas5 = new TCanvas("canvas5","canvas5",500,500); |
125 | canvas5->cd(); |
125 | canvas5->cd(); |
126 | for(int i=chYstart; i<=chYend; i++) { |
126 | for(int i=chYstart; i<=chYend; i++) { |
127 | TH1F* h_hitsy; |
127 | TH1F* h_hitsy; |
128 | char hname[128]; |
128 | char hname[128]; |
129 | int chPosition = map[chXstart][i]; |
129 | int chPosition = map[chXstart][i]; |
130 | sprintf(hname, "hnhitsy%d", chPosition); |
130 | sprintf(hname, "hnhitsy%d", chPosition); |
131 | h_hitsy = (TH1F*)rootfile->Get(hname); |
131 | h_hitsy = (TH1F*)rootfile->Get(hname); |
132 | h_hitsy->SetTitle("Scan Y;y [mm]; Entries"); |
132 | h_hitsy->SetTitle("Scan Y;y [mm]; Entries"); |
133 | h_hitsy->GetYaxis()->SetTitleOffset(1.3); |
133 | h_hitsy->GetYaxis()->SetTitleOffset(1.3); |
134 | h_hitsy->SetStats(0); |
134 | h_hitsy->SetStats(0); |
135 | if (i == chYstart) |
135 | if (i == chYstart) |
136 | h_hitsy->Draw(); |
136 | h_hitsy->Draw(); |
137 | else { |
137 | else { |
138 | h_hitsy->SetLineColor(i+1); |
138 | h_hitsy->SetLineColor(i+1); |
139 | h_hitsy->Draw("same"); |
139 | h_hitsy->Draw("same"); |
140 | } |
140 | } |
141 | } |
141 | } |
142 | //sprintf(fullname, "ps/%s_Yshare.eps", filename); |
142 | //sprintf(fullname, "ps/%s_Yshare.eps", filename); |
143 | //canvas5->SaveAs(fullname); |
143 | //canvas5->SaveAs(fullname); |
144 | } |
144 | } |
145 | 145 | ||
146 | /** Draw one channel in x and y |
146 | /** Draw one channel in x and y |
147 | * Without and with the noise subtraction |
147 | * Without and with the noise subtraction |
148 | */ |
148 | */ |
149 | if (strstr(plopt, "noise") != NULL) { |
149 | if (strstr(plopt, "noise") != NULL) { |
150 | TCanvas *canvas51 = new TCanvas("canvas51","canvas51",500,500); |
150 | TCanvas *canvas51 = new TCanvas("canvas51","canvas51",500,500); |
151 | canvas51->cd(); |
151 | canvas51->cd(); |
152 | 152 | ||
153 | // Get the filled histogram from scan |
153 | // Get the filled histogram from scan |
154 | TH1F* h_x; |
154 | TH1F* h_x; |
155 | char name[32]; |
155 | char name[32]; |
156 | sprintf(name, "hnhitsx%d", chXstart); |
156 | sprintf(name, "hnhitsx%d", chXstart); |
157 | h_x = (TH1F*)rootfile->Get(name); |
157 | h_x = (TH1F*)rootfile->Get(name); |
158 | h_x->DrawCopy(); |
158 | h_x->DrawCopy(); |
159 | 159 | ||
160 | TH1F* h_y; |
160 | TH1F* h_y; |
161 | sprintf(name, "hnhitsy%d", chXstart); |
161 | sprintf(name, "hnhitsy%d", chXstart); |
162 | h_y = (TH1F*)rootfile->Get(name); |
162 | h_y = (TH1F*)rootfile->Get(name); |
163 | 163 | ||
164 | // Create and fill corrected histogram |
164 | // Create and fill corrected histogram |
165 | Int_t binsX = h_x->GetXaxis()->GetNbins(); |
165 | Int_t binsX = h_x->GetXaxis()->GetNbins(); |
166 | if(debug) printf("nBins: %d\n", binsX); |
166 | if(debug) printf("nBins: %d\n", binsX); |
167 | //Double_t xLowUser = h_x->GetXaxis()->GetBinLowEdge(1); |
167 | //Double_t xLowUser = h_x->GetXaxis()->GetBinLowEdge(1); |
168 | //Double_t xUpUser = h_x->GetXaxis()->GetBinUpEdge(binsX); |
168 | //Double_t xUpUser = h_x->GetXaxis()->GetBinUpEdge(binsX); |
169 | Double_t xLowUser = h_x->GetXaxis()->GetXmin(); |
169 | Double_t xLowUser = h_x->GetXaxis()->GetXmin(); |
170 | Double_t xUpUser = h_x->GetXaxis()->GetXmax(); |
170 | Double_t xUpUser = h_x->GetXaxis()->GetXmax(); |
171 | TH1F* h_xCorrected = new TH1F("h_xCorrected",";;",binsX, xLowUser, xUpUser); |
171 | TH1F* h_xCorrected = new TH1F("h_xCorrected",";;",binsX, xLowUser, xUpUser); |
172 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
172 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
173 | double noise = h_x->GetBinContent(2); |
173 | double noise = h_x->GetBinContent(2); |
174 | for (int jk = 0; jk < binsX; jk++) { |
174 | for (int jk = 0; jk < binsX; jk++) { |
175 | double signal = h_x->GetBinContent(jk+1); |
175 | double signal = h_x->GetBinContent(jk+1); |
176 | signal -= noise; |
176 | signal -= noise; |
177 | double x = h_x->GetXaxis()->GetBinCenter(jk+1); |
177 | double x = h_x->GetXaxis()->GetBinCenter(jk+1); |
178 | h_xCorrected->Fill(x, signal); |
178 | h_xCorrected->Fill(x, signal); |
179 | } |
179 | } |
180 | h_xCorrected->SetLineColor(2); |
180 | h_xCorrected->SetLineColor(2); |
181 | h_xCorrected->DrawCopy("same"); |
181 | h_xCorrected->DrawCopy("same"); |
182 | 182 | ||
183 | //Create and fill corrected histogram |
183 | //Create and fill corrected histogram |
184 | Int_t binsY = h_y->GetXaxis()->GetNbins(); |
184 | Int_t binsY = h_y->GetXaxis()->GetNbins(); |
185 | Double_t yLowUser = h_y->GetXaxis()->GetXmin(); |
185 | Double_t yLowUser = h_y->GetXaxis()->GetXmin(); |
186 | Double_t yUpUser = h_y->GetXaxis()->GetXmax(); |
186 | Double_t yUpUser = h_y->GetXaxis()->GetXmax(); |
187 | TH1F* h_yCorrected = new TH1F("h_yCorrected",";;",binsY, yLowUser, yUpUser); |
187 | TH1F* h_yCorrected = new TH1F("h_yCorrected",";;",binsY, yLowUser, yUpUser); |
188 | noise = h_y->GetBinContent(2); |
188 | noise = h_y->GetBinContent(2); |
189 | for (int jk = 0; jk < binsY; jk++) { |
189 | for (int jk = 0; jk < binsY; jk++) { |
190 | double signal = h_y->GetBinContent(jk+1); |
190 | double signal = h_y->GetBinContent(jk+1); |
191 | signal -= noise; |
191 | signal -= noise; |
192 | double x = h_y->GetXaxis()->GetBinCenter(jk+1); |
192 | double x = h_y->GetXaxis()->GetBinCenter(jk+1); |
193 | h_yCorrected->Fill(x, signal); |
193 | h_yCorrected->Fill(x, signal); |
194 | } |
194 | } |
195 | TCanvas *canvas52 = new TCanvas("canvas52","canvas52",500,500); |
195 | TCanvas *canvas52 = new TCanvas("canvas52","canvas52",500,500); |
196 | canvas52->cd(); |
196 | canvas52->cd(); |
197 | h_y->DrawCopy(); |
197 | h_y->DrawCopy(); |
198 | h_yCorrected->SetLineColor(2); |
198 | h_yCorrected->SetLineColor(2); |
199 | h_yCorrected->DrawCopy("same"); |
199 | h_yCorrected->DrawCopy("same"); |
200 | } |
200 | } |
201 | 201 | ||
202 | /** Draws the signal from 8 channels in x-row |
202 | /** Draws the signal from 8 channels in x-row |
203 | * for one specific y bin, so the background and cross-talk |
203 | * for one specific y bin, so the background and cross-talk |
204 | * can be estimated. |
204 | * can be estimated. |
205 | * Draws also a 2d scan of these channels. |
205 | * Draws also a 2d scan of these channels. |
206 | */ |
206 | */ |
207 | if (strstr(plopt, "line") != NULL) { |
207 | if (strstr(plopt, "line") != NULL) { |
208 | TCanvas* canvas6 = new TCanvas("canvas6","canvas6",500,500); |
208 | TCanvas* canvas6 = new TCanvas("canvas6","canvas6",500,500); |
209 | canvas6->cd(0); |
209 | canvas6->cd(0); |
210 | gStyle->SetOptStat(0); |
210 | gStyle->SetOptStat(0); |
211 | 211 | ||
212 | TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
212 | TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
213 | Int_t binsX = h0->GetXaxis()->GetNbins(); |
213 | Int_t binsX = h0->GetXaxis()->GetNbins(); |
214 | Int_t minX = h0->GetXaxis()->GetFirst(); |
214 | Int_t minX = h0->GetXaxis()->GetFirst(); |
215 | Int_t maxX = h0->GetXaxis()->GetLast()+1; |
215 | Int_t maxX = h0->GetXaxis()->GetLast()+1; |
216 | Int_t binsY = h0->GetYaxis()->GetNbins(); |
216 | Int_t binsY = h0->GetYaxis()->GetNbins(); |
217 | Int_t minY = h0->GetYaxis()->GetFirst(); |
217 | Int_t minY = h0->GetYaxis()->GetFirst(); |
218 | Int_t maxY = h0->GetYaxis()->GetLast()+1; |
218 | Int_t maxY = h0->GetYaxis()->GetLast()+1; |
219 | Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
219 | Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
220 | Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
220 | Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
221 | Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
221 | Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
222 | Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
222 | Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
223 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
223 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
224 | 224 | ||
225 | //! 1-dimension position in x vs. hits |
225 | //! 1-dimension position in x vs. hits |
226 | TH2F* h[8]; |
226 | TH2F* h[8]; |
227 | TH1F* h_line[8]; |
227 | TH1F* h_line[8]; |
228 | for(int j=0; j<8; j++) { |
228 | for(int j=0; j<8; j++) { |
229 | h_line[j] = new TH1F("h_line", "h_line", binsX, xLowUser, xUpUser); |
229 | h_line[j] = new TH1F("h_line", "h_line", binsX, xLowUser, xUpUser); |
230 | } |
230 | } |
231 | 231 | ||
232 | for(int j=chXstart; j<=chXend; j++) { |
232 | for(int j=chXstart; j<=chXend; j++) { |
233 | int chPosition = map[j][chYstart]; |
233 | int chPosition = map[j][chYstart]; |
234 | char hname[128]; |
234 | char hname[128]; |
235 | sprintf(hname, "h2d%d", chPosition); |
235 | sprintf(hname, "h2d%d", chPosition); |
236 | int histogram = j; |
236 | int histogram = j; |
237 | h[histogram] = (TH2F *) rootfile->Get(hname); |
237 | h[histogram] = (TH2F *) rootfile->Get(hname); |
238 | int noise = getNoise(h[histogram], 1, 160); |
238 | int noise = getNoise(h[histogram], 1, 160); |
239 | for(int k=minX; k<=maxX; k++) { |
239 | for(int k=minX; k<=maxX; k++) { |
240 | int l=chYstart*20+12; |
240 | int l=chYstart*20+12; |
241 | //for(int l=12; l<=16; l++) { |
241 | //for(int l=12; l<=16; l++) { |
242 | double signal = h[histogram]->GetBinContent(k,l); |
242 | double signal = h[histogram]->GetBinContent(k,l); |
243 | //signal -= noise; |
243 | //signal -= noise; |
244 | //signal /= 5*10000.0; |
244 | //signal /= 5*10000.0; |
245 | double eta = -log(1 - signal); |
245 | double eta = -log(1 - signal); |
246 | double x = xLowUser + k*(xUpUser-xLowUser)/double(binsX); |
246 | double x = xLowUser + k*(xUpUser-xLowUser)/double(binsX); |
247 | //double y = l*(yUpUser-yLowUser)/double(binsY); |
247 | //double y = l*(yUpUser-yLowUser)/double(binsY); |
248 | h_line[j]->Fill(x-c_xOffset, signal); |
248 | h_line[j]->Fill(x-c_xOffset, signal); |
249 | //} |
249 | //} |
250 | } |
250 | } |
251 | if (j == chXstart) { |
251 | if (j == chXstart) { |
252 | h_line[j]->SetTitle("SiPM#2 w/o noise subtraction;x[mm];Hits"); |
252 | h_line[j]->SetTitle("SiPM#2 w/o noise subtraction;x[mm];Hits"); |
253 | //h_line[j]->GetYaxis()->SetRangeUser(-0.05, 0.3); |
253 | //h_line[j]->GetYaxis()->SetRangeUser(-0.05, 0.3); |
254 | //h_line[j]->GetYaxis()->SetRangeUser(-50, 2500); |
254 | //h_line[j]->GetYaxis()->SetRangeUser(-50, 2500); |
255 | h_line[j]->Draw(""); |
255 | h_line[j]->Draw(""); |
256 | } |
256 | } |
257 | else { |
257 | else { |
258 | h_line[j]->SetLineColor(j+1); |
258 | h_line[j]->SetLineColor(j+1); |
259 | h_line[j]->Draw("same"); |
259 | h_line[j]->Draw("same"); |
260 | } |
260 | } |
261 | } |
261 | } |
262 | 262 | ||
263 | 263 | ||
264 | //! 2d scan |
264 | //! 2d scan |
265 | TCanvas* canvas61 = new TCanvas("canvas61","canvas61",8*200,300); |
265 | TCanvas* canvas61 = new TCanvas("canvas61","canvas61",8*200,300); |
266 | canvas61->cd(); |
266 | canvas61->cd(); |
267 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
267 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
268 | for(int i=chXstart; i<=chXend; i++) { |
268 | for(int i=chXstart; i<=chXend; i++) { |
269 | //int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
269 | //int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
270 | //int canvasPosition = nX*(chYend-i)+chXstart+1; |
270 | //int canvasPosition = nX*(chYend-i)+chXstart+1; |
271 | //if (debug) printf("canvas %d\n",canvasPosition); |
271 | //if (debug) printf("canvas %d\n",canvasPosition); |
272 | int chPosition = map[i][chYstart]; |
272 | int chPosition = map[i][chYstart]; |
273 | char hname[128]; |
273 | char hname[128]; |
274 | sprintf(hname, "h2d%d", chPosition); |
274 | sprintf(hname, "h2d%d", chPosition); |
275 | int histogram = i; |
275 | int histogram = i; |
276 | h[histogram] = (TH2F *) rootfile->Get(hname); |
276 | h[histogram] = (TH2F *) rootfile->Get(hname); |
277 | int noise = getNoise(h[histogram], 1, 100); |
277 | int noise = getNoise(h[histogram], 1, 100); |
278 | 278 | ||
279 | for(int k=minX; k<=maxX; k++) { |
279 | for(int k=minX; k<=maxX; k++) { |
280 | for(int l=minY; l<=maxY; l++) { |
280 | for(int l=minY; l<=maxY; l++) { |
281 | int signal = h[histogram]->GetBinContent(k,l); // detected |
281 | int signal = h[histogram]->GetBinContent(k,l); // detected |
282 | //p /= 10000.; |
282 | //p /= 10000.; |
283 | //double p0 = 1.0 - p; // events with zero photons |
283 | //double p0 = 1.0 - p; // events with zero photons |
284 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
284 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
285 | //double eta = -log(p0); |
285 | //double eta = -log(p0); |
286 | //printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
286 | //printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
287 | //double signal = ((p - noise) > 0.1) ? (p-noise) : 0.1; |
287 | //double signal = ((p - noise) > 0.1) ? (p-noise) : 0.1; |
288 | double p = signal - noise; |
288 | double p = signal - noise; |
289 | p /= 10000.0; |
289 | p /= 10000.0; |
290 | double eta = -log(1 - p); |
290 | double eta = -log(1 - p); |
291 | //double x = k*(xUpUser-xLowUser)/double(binsX); |
291 | //double x = k*(xUpUser-xLowUser)/double(binsX); |
292 | //double y = l*(yUpUser-yLowUser)/double(binsY); |
292 | //double y = l*(yUpUser-yLowUser)/double(binsY); |
293 | double x = h[histogram]->GetXaxis()->GetBinCenter(k); |
293 | double x = h[histogram]->GetXaxis()->GetBinCenter(k); |
294 | double y = h[histogram]->GetYaxis()->GetBinCenter(l); |
294 | double y = h[histogram]->GetYaxis()->GetBinCenter(l); |
295 | h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
295 | h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
296 | } |
296 | } |
297 | } |
297 | } |
298 | } |
298 | } |
299 | h_corrected->SetTitle("SiPM#2 n_pe = - ln(P0);x[mm];y[mm]"); |
299 | h_corrected->SetTitle("SiPM#2 n_pe = - ln(P0);x[mm];y[mm]"); |
300 | h_corrected->GetZaxis()->SetRangeUser(-0.05,0.30); |
300 | h_corrected->GetZaxis()->SetRangeUser(-0.05,0.30); |
301 | h_corrected->Draw("colz"); |
301 | h_corrected->Draw("colz"); |
302 | 302 | ||
303 | } |
303 | } |
304 | 304 | ||
305 | /** Draws the sum of the channels |
305 | /** Draws the sum of the channels |
306 | * Each channel is a 2d plot |
306 | * Each channel is a 2d plot |
307 | * Intended for the study of 1 channel |
307 | * Intended for the study of 1 channel |
308 | * |
308 | * |
309 | * Input: scanned channel coordinates (x,y) |
309 | * Input: scanned channel coordinates (x,y) |
310 | */ |
310 | */ |
311 | if (strstr(plopt, "2d") != NULL) { |
311 | if (strstr(plopt, "2d") != NULL) { |
312 | 312 | ||
313 | int nX = 3; |
313 | int nX = 3; |
314 | int nY = 3; |
314 | int nY = 3; |
315 | TCanvas* canvas7 = new TCanvas("canvas7","canvas7", nX*400,nY*400); |
315 | TCanvas* canvas7 = new TCanvas("canvas7","canvas7", nX*400,nY*400); |
316 | printf("nx %d ny %d\n",nX,nY); |
316 | printf("nx %d ny %d\n",nX,nY); |
317 | canvas7->Divide(nX,nY); |
317 | canvas7->Divide(nX,nY); |
318 | for(int i=chXend-1; i<=chXend+1; i++) { |
318 | for(int i=chXend-1; i<=chXend+1; i++) { |
319 | for(int j=chXstart-1; j<=chXstart+1; j++) { |
319 | for(int j=chXstart-1; j<=chXstart+1; j++) { |
320 | //int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
320 | //int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
321 | int canvasPosition = nX*(chXend+1-i) + (j-chXstart+1) +1; |
321 | int canvasPosition = nX*(chXend+1-i) + (j-chXstart+1) +1; |
322 | if (debug) printf("canvas %d\n",canvasPosition); |
322 | if (debug) printf("canvas %d\n",canvasPosition); |
323 | canvas7->cd(canvasPosition); |
323 | canvas7->cd(canvasPosition); |
324 | char hname[128]; |
324 | char hname[128]; |
325 | int chPosition = map[j][i]; |
325 | int chPosition = map[j][i]; |
326 | sprintf(hname, "h2d%d", chPosition); |
326 | sprintf(hname, "h2d%d", chPosition); |
327 | TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
327 | TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
328 | h_2d->Draw("colz"); |
328 | h_2d->Draw("colz"); |
329 | } //x |
329 | } //x |
330 | } |
330 | } |
331 | 331 | ||
332 | // Number of photoelectrons - Poissonian correction |
332 | // Number of photoelectrons - Poissonian correction |
333 | TCanvas* canvas8 = new TCanvas("canvas8","canvas8", 1000,1000); |
333 | TCanvas* canvas8 = new TCanvas("canvas8","canvas8", 1000,1000); |
334 | canvas8->cd(); |
334 | canvas8->cd(); |
335 | gStyle->SetOptStat(0); |
335 | gStyle->SetOptStat(0); |
336 | char hname[128]; |
336 | char hname[128]; |
337 | int chPosition = map[chXstart][chXend]; |
337 | int chPosition = map[chXstart][chXend]; |
338 | sprintf(hname, "h2d%d", chPosition); |
338 | sprintf(hname, "h2d%d", chPosition); |
339 | TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
339 | TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
340 | 340 | ||
341 | Int_t binsX = h_2d->GetXaxis()->GetNbins(); |
341 | Int_t binsX = h_2d->GetXaxis()->GetNbins(); |
342 | Int_t minX = h_2d->GetXaxis()->GetFirst(); |
342 | Int_t minX = h_2d->GetXaxis()->GetFirst(); |
343 | Int_t maxX = h_2d->GetXaxis()->GetLast()+1; |
343 | Int_t maxX = h_2d->GetXaxis()->GetLast()+1; |
344 | Int_t binsY = h_2d->GetYaxis()->GetNbins(); |
344 | Int_t binsY = h_2d->GetYaxis()->GetNbins(); |
345 | Int_t minY = h_2d->GetYaxis()->GetFirst(); |
345 | Int_t minY = h_2d->GetYaxis()->GetFirst(); |
346 | Int_t maxY = h_2d->GetYaxis()->GetLast()+1; |
346 | Int_t maxY = h_2d->GetYaxis()->GetLast()+1; |
347 | Double_t xLowUser = h_2d->GetXaxis()->GetXmin(); |
347 | Double_t xLowUser = h_2d->GetXaxis()->GetXmin(); |
348 | Double_t xUpUser = h_2d->GetXaxis()->GetXmax(); |
348 | Double_t xUpUser = h_2d->GetXaxis()->GetXmax(); |
349 | Double_t yLowUser = h_2d->GetYaxis()->GetXmin(); |
349 | Double_t yLowUser = h_2d->GetYaxis()->GetXmin(); |
350 | Double_t yUpUser = h_2d->GetYaxis()->GetXmax(); |
350 | Double_t yUpUser = h_2d->GetYaxis()->GetXmax(); |
351 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
351 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
352 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, xLowUser, xUpUser, binsY, yLowUser, yUpUser); |
352 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, xLowUser, xUpUser, binsY, yLowUser, yUpUser); |
353 | //TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, minX, maxX, binsY, minY, maxY); |
353 | //TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, minX, maxX, binsY, minY, maxY); |
354 | 354 | ||
355 | double noise = getNoise(h_2d, 1, 70); |
355 | double noise = getNoise(h_2d, 1, 70); |
356 | if(debug) printf("Noise = %f\n", noise); |
356 | if(debug) printf("Noise = %f\n", noise); |
357 | for(int k=minX; k<=maxX; k++) { |
357 | for(int k=minX; k<=maxX; k++) { |
358 | for(int j=minY; j<=maxY; j++) { |
358 | for(int j=minY; j<=maxY; j++) { |
359 | double signal = h_2d->GetBinContent(k,j); // detected |
359 | double signal = h_2d->GetBinContent(k,j); // detected |
360 | //double p = ((signal - noise) > 1) ? (signal-noise) : 1; |
360 | //double p = ((signal - noise) > 1) ? (signal-noise) : 1; |
361 | double p = signal - noise; |
361 | double p = signal - noise; |
362 | p /= 10000.; |
362 | p /= 10000.; |
363 | double p0 = 1.0 - p; // events with zero photons |
363 | double p0 = 1.0 - p; // events with zero photons |
364 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
364 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
365 | double eta = -log(p0); // constant of the poissonian statistics |
365 | double eta = -log(p0); // constant of the poissonian statistics |
366 | if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
366 | if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
367 | //double x = xLowUser + k*(xUpUser - xLowUser) / double(binsX); |
367 | //double x = xLowUser + k*(xUpUser - xLowUser) / double(binsX); |
368 | double x = h_2d->GetXaxis()->GetBinCenter(k); |
368 | double x = h_2d->GetXaxis()->GetBinCenter(k); |
369 | //double y = yLowUser + j*(yUpUser-yLowUser)/double(binsY); |
369 | //double y = yLowUser + j*(yUpUser-yLowUser)/double(binsY); |
370 | double y = h_2d->GetYaxis()->GetBinCenter(j); |
370 | double y = h_2d->GetYaxis()->GetBinCenter(j); |
371 | if (debug) printf("x=%f y=%f\n", x, y); |
371 | if (debug) printf("x=%f y=%f\n", x, y); |
372 | h_corrected->Fill(x + c_xOffset, y + c_yOffset, eta); |
372 | h_corrected->Fill(x + c_xOffset, y + c_yOffset, eta); |
373 | } |
373 | } |
374 | } |
374 | } |
375 |
|
375 | h_corrected->SetTitle("n_pe = - ln P(0);x[mm];y[mm]"); |
376 | h_corrected |
376 | //h_corrected->SetTitle(";x[mm];y[mm]"); |
- | 377 | gStyle->SetPalette(52,0); |
|
377 | h_corrected->GetZaxis()->SetRangeUser(-0. |
378 | h_corrected->GetZaxis()->SetRangeUser(-0.05,0.7); |
378 | h_corrected->SetContour(50); |
379 | h_corrected->SetContour(50); |
379 | h_corrected->Draw("colz"); |
380 | h_corrected->Draw("colz"); |
380 | 381 | ||
381 | // collection efficiency |
382 | // collection efficiency |
382 | int nPoints =0; |
383 | int nPoints =0; |
383 | double efficiency=0; |
384 | double efficiency=0; |
384 | for (int i=18; i<=58; i++) { |
385 | for (int i=18; i<=58; i++) { |
385 | for (int j=19; j<=59; j++) { |
386 | for (int j=19; j<=59; j++) { |
386 | double signal = h_corrected->GetBinContent(i,j); |
387 | double signal = h_corrected->GetBinContent(i,j); |
387 | if(debug) printf("signal %f\n",signal); |
388 | if(debug) printf("signal %f\n",signal); |
388 | efficiency += signal; |
389 | efficiency += signal; |
389 | nPoints++; |
390 | nPoints++; |
390 | } |
391 | } |
391 | } |
392 | } |
392 | printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
393 | printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
393 | } |
394 | } |
394 | 395 | ||
395 | /** Draws the sum of channel signals |
396 | /** Draws the sum of channel signals |
396 | * Each channel is a 2d ('h2d') histogram |
397 | * Each channel is a 2d ('h2d') histogram |
397 | * Suitable for 8x8 chs scan |
398 | * Suitable for 8x8 chs scan |
398 | */ |
399 | */ |
399 | if( strstr(plopt, "sum") != NULL ) { |
400 | if( strstr(plopt, "sum") != NULL ) { |
400 | int nX = chXend - chXstart + 1; |
401 | int nX = chXend - chXstart + 1; |
401 | int nY = chYend - chYstart + 1; |
402 | int nY = chYend - chYstart + 1; |
402 | TCanvas* canvas12 = new TCanvas("canvas12","c2",8*200, 8*200); |
403 | TCanvas* canvas12 = new TCanvas("canvas12","c2",8*200, 8*200); |
403 | canvas12->cd(); |
404 | canvas12->cd(); |
404 | gStyle->SetOptStat(0); |
405 | gStyle->SetOptStat(0); |
405 | 406 | ||
406 | // final histogram parameters |
407 | // final histogram parameters |
407 | TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
408 | TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
408 | Int_t binsX = h0->GetXaxis()->GetNbins(); |
409 | Int_t binsX = h0->GetXaxis()->GetNbins(); |
409 | Int_t minX = h0->GetXaxis()->GetFirst(); |
410 | Int_t minX = h0->GetXaxis()->GetFirst(); |
410 | Int_t maxX = h0->GetXaxis()->GetLast()+1; |
411 | Int_t maxX = h0->GetXaxis()->GetLast()+1; |
411 | Int_t binsY = h0->GetYaxis()->GetNbins(); |
412 | Int_t binsY = h0->GetYaxis()->GetNbins(); |
412 | Int_t minY = h0->GetYaxis()->GetFirst(); |
413 | Int_t minY = h0->GetYaxis()->GetFirst(); |
413 | Int_t maxY = h0->GetYaxis()->GetLast()+1; |
414 | Int_t maxY = h0->GetYaxis()->GetLast()+1; |
414 | Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
415 | Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
415 | Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
416 | Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
416 | Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
417 | Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
417 | Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
418 | Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
418 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
419 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
419 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
420 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
420 | TH2F* h[nX*nY]; |
421 | TH2F* h[nX*nY]; |
421 | 422 | ||
422 | // 2d histogram noise subtraction and poisson scaling |
423 | // 2d histogram noise subtraction and poisson scaling |
423 | for(int i=chYstart; i<=chYend; i++) { |
424 | for(int i=chYstart; i<=chYend; i++) { |
424 | for(int j=chXstart; j<=chXend; j++) { |
425 | for(int j=chXstart; j<=chXend; j++) { |
425 | int chPosition = map[j][i]; |
426 | int chPosition = map[j][i]; |
426 | char hname[128]; |
427 | char hname[128]; |
427 | sprintf(hname, "h2d%d", chPosition); |
428 | sprintf(hname, "h2d%d", chPosition); |
428 | int histogram = nX*(i-chYstart)+(j-chXstart); |
429 | int histogram = nX*(i-chYstart)+(j-chXstart); |
429 | h[histogram] = (TH2F *) rootfile->Get(hname); |
430 | h[histogram] = (TH2F *) rootfile->Get(hname); |
430 | int noise = getNoise(h[histogram], 1, 170); |
431 | int noise = getNoise(h[histogram], 1, 170); |
431 | if (debug) printf("noise: %d\n",noise); |
432 | if (debug) printf("noise: %d\n",noise); |
432 | for(int k=minX; k<=maxX; k++) { |
433 | for(int k=minX; k<=maxX; k++) { |
433 | for(int l=minY; l<=maxY; l++) { |
434 | for(int l=minY; l<=maxY; l++) { |
434 | int signal = h[histogram]->GetBinContent(k,l); // detected |
435 | int signal = h[histogram]->GetBinContent(k,l); // detected |
435 | //double p = ((signal - noise) > 0.1) ? (signal-noise) : 0.1; |
436 | //double p = ((signal - noise) > 0.1) ? (signal-noise) : 0.1; |
436 | double p = signal - noise; |
437 | double p = signal - noise; |
437 | p /= 10000.; |
438 | p /= 10000.; |
438 | double p0 = 1.0 - p; // events with zero photons |
439 | double p0 = 1.0 - p; // events with zero photons |
439 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
440 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
440 | double eta = -log(p0); |
441 | double eta = -log(p0); |
441 | //printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
442 | //printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
442 | double x = k*(xUpUser-xLowUser)/double(binsX); |
443 | double x = k*(xUpUser-xLowUser)/double(binsX); |
443 | double y = l*(yUpUser-yLowUser)/double(binsY); |
444 | double y = l*(yUpUser-yLowUser)/double(binsY); |
444 | h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
445 | h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
445 | } |
446 | } |
446 | } |
447 | } |
447 | 448 | ||
448 | } |
449 | } |
449 | } |
450 | } |
450 | 451 | ||
451 | h_corrected->SetTitle("SiPM#2 n_p.e.;x[mm];y[mm]"); |
452 | h_corrected->SetTitle("SiPM#2 n_p.e.;x[mm];y[mm]"); |
452 | h_corrected->GetZaxis()->SetRangeUser(-0.05,.30); |
453 | h_corrected->GetZaxis()->SetRangeUser(-0.05,.30); |
453 | h_corrected->Draw("colz"); |
454 | h_corrected->Draw("colz"); |
454 | 455 | ||
455 | TCanvas* canvas13 = new TCanvas("canvas13","canvas13",600,300); |
456 | TCanvas* canvas13 = new TCanvas("canvas13","canvas13",600,300); |
456 | canvas13->Divide(2); |
457 | canvas13->Divide(2); |
457 | canvas13->cd(1); |
458 | canvas13->cd(1); |
458 | h[16]->Draw("colz"); |
459 | h[16]->Draw("colz"); |
459 | canvas13->cd(2); |
460 | canvas13->cd(2); |
460 | h[8]->Draw("colz"); |
461 | h[8]->Draw("colz"); |
461 | } |
462 | } |
462 | 463 | ||
463 | /** Draws the beam profile and fits it with error function |
464 | /** Draws the beam profile and fits it with error function |
464 | * on some background function |
465 | * on some background function |
465 | */ |
466 | */ |
466 | if (strstr(plopt, "beam") != NULL) { |
467 | if (strstr(plopt, "beam") != NULL) { |
467 | 468 | ||
468 | TCanvas* canvas9 = new TCanvas("canvas9","canvas9", 500,500); |
469 | TCanvas* canvas9 = new TCanvas("canvas9","canvas9", 500,500); |
469 | canvas9->cd(); |
470 | canvas9->cd(); |
470 | char hname[128]; |
471 | char hname[128]; |
471 | sprintf(hname, "hnhitsx%d", 36); |
472 | sprintf(hname, "hnhitsx%d", 36); |
472 | TH1F* h_laser = (TH1F*)rootfile->Get(hname); |
473 | TH1F* h_laser = (TH1F*)rootfile->Get(hname); |
473 | h_laser->Draw(); |
474 | h_laser->Draw(); |
474 | h_laser->SetStats(1); |
475 | h_laser->SetStats(1); |
475 | 476 | ||
476 | TF1* err = new TF1("err","[0]+[1]*TMath::Erf((x-[2])/[3])",17.02,17.30); |
477 | TF1* err = new TF1("err","[0]+[1]*TMath::Erf((x-[2])/[3])",17.02,17.30); |
477 | err->SetParameter(0,2500); |
478 | err->SetParameter(0,2500); |
478 | err->SetParameter(1, h_laser->GetMaximum()); |
479 | err->SetParameter(1, h_laser->GetMaximum()); |
479 | err->SetParameter(2, h_laser->GetBinCenter(h_laser->GetMaximumBin())); |
480 | err->SetParameter(2, h_laser->GetBinCenter(h_laser->GetMaximumBin())); |
480 | err->SetParameter(3, 0.001); |
481 | err->SetParameter(3, 0.001); |
481 | h_laser->Fit(err,"qr"); |
482 | h_laser->Fit(err,"qr"); |
482 | h_laser->Fit(err,"lr"); |
483 | h_laser->Fit(err,"lr"); |
483 | double sigma = err->GetParameter(3); |
484 | double sigma = err->GetParameter(3); |
484 | printf("sigma = %2.0f um, FWHM = %2.0f um\n", sigma*1000, 2.35*sigma*1000); |
485 | printf("sigma = %2.0f um, FWHM = %2.0f um\n", sigma*1000, 2.35*sigma*1000); |
485 | } |
486 | } |
486 | 487 | ||
487 | if (strstr(plopt, "map") != NULL) { |
488 | if (strstr(plopt, "map") != NULL) { |
488 | TCanvas* canvas1 = new TCanvas("canvas1","canvas1",1000,1000); |
489 | TCanvas* canvas1 = new TCanvas("canvas1","canvas1",1000,1000); |
489 | canvas1->cd(); |
490 | canvas1->cd(); |
490 | TH2I* h_map = new TH2I("h_map","h_map",8,-0.5,7.5,8,-0.5,7.5); |
491 | TH2I* h_map = new TH2I("h_map","h_map",8,-0.5,7.5,8,-0.5,7.5); |
491 | for (int i=7; i>=0; i--) { |
492 | for (int i=7; i>=0; i--) { |
492 | //for (int j=7; j>=0; j--) printf("(%d, %d) ", j,i); |
493 | //for (int j=7; j>=0; j--) printf("(%d, %d) ", j,i); |
493 | for (int j=7; j>=0; j--) { |
494 | for (int j=7; j>=0; j--) { |
494 | printf("%2d (%d %d)\n", map[j][i], j*10080+5040, i*10080+5040); |
495 | printf("%2d (%d %d)\n", map[j][i], j*10080+5040, i*10080+5040); |
495 | h_map->Fill(j,i, map[j][i]); |
496 | h_map->Fill(j,i, map[j][i]); |
496 | } |
497 | } |
497 | printf("\n"); |
498 | printf("\n"); |
498 | } |
499 | } |
499 | gStyle->SetOptStat(0); |
500 | gStyle->SetOptStat(0); |
500 | h_map->Draw("text"); |
501 | h_map->Draw("text"); |
501 | } |
502 | } |
502 | 503 | ||
503 | return(0); |
504 | return(0); |
504 | } |
505 | } |
505 | 506 | ||
506 | /** Function calculates the noise from one channel |
507 | /** Function calculates the noise from one channel |
507 | * it runs through the bins along x and returns the average value |
508 | * it runs through the bins along x and returns the average value |
508 | */ |
509 | */ |
509 | double getNoise(TH2F* histogram, int yStart, int yEnd) |
510 | double getNoise(TH2F* histogram, int yStart, int yEnd) |
510 | { |
511 | { |
511 | double noise=0; |
512 | double noise=0; |
512 | int count=0; |
513 | int count=0; |
513 | for(int j=yStart; j<yEnd; j++) { |
514 | for(int j=yStart; j<yEnd; j++) { |
514 | double value = histogram->GetBinContent(j,2); |
515 | double value = histogram->GetBinContent(j,2); |
515 | //if (noise < value) noise = value; |
516 | //if (noise < value) noise = value; |
516 | noise += value; |
517 | noise += value; |
517 | count++; |
518 | count++; |
518 | } |
519 | } |
519 | return (noise/double(count)); |
520 | return (noise/double(count)); |
520 | } |
521 | } |
521 | 522 | ||
522 | 523 |