Rev 50 | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed
| Rev 50 | Rev 91 | ||
|---|---|---|---|
| Line -... | Line 1... | ||
| - | 1 | // ***************************************************************************** |
|
| - | 2 | // Jozef Stefan Institute |
|
| - | 3 | // Dino Tahirovic |
|
| - | 4 | // 2012-2014 |
|
| - | 5 | // |
|
| - | 6 | // Options: |
|
| - | 7 | // 2d - 1 channel analysis |
|
| - | 8 | // s - save an image (png) after drawinv the canvas |
|
| - | 9 | // p - save every canvas as page in pdf |
|
| - | 10 | // |
|
| - | 11 | // ***************************************************************************** |
|
| - | 12 | ||
| 1 | #include "TROOT.h" |
13 | #include "TROOT.h" |
| 2 | #include "TFile.h" |
14 | #include "TFile.h" |
| 3 | #include "TBenchmark.h" |
15 | #include "TBenchmark.h" |
| 4 | #include "TH1F.h" |
16 | #include "TH1F.h" |
| 5 | #include "TH2F.h" |
17 | #include "TH2F.h" |
| Line 350... | Line 362... | ||
| 350 | if ( !((i == parameter2) and (j == parameter1)) ) h_2d->GetZaxis()->SetRangeUser(0,200); |
362 | if ( !((i == parameter2) and (j == parameter1)) ) h_2d->GetZaxis()->SetRangeUser(0,200); |
| 351 | } //x |
363 | } //x |
| 352 | } |
364 | } |
| 353 | 365 | ||
| 354 | // Number of photoelectrons - Poissonian correction |
366 | // Number of photoelectrons - Poissonian correction |
| 355 | TCanvas* canvas8 = new TCanvas("canvas8","canvas8", |
367 | TCanvas* canvas8 = new TCanvas("canvas8","canvas8", 1600,1600); |
| 356 | canvas8->cd(); |
368 | canvas8->cd(); |
| 357 | gStyle->SetOptStat(0); |
369 | gStyle->SetOptStat(0); |
| 358 | char hname[128]; |
370 | char hname[128]; |
| 359 | int chPosition = map[parameter1][parameter2]; |
371 | int chPosition = map[parameter1][parameter2]; |
| 360 | sprintf(hname, "h2d%d", chPosition); |
372 | sprintf(hname, "h2d%d", chPosition); |
| Line 376... | Line 388... | ||
| 376 | double diffY = yUpUser - yLowUser; |
388 | double diffY = yUpUser - yLowUser; |
| 377 | //TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, 3.00, 3.00+diffX, binsY, 8.00, 8.00+diffY); |
389 | //TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, 3.00, 3.00+diffX, binsY, 8.00, 8.00+diffY); |
| 378 | //TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, minX, maxX, binsY, minY, maxY); |
390 | //TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, minX, maxX, binsY, minY, maxY); |
| 379 | 391 | ||
| 380 | double noise = getNoise(h_2d, 1, 70); |
392 | double noise = getNoise(h_2d, 1, 70); |
| 381 | noise = |
393 | noise = 61; |
| 382 |
|
394 | printf("[2d scan] Average noise = %f\n", noise); |
| 383 | for(int k=minX; k<=maxX; k++) { |
395 | for(int k=minX; k<=maxX; k++) { |
| 384 | for(int j=minY; j<=maxY; j++) { |
396 | for(int j=minY; j<=maxY; j++) { |
| 385 | double signal = h_2d->GetBinContent(k,j); // detected |
397 | double signal = h_2d->GetBinContent(k,j); // detected |
| 386 | //double p = ((signal - noise) > 1) ? (signal-noise) : 1; |
- | |
| 387 | double p = signal - noise; |
398 | double p = signal - noise; |
| 388 | p /= 10000.; |
399 | p /= 10000.; |
| 389 | double p0 = 1.0 - p; // events with zero photons |
400 | double p0 = 1.0 - p; // events with zero photons |
| 390 | //double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
- | |
| 391 | double eta = -log(p0); // constant of the poissonian statistics |
401 | double eta = -log(p0); // constant of the poissonian statistics |
| 392 | if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
402 | if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
| 393 | //double x = xLowUser + k*(xUpUser - xLowUser) / double(binsX); |
- | |
| 394 | double x = h_2d->GetXaxis()->GetBinCenter(k); |
403 | double x = h_2d->GetXaxis()->GetBinCenter(k); |
| 395 | //double y = yLowUser + j*(yUpUser-yLowUser)/double(binsY); |
- | |
| 396 | double y = h_2d->GetYaxis()->GetBinCenter(j); |
404 | double y = h_2d->GetYaxis()->GetBinCenter(j); |
| 397 | if (debug) printf("x=%f y=%f\n", x, y); |
405 | if (debug) printf("x=%f y=%f\n", x, y); |
| 398 | h_corrected->Fill(x, y, |
406 | h_corrected->Fill(x, y, eta); |
| 399 | } |
407 | } |
| 400 | } |
408 | } |
| 401 | //h_corrected->SetTitle("n_pe = - ln P(0);x[mm];y[mm]"); |
409 | //h_corrected->SetTitle("n_pe = - ln P(0);x[mm];y[mm]"); |
| 402 | h_corrected->SetTitle(";x[mm];y[mm]"); |
410 | h_corrected->SetTitle(";x[mm];y[mm]"); |
| 403 | //gStyle->SetPalette(52,0); // black and white for print |
411 | //gStyle->SetPalette(52,0); // black and white for print |
| 404 | //gStyle->SetPalette(1); //black and white 2nd option |
412 | //gStyle->SetPalette(1); //black and white 2nd option |
| 405 | //SetGS(); // inverse |
413 | //SetGS(); // inverse |
| 406 | h_corrected |
414 | //h_corrected->GetZaxis()->SetRangeUser(-0.01,0.12); |
| 407 | h_corrected->SetContour( |
415 | h_corrected->SetContour(50); |
| 408 | h_corrected->Draw("colz"); |
416 | h_corrected->Draw("colz"); |
| 409 | 417 | ||
| 410 | // collection efficiency |
418 | // collection efficiency |
| 411 | int nPoints =0; |
419 | int nPoints =0; |
| 412 | double efficiency=0; |
420 | double efficiency=0; |
| 413 | for (int i= |
421 | for (int i=minX; i<=maxX; i++) { |
| 414 | for (int j= |
422 | for (int j=minY; j<=maxY; j++) { |
| 415 | double signal = h_corrected->GetBinContent(i,j); |
423 | double signal = h_corrected->GetBinContent(i,j); |
| 416 | if(debug) printf("signal %f\n",signal); |
424 | if(debug) printf("signal %f\n",signal); |
| 417 | efficiency += signal; |
425 | if (signal > 0.02) efficiency += signal; |
| 418 | nPoints++; |
426 | nPoints++; |
| 419 | } |
427 | } |
| 420 | } |
428 | } |
| 421 | printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
429 | printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
| 422 | 430 | ||
| 423 | if(strstr(plopt, " |
431 | if(strstr(plopt, "s") != NULL) { |
| - | 432 | char photo[128]; |
|
| - | 433 | sprintf(photo, "ps/%s.png", filename); |
|
| 424 | canvas8 |
434 | //canvas8->Print("result.pdf"); |
| - | 435 | canvas8->SaveAs(photo); |
|
| 425 | } |
436 | } |
| 426 | 437 | ||
| 427 | } |
438 | } |
| 428 | 439 | ||
| 429 | /** Draws the sum of channel signals |
440 | /** Draws the sum of channel signals |
| Line 549... | Line 560... | ||
| 549 | double getNoise(TH2F* histogram, int yStart, int yEnd) |
560 | double getNoise(TH2F* histogram, int yStart, int yEnd) |
| 550 | { |
561 | { |
| 551 | double noise=0; |
562 | double noise=0; |
| 552 | int count=0; |
563 | int count=0; |
| 553 | for(int j=yStart; j<yEnd; j++) { |
564 | for(int j=yStart; j<yEnd; j++) { |
| 554 | double value = histogram->GetBinContent(j, |
565 | double value = histogram->GetBinContent(j,1); |
| 555 | //if (noise < value) noise = value; |
566 | //if (noise < value) noise = value; |
| 556 | noise += value; |
567 | noise += value; |
| 557 | count++; |
568 | count++; |
| 558 | } |
569 | } |
| 559 | return (noise/double(count)); |
570 | return (noise/double(count)); |