Rev 50 | Details | Compare with Previous | Last modification | View Log | RSS feed
| Rev | Author | Line No. | Line |
|---|---|---|---|
| 35 | f9daq | 1 | #include "TROOT.h" |
| 2 | #include "TFile.h" |
||
| 3 | #include "TBenchmark.h" |
||
| 4 | #include "TH1F.h" |
||
| 5 | #include "TH2F.h" |
||
| 6 | #include "TH3F.h" |
||
| 7 | #include "TCanvas.h" |
||
| 8 | #include "TStyle.h" |
||
| 9 | #include "TPad.h" |
||
| 10 | #include "TF1.h" |
||
| 11 | #include "TGraph.h" |
||
| 12 | #include "TSpectrum.h" |
||
| 13 | #include "stdio.h" |
||
| 38 | f9daq | 14 | #include "THStack.h" |
| 40 | f9daq | 15 | #include "TPaveText.h" |
| 35 | f9daq | 16 | |
| 17 | #include "include/RTUtil.h" |
||
| 18 | |||
| 19 | double getNoise(TH2F*, int, int); |
||
| 20 | |||
| 40 | f9daq | 21 | int biasScan(char filename[256] = "test", char plopt[256]="th", int chXstart=0, int chXend=7, double par1=0.0, double par2=0.0, bool debug = false) |
| 35 | f9daq | 22 | { |
| 23 | const int c_nChannels = 64; |
||
| 40 | f9daq | 24 | const double c_xOffset = 0; // mm |
| 25 | const double c_yOffset = 0; |
||
| 35 | f9daq | 26 | |
| 27 | int map[8][8]={{32,34,53,55,40,42,61,63}, |
||
| 28 | {48,50,37,39,56,58,45,47}, |
||
| 29 | {33,35,52,54,41,43,60,62}, |
||
| 30 | {49,51,36,38,57,59,44,46}, |
||
| 31 | {17,19,4,6,25,27,12,14}, |
||
| 32 | {1,3,20,22,9,11,28,30}, |
||
| 33 | {16,18,5,7,24,26,13,15}, |
||
| 34 | {0,2,21,23,8,10,29,31} |
||
| 35 | }; |
||
| 36 | |||
| 38 | f9daq | 37 | // Set draw style |
| 35 | f9daq | 38 | RTSetStyle(gStyle); |
| 39 | |||
| 40 | // open the file with histograms |
||
| 41 | char fnameroot[256]; |
||
| 42 | TFile* rootfile; |
||
| 43 | sprintf(fnameroot, "root/%s.root", filename); |
||
| 44 | rootfile = (TFile *) gROOT->FindObject(filename); |
||
| 45 | if(rootfile==NULL) rootfile = new TFile(fnameroot); |
||
| 46 | if(rootfile==NULL) { |
||
| 47 | printf("Cannot open root file %s!!!\n",fnameroot); |
||
| 48 | return(0); |
||
| 49 | } |
||
| 50 | |||
| 51 | if(strstr(plopt, "tdc") != NULL) { |
||
| 52 | TCanvas *canvas21 = new TCanvas("canvas21","canvas21",1600,800); |
||
| 53 | TH3F* h0 = (TH3F*) rootfile->Get("h_correctedTDC"); |
||
| 54 | Int_t binsX = h0->GetXaxis()->GetNbins(); |
||
| 55 | Int_t minX = h0->GetXaxis()->GetFirst(); |
||
| 56 | Int_t maxX = h0->GetXaxis()->GetLast()+1; |
||
| 57 | Int_t binsY = h0->GetYaxis()->GetNbins(); |
||
| 58 | Int_t minY = h0->GetYaxis()->GetFirst(); |
||
| 59 | Int_t maxY = h0->GetYaxis()->GetLast()+1; |
||
| 60 | Int_t binsZ = h0->GetZaxis()->GetNbins(); |
||
| 61 | Int_t minZ = h0->GetZaxis()->GetFirst(); |
||
| 62 | Int_t maxZ = h0->GetZaxis()->GetLast()+1; |
||
| 63 | Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
||
| 64 | Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
||
| 65 | Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
||
| 66 | Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
||
| 67 | |||
| 68 | canvas21->Divide(2); |
||
| 69 | //h0->Draw(); |
||
| 70 | canvas21->cd(1); |
||
| 71 | //h0->Project3D("xy")->Draw("colz"); |
||
| 72 | //canvas21->cd(2); |
||
| 73 | //h0->Project3D("yz")->Draw("colz"); |
||
| 74 | //canvas21->cd(3); |
||
| 75 | TH2D* h_correctedTDC = (TH2D*) h0->Project3D("xz"); |
||
| 76 | //h0->GetZaxis()->SetRangeUser(-5,5); |
||
| 77 | //h0->Project3D("xzo")->Draw("colz"); |
||
| 40 | f9daq | 78 | //h_correctedTDC->SetTitle("; t [ns]; Channel"); |
| 35 | f9daq | 79 | h_correctedTDC->Draw("colz"); |
| 80 | /* |
||
| 81 | TH1F* tdc1 = new TH1F("tdc1",";TDC [ns];Events",binsZ, minZ, maxZ); |
||
| 82 | //for(int j=minY; j<maxY; j++) { |
||
| 83 | for(int k=minZ; k<maxZ; k++) { |
||
| 84 | double signal = h0->GetBinContent(1,1,k); |
||
| 85 | tdc1->Fill(k,signal); |
||
| 86 | //} |
||
| 87 | } |
||
| 88 | tdc1->Draw(); |
||
| 89 | */ |
||
| 90 | |||
| 91 | TH2D* h_timeWalk = (TH2D*) h0->Project3D("zy"); |
||
| 92 | canvas21->cd(2); |
||
| 40 | f9daq | 93 | //h_timeWalk->SetTitle(";Threshold [V]; t [ns]"); |
| 35 | f9daq | 94 | h_timeWalk->Draw("colz"); |
| 95 | } |
||
| 96 | |||
| 36 | f9daq | 97 | if (strstr(plopt, "bias") != NULL) { |
| 98 | TCanvas* canvas01 = new TCanvas("canvas01","",800, 800); |
||
| 99 | canvas01->cd(); |
||
| 100 | TH2F* h_bias = (TH2F*) rootfile->Get("h_bias"); |
||
| 101 | h_bias->Draw("colz"); |
||
| 50 | f9daq | 102 | |
| 103 | TCanvas* canvas02 = new TCanvas("canvas02","",800, 800); |
||
| 104 | //canvas02->Divide(2,2); |
||
| 105 | canvas02->cd(); |
||
| 106 | TH3F* h_biasThreshold = (TH3F*) rootfile->Get("h_biasThreshold"); |
||
| 107 | //h_biasThreshold->Draw(); |
||
| 108 | //canvas02->cd(2); |
||
| 109 | //h_biasThreshold->GetXaxis()->SetRange(1,5); |
||
| 110 | h_biasThreshold->Project3D("zy"); |
||
| 91 | f9daq | 111 | h_biasThreshold->SetTitle(""); |
| 50 | f9daq | 112 | h_biasThreshold_zy->DrawCopy("colz"); |
| 113 | |||
| 114 | TCanvas* canvas03 = new TCanvas("canvas03","",1200, 1200); |
||
| 115 | canvas03->Divide(4,4); |
||
| 116 | int x1=0; |
||
| 117 | int x2=3; |
||
| 118 | int y1=0; |
||
| 119 | int y2=3; |
||
| 120 | for(int j=x1; j<=x2; j++){ |
||
| 121 | for(int k=y1; k<=y2; k++){ |
||
| 122 | int channel = map[j][k]; |
||
| 123 | int canvasPosition = 4*(4-k-1) + 4 - j; |
||
| 124 | canvas03->cd(canvasPosition); |
||
| 125 | h_biasThreshold->GetXaxis()->SetRange(channel+1,channel+1); |
||
| 126 | char histoName[128]; |
||
| 127 | sprintf(histoName, "Ch %d", channel); |
||
| 128 | h_biasThreshold->Project3D("zy"); |
||
| 129 | //sprintf(option, "h_biasThresholsd_zy%f", channel); |
||
| 130 | h_biasThreshold_zy->SetTitle(histoName); |
||
| 131 | h_biasThreshold_zy->DrawCopy("colz"); |
||
| 132 | |||
| 133 | } |
||
| 134 | } |
||
| 135 | TCanvas* canvas04 = new TCanvas("canvas04","",1200, 1200); |
||
| 136 | canvas04->Divide(4,4); |
||
| 137 | int x1=0; |
||
| 138 | int x2=3; |
||
| 139 | int y1=4; |
||
| 140 | int y2=7; |
||
| 141 | for(int j=x1; j<=x2; j++){ |
||
| 142 | for(int k=y1; k<=y2; k++){ |
||
| 143 | int channel = map[j][k]; |
||
| 144 | int canvasPosition = 4*(8-k-1) + 4 - j; |
||
| 145 | canvas04->cd(canvasPosition); |
||
| 146 | h_biasThreshold->GetXaxis()->SetRange(channel+1,channel+1); |
||
| 147 | char histoName[128]; |
||
| 148 | sprintf(histoName, "Ch %d", channel); |
||
| 149 | h_biasThreshold->Project3D("zy"); |
||
| 150 | //sprintf(option, "h_biasThresholsd_zy%f", channel); |
||
| 151 | h_biasThreshold_zy->SetTitle(histoName); |
||
| 152 | h_biasThreshold_zy->DrawCopy("colz"); |
||
| 153 | |||
| 154 | } |
||
| 155 | } |
||
| 156 | TCanvas* canvas05 = new TCanvas("canvas05","",1200, 1200); |
||
| 157 | canvas05->Divide(4,4); |
||
| 158 | int x1=4; |
||
| 159 | int x2=7; |
||
| 160 | int y1=0; |
||
| 161 | int y2=3; |
||
| 162 | for(int j=x1; j<=x2; j++){ |
||
| 163 | for(int k=y1; k<=y2; k++){ |
||
| 164 | int channel = map[j][k]; |
||
| 165 | int canvasPosition = 4*(4-k-1) + 8 - j; |
||
| 166 | canvas05->cd(canvasPosition); |
||
| 167 | h_biasThreshold->GetXaxis()->SetRange(channel+1,channel+1); |
||
| 168 | char histoName[128]; |
||
| 169 | sprintf(histoName, "Ch %d", channel); |
||
| 170 | h_biasThreshold->Project3D("zy"); |
||
| 171 | //sprintf(option, "h_biasThresholsd_zy%f", channel); |
||
| 172 | h_biasThreshold_zy->SetTitle(histoName); |
||
| 173 | h_biasThreshold_zy->DrawCopy("colz"); |
||
| 174 | |||
| 175 | } |
||
| 176 | } |
||
| 177 | |||
| 178 | TCanvas* canvas06 = new TCanvas("canvas06","",1200, 1200); |
||
| 179 | canvas06->Divide(4,4); |
||
| 180 | int x1=4; |
||
| 181 | int x2=7; |
||
| 182 | int y1=4; |
||
| 183 | int y2=7; |
||
| 184 | for(int j=x1; j<=x2; j++){ |
||
| 185 | for(int k=y1; k<=y2; k++){ |
||
| 186 | int channel = map[j][k]; |
||
| 187 | int canvasPosition = 4*(8-k-1) + 8 - j; |
||
| 188 | canvas06->cd(canvasPosition); |
||
| 189 | h_biasThreshold->GetXaxis()->SetRange(channel+1,channel+1); |
||
| 190 | char histoName[128]; |
||
| 191 | sprintf(histoName, "Ch %d", channel); |
||
| 192 | h_biasThreshold->Project3D("zy"); |
||
| 193 | //sprintf(option, "h_biasThresholsd_zy%f", channel); |
||
| 194 | h_biasThreshold_zy->SetTitle(histoName); |
||
| 195 | h_biasThreshold_zy->DrawCopy("colz"); |
||
| 196 | |||
| 197 | } |
||
| 198 | } |
||
| 199 | |||
| 36 | f9daq | 200 | } |
| 201 | |||
| 202 | if (strstr(plopt, "th") != NULL) { |
||
| 203 | TCanvas* canvas1 = new TCanvas("canvas1","",800,800); |
||
| 204 | canvas1->cd(); |
||
| 205 | //gPad->SetLogz(); |
||
| 206 | //canvas1->SetLogz(); |
||
| 207 | TH2F* h_threshold = (TH2F*) rootfile->Get("h_threshold"); |
||
| 208 | h_threshold->SetTitleOffset(1.5,"y"); |
||
| 209 | //h_threshold->GetZaxis()->SetRangeUser(0,1000); |
||
| 210 | h_threshold->SetContour(20); |
||
| 211 | |||
| 212 | h_threshold->Draw("colz"); |
||
| 213 | |||
| 214 | TCanvas* canvas2 = new TCanvas("canvas2","",1600,800); |
||
| 215 | canvas2->Divide(2); |
||
| 216 | TH1D* h_projection1 = h_threshold->ProjectionY("Ch 37",38,38); |
||
| 217 | canvas2->cd(1); |
||
| 218 | //h_projection1->GetYaxis()->SetRangeUser(0,1000); |
||
| 219 | //gPad->SetLogy(); |
||
| 220 | h_projection1->Draw(); |
||
| 221 | |||
| 222 | TH1D* h_projection2 = h_threshold->ProjectionY("Ch 38",39,39); |
||
| 223 | canvas2->cd(2); |
||
| 224 | //gPad->SetLogy(); |
||
| 91 | f9daq | 225 | h_projection2->Draw(); |
| 226 | |||
| 227 | TCanvas* canvas21 = new TCanvas("canvas21","",1200,1200); |
||
| 228 | canvas21->Divide(3,3); |
||
| 229 | TH3F* h_biasThreshold = (TH3F*) rootfile->Get("h_biasThreshold"); |
||
| 230 | //h_biasThreshold->Draw(); |
||
| 231 | //canvas02->cd(2); |
||
| 232 | //h_biasThreshold->GetXaxis()->SetRange(1,5); |
||
| 233 | h_biasThreshold->Project3D("zy"); |
||
| 234 | canvas21->cd(1); |
||
| 235 | h_biasThreshold_zy->DrawCopy("colz"); |
||
| 236 | canvas21->cd(2); |
||
| 237 | for (int jk=0; jk<9; jk++){ |
||
| 238 | h_biasThreshold_zy->ProjectionY("_py",jk*2+1,jk*2+1); |
||
| 239 | canvas21->cd(jk+1); |
||
| 240 | char histoName[128]; |
||
| 241 | double bias = h_biasThreshold_zy->GetXaxis()->GetBinCenter(jk*2+1); |
||
| 242 | sprintf(histoName, "Bias = %4.2f V", bias); |
||
| 243 | h_biasThreshold_zy_py->SetTitle(histoName); |
||
| 244 | h_biasThreshold_zy_py->DrawCopy(); |
||
| 245 | } |
||
| 246 | |||
| 36 | f9daq | 247 | } |
| 248 | |||
| 35 | f9daq | 249 | if(strstr(plopt,"16") != NULL) |
| 250 | { |
||
| 36 | f9daq | 251 | TH2F* h_threshold = (TH2F*) rootfile->Get("h_bias"); |
| 35 | f9daq | 252 | TCanvas* canvas = new TCanvas("canvas","",1600,1600); |
| 38 | f9daq | 253 | TCanvas* canvas11 = new TCanvas("canvas11","",800,800); |
| 35 | f9daq | 254 | canvas->cd(); |
| 255 | TVirtualPad *main = new TPad("main","main",0,0,1,1,10,1); |
||
| 256 | main->Draw(); |
||
| 38 | f9daq | 257 | main->cd(); |
| 35 | f9daq | 258 | main->Divide(4,4); |
| 259 | for(int i=chXstart; i<chXend; i++) { |
||
| 40 | f9daq | 260 | for(int j=(int)par1; j<(int)par2; j++) { |
| 35 | f9daq | 261 | int channel = map[i][j]; |
| 262 | TH1D* h_projection = h_threshold->ProjectionY("",channel+1,channel+1); |
||
| 40 | f9daq | 263 | int canvasPosition = i-chXend+4*((int)par2-j)+1; |
| 35 | f9daq | 264 | printf(" %d ", canvasPosition); |
| 265 | main->cd(canvasPosition); |
||
| 38 | f9daq | 266 | gPad->SetMargin(0.08, 0.08, 0.08, 0.08); |
| 35 | f9daq | 267 | char name[128]; |
| 268 | sprintf(name,"Channel %d",channel); |
||
| 269 | h_projection->SetTitle(name); |
||
| 38 | f9daq | 270 | h_projection->GetYaxis()->SetRangeUser(0,10000); |
| 271 | h_projection->DrawCopy("pe1x0"); |
||
| 272 | canvas11->cd(); |
||
| 273 | if (canvasPosition == 1) { |
||
| 274 | h_projection->SetMarkerStyle(kFullDotMedium); |
||
| 275 | h_projection->GetYaxis()->SetRangeUser(0,10000); |
||
| 276 | h_projection->DrawCopy("Pl"); |
||
| 277 | } |
||
| 278 | else { |
||
| 279 | h_projection->SetMarkerStyle(kFullDotMedium); |
||
| 280 | h_projection->SetLineColor(canvasPosition); |
||
| 281 | h_projection->SetMarkerColor(canvasPosition); |
||
| 282 | h_projection->SetTitle(""); |
||
| 283 | h_projection->DrawCopy("plsame"); |
||
| 284 | } |
||
| 35 | f9daq | 285 | } |
| 38 | f9daq | 286 | |
| 35 | f9daq | 287 | } |
| 288 | printf("\n"); |
||
| 289 | } |
||
| 290 | |||
| 36 | f9daq | 291 | if(strstr(plopt, "1ch") != NULL) |
| 292 | { |
||
| 293 | TH2F* h_bias = (TH2F*) rootfile->Get("h_bias"); |
||
| 38 | f9daq | 294 | TCanvas* canvas = new TCanvas("canvas","",800,800); |
| 36 | f9daq | 295 | canvas->cd(); |
| 296 | int channel = map[chXstart][chXend]; |
||
| 297 | TH1D* h_projection = h_bias->ProjectionY("",channel+1,channel+1); |
||
| 298 | char name[128]; |
||
| 299 | sprintf(name,"Channel %d",channel); |
||
| 300 | h_projection->SetTitle(name); |
||
| 38 | f9daq | 301 | h_projection->SetMarkerStyle(kFullDotMedium); |
| 302 | //h_projection->SetMarkerSize(8); |
||
| 303 | gPad->SetMargin(0.1,0.1,0.1,0.1); |
||
| 40 | f9daq | 304 | |
| 305 | |||
| 306 | TF1* f_linear = new TF1("f_linear", "[0] + [1]*x", par1, par2); |
||
| 307 | f_linear->SetParameter(0, -1); |
||
| 308 | f_linear->SetParameter(1, 0.1); |
||
| 309 | f_linear->SetParNames("p0", "p1"); |
||
| 310 | h_projection->Fit(f_linear,"q"); |
||
| 311 | h_projection->Fit(f_linear,"r"); |
||
| 38 | f9daq | 312 | h_projection->Draw("pe1x0"); |
| 40 | f9daq | 313 | |
| 314 | double p0 = f_linear->GetParameter(0); |
||
| 315 | double p1 = f_linear->GetParameter(1); |
||
| 316 | printf("p0 %f p1 %f V_b = %f [V]\n", p0, p1, -p0/p1); |
||
| 317 | TPaveText *text = new TPaveText(0.15,0.5,0.35,0.56); |
||
| 318 | text->SetFillColor(0); |
||
| 319 | text->SetBorderSize(0.1); |
||
| 320 | text->SetTextSize(0.030); |
||
| 321 | char string[128]; |
||
| 322 | sprintf(string, "V_b = %2.2f", -p0/p1); |
||
| 323 | text->AddText(string); |
||
| 324 | text->Draw("same"); |
||
| 325 | text->Paint(); |
||
| 36 | f9daq | 326 | } |
| 327 | |||
| 35 | f9daq | 328 | if( strstr(plopt, "all") != NULL ) { |
| 329 | TCanvas *canvas2 = new TCanvas("canvas2","Hits x;;",2000,2000); |
||
| 330 | TCanvas *canvas3 = new TCanvas("canvas3","Hits y;;",2000,2000); |
||
| 331 | canvas2->Divide(8,8); |
||
| 332 | canvas3->Divide(8,8); |
||
| 333 | TH1F* h_hitsx; |
||
| 334 | TH1F* h_hitsy; |
||
| 335 | for(int i=0; i<c_nChannels; i++) { |
||
| 336 | canvas2->cd(i+1); |
||
| 337 | char hname[128]; |
||
| 338 | sprintf(hname, "hnhitsx%d", i); |
||
| 339 | h_hitsx = (TH1F*)rootfile->Get(hname); |
||
| 340 | h_hitsx->Draw(); |
||
| 341 | canvas3->cd(i+1); |
||
| 342 | sprintf(hname, "hnhitsy%d", i); |
||
| 343 | h_hitsy = (TH1F*)rootfile->Get(hname); |
||
| 344 | h_hitsy->Draw(); |
||
| 345 | } |
||
| 346 | } |
||
| 347 | |||
| 348 | if( strstr(plopt, "x") != NULL ) { |
||
| 349 | TCanvas *canvas10 = new TCanvas("canvas10","Ch x;;",500,500); |
||
| 350 | TH1F* h_hitsx; |
||
| 351 | canvas10->cd(); |
||
| 352 | char hname[128]; |
||
| 353 | sprintf(hname, "hnhitsx%d", chXstart); |
||
| 354 | h_hitsx = (TH1F*)rootfile->Get(hname); |
||
| 355 | h_hitsx->Draw("colz"); |
||
| 356 | } |
||
| 357 | |||
| 358 | if( strstr(plopt, "y") != NULL ) { |
||
| 359 | TCanvas *canvas11 = new TCanvas("canvas11","Ch x;;",500,500); |
||
| 360 | TH1F* h_hitsy; |
||
| 361 | canvas11->cd(); |
||
| 362 | char hname[128]; |
||
| 363 | sprintf(hname, "hnhitsy%d", chXstart); |
||
| 364 | h_hitsy = (TH1F*)rootfile->Get(hname); |
||
| 365 | h_hitsy->Draw(); |
||
| 366 | } |
||
| 40 | f9daq | 367 | |
| 35 | f9daq | 368 | |
| 369 | /** Draws the sum of the channels |
||
| 370 | * Each channel is a 2d plot |
||
| 371 | * Intended for the study of 1 channel |
||
| 372 | */ |
||
| 373 | if (strstr(plopt, "2d") != NULL) { |
||
| 374 | |||
| 375 | int nX = chXend - chXstart + 1; |
||
| 40 | f9daq | 376 | int nY = (int)par2 - (int)par1 + 1; |
| 35 | f9daq | 377 | TCanvas* canvas7 = new TCanvas("canvas7","canvas7", nX*400,nY*400); |
| 378 | printf("nx %d ny %d\n",nX,nY); |
||
| 379 | canvas7->Divide(nX,nY); |
||
| 40 | f9daq | 380 | for(int i=(int)par1; i<=(int)par2; i++) { |
| 35 | f9daq | 381 | for(int j=chXstart; j<=chXend; j++) { |
| 382 | //int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
||
| 40 | f9daq | 383 | int canvasPosition = nX*((int)par2-i)+(j-chXstart)+1; |
| 35 | f9daq | 384 | if (debug) printf("canvas %d\n",canvasPosition); |
| 385 | canvas7->cd(canvasPosition); |
||
| 386 | char hname[128]; |
||
| 387 | int chPosition = map[j][i]; |
||
| 388 | sprintf(hname, "h2d%d", chPosition); |
||
| 389 | TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
||
| 390 | h_2d->Draw("colz"); |
||
| 391 | } //x |
||
| 392 | } |
||
| 393 | |||
| 394 | // Number of photoelectrons - Poissonian correction |
||
| 395 | TCanvas* canvas8 = new TCanvas("canvas8","canvas8", 1000,1000); |
||
| 396 | canvas8->cd(); |
||
| 397 | gStyle->SetOptStat(0); |
||
| 398 | char hname[128]; |
||
| 399 | int chPosition = map[1][2]; |
||
| 400 | sprintf(hname, "h2d%d", chPosition); |
||
| 401 | TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
||
| 402 | |||
| 403 | Int_t binsX = h_2d->GetXaxis()->GetNbins(); |
||
| 404 | Int_t minX = h_2d->GetXaxis()->GetFirst(); |
||
| 405 | Int_t maxX = h_2d->GetXaxis()->GetLast()+1; |
||
| 406 | Int_t binsY = h_2d->GetYaxis()->GetNbins(); |
||
| 407 | Int_t minY = h_2d->GetYaxis()->GetFirst(); |
||
| 408 | Int_t maxY = h_2d->GetYaxis()->GetLast()+1; |
||
| 409 | Double_t xLowUser = h_2d->GetXaxis()->GetBinLowEdge(minX); |
||
| 410 | Double_t xUpUser = h_2d->GetXaxis()->GetBinUpEdge(maxX); |
||
| 411 | Double_t yLowUser = h_2d->GetYaxis()->GetBinLowEdge(minY); |
||
| 412 | Double_t yUpUser = h_2d->GetYaxis()->GetBinUpEdge(maxY); |
||
| 413 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
||
| 414 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, xLowUser, xUpUser, binsY, yLowUser, yUpUser); |
||
| 415 | |||
| 416 | double noise = getNoise(h_2d, 1, 89); |
||
| 417 | if(debug) printf("Noise = %f\n", noise); |
||
| 418 | for(int k=minX; k<=maxX; k++) { |
||
| 419 | for(int j=minY; j<=maxY; j++) { |
||
| 420 | double signal = h_2d->GetBinContent(k,j); // detected |
||
| 421 | //double p = ((signal - noise) > 1) ? (signal-noise) : 1; |
||
| 422 | double p = signal - noise; |
||
| 423 | p /= 10000.; |
||
| 424 | double p0 = 1.0 - p; // events with zero photons |
||
| 425 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
||
| 426 | double eta = -log(p0); // constant of the poissonian statistics |
||
| 427 | if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
||
| 428 | double x = k*(xUpUser-xLowUser)/double(binsX); |
||
| 429 | double y = j*(yUpUser-yLowUser)/double(binsY); |
||
| 430 | h_corrected->Fill(x+3,y+8, eta); |
||
| 431 | } |
||
| 432 | } |
||
| 433 | h_corrected->SetTitle("n_pe = - ln P(0);x[mm];y[mm]"); |
||
| 434 | h_corrected->GetZaxis()->SetRangeUser(-0.05,0.3); |
||
| 435 | h_corrected->Draw("colz"); |
||
| 436 | |||
| 437 | // collection efficiency |
||
| 438 | int nPoints =0; |
||
| 439 | double efficiency=0; |
||
| 440 | for (int i=18; i<=58; i++) { |
||
| 441 | for (int j=19; j<=59; j++) { |
||
| 442 | double signal = h_corrected->GetBinContent(i,j); |
||
| 443 | if(debug) printf("signal %f\n",signal); |
||
| 444 | efficiency += signal; |
||
| 445 | nPoints++; |
||
| 446 | } |
||
| 447 | } |
||
| 448 | printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
||
| 449 | } |
||
| 450 | |||
| 451 | /** Draws the sum of channel signals |
||
| 452 | * Each channel is a 2d ('h2d') histogram |
||
| 453 | * Suitable for 8x8 chs scan |
||
| 454 | */ |
||
| 455 | if( strstr(plopt, "sum") != NULL ) { |
||
| 40 | f9daq | 456 | int nX = 7 - 0 + 1; |
| 457 | int nY = 7 - 0 + 1; |
||
| 35 | f9daq | 458 | TCanvas* canvas12 = new TCanvas("canvas12","c2",8*200, 8*200); |
| 459 | canvas12->cd(); |
||
| 460 | gStyle->SetOptStat(0); |
||
| 461 | |||
| 462 | // final histogram parameters |
||
| 463 | TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
||
| 464 | Int_t binsX = h0->GetXaxis()->GetNbins(); |
||
| 465 | Int_t minX = h0->GetXaxis()->GetFirst(); |
||
| 466 | Int_t maxX = h0->GetXaxis()->GetLast()+1; |
||
| 467 | Int_t binsY = h0->GetYaxis()->GetNbins(); |
||
| 468 | Int_t minY = h0->GetYaxis()->GetFirst(); |
||
| 469 | Int_t maxY = h0->GetYaxis()->GetLast()+1; |
||
| 470 | Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
||
| 471 | Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
||
| 472 | Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
||
| 473 | Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
||
| 474 | if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
||
| 475 | TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
||
| 476 | TH2F* h[9]; |
||
| 477 | |||
| 478 | // 2d histogram noise subtraction and poisson scaling |
||
| 40 | f9daq | 479 | for(int i=0; i<=7; i++) { |
| 35 | f9daq | 480 | for(int j=chXstart; j<=chXend; j++) { |
| 481 | int chPosition = map[j][i]; |
||
| 482 | char hname[128]; |
||
| 483 | sprintf(hname, "h2d%d", chPosition); |
||
| 40 | f9daq | 484 | int histogram = nX*(i-0)+(j-chXstart); |
| 35 | f9daq | 485 | h[histogram] = (TH2F *) rootfile->Get(hname); |
| 486 | int noise = getNoise(h[histogram], 1, 170); |
||
| 487 | if (debug) printf("noise: %d\n",noise); |
||
| 488 | for(int k=minX; k<=maxX; k++) { |
||
| 489 | for(int l=minY; l<=maxY; l++) { |
||
| 490 | int signal = h[histogram]->GetBinContent(k,l); // detected |
||
| 491 | //double p = ((signal - noise) > 0.1) ? (signal-noise) : 0.1; |
||
| 492 | double p = signal - noise; |
||
| 493 | p /= 10000.; |
||
| 494 | double p0 = 1.0 - p; // events with zero photons |
||
| 495 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
||
| 496 | double eta = -log(p0); |
||
| 497 | //printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
||
| 498 | double x = k*(xUpUser-xLowUser)/double(binsX); |
||
| 499 | double y = l*(yUpUser-yLowUser)/double(binsY); |
||
| 500 | h_corrected->Fill(x-c_xOffset, y-c_yOffset, signal); |
||
| 501 | } |
||
| 502 | } |
||
| 503 | |||
| 504 | } |
||
| 505 | } |
||
| 506 | |||
| 507 | h_corrected->SetTitle("SiPM#2 n_p.e.;x[mm];y[mm]"); |
||
| 508 | //h_corrected->GetZaxis()->SetRangeUser(-0.05,.30); |
||
| 509 | h_corrected->Draw("colz"); |
||
| 510 | |||
| 511 | TCanvas* canvas13 = new TCanvas("canvas13","canvas13",600,300); |
||
| 512 | canvas13->Divide(2); |
||
| 513 | canvas13->cd(1); |
||
| 514 | //h[16]->Draw("colz"); |
||
| 515 | canvas13->cd(2); |
||
| 516 | //h[8]->Draw("colz"); |
||
| 517 | } |
||
| 518 | |||
| 519 | return(0); |
||
| 520 | } |
||
| 521 | |||
| 522 | /** Function calculates the noise from one channel |
||
| 523 | * it runs through the bins along x and returns the average value |
||
| 524 | */ |
||
| 525 | double getNoise(TH2F* histogram, int yStart, int yEnd) |
||
| 526 | { |
||
| 527 | double noise=0; |
||
| 528 | int count=0; |
||
| 529 | for(int j=yStart; j<yEnd; j++) { |
||
| 530 | double value = histogram->GetBinContent(j,2); |
||
| 531 | //if (noise < value) noise = value; |
||
| 532 | noise += value; |
||
| 533 | count++; |
||
| 534 | } |
||
| 535 | return (noise/double(count)); |
||
| 536 | } |
||
| 537 |