/sipmScan/examples/thresholdScan.C |
---|
0,0 → 1,513 |
#include "TROOT.h" |
#include "TFile.h" |
#include "TBenchmark.h" |
#include "TH1F.h" |
#include "TH2F.h" |
#include "TH3F.h" |
#include "TCanvas.h" |
#include "TStyle.h" |
#include "TPad.h" |
#include "TF1.h" |
#include "TGraph.h" |
#include "TSpectrum.h" |
#include "stdio.h" |
#include "include/RTUtil.h" |
double getNoise(TH2F*, int, int); |
int thresholdScan(char filename[256] = "test", char plopt[256]="th", int chXstart=0, int chXend=7, int chYstart=0, int chYend=7, bool debug = false) |
{ |
const int c_nChannels = 64; |
const double c_xOffset = 1; // mm |
const double c_yOffset = 0.7; |
int map[8][8]={{32,34,53,55,40,42,61,63}, |
{48,50,37,39,56,58,45,47}, |
{33,35,52,54,41,43,60,62}, |
{49,51,36,38,57,59,44,46}, |
{17,19,4,6,25,27,12,14}, |
{1,3,20,22,9,11,28,30}, |
{16,18,5,7,24,26,13,15}, |
{0,2,21,23,8,10,29,31} |
}; |
RTSetStyle(gStyle); |
// open the file with histograms |
char fnameroot[256]; |
TFile* rootfile; |
sprintf(fnameroot, "root/%s.root", filename); |
rootfile = (TFile *) gROOT->FindObject(filename); |
if(rootfile==NULL) rootfile = new TFile(fnameroot); |
if(rootfile==NULL) { |
printf("Cannot open root file %s!!!\n",fnameroot); |
return(0); |
} |
if (strstr(plopt, "th") != NULL) { |
TCanvas* canvas1 = new TCanvas("canvas1","",800,800); |
canvas1->cd(); |
//gPad->SetLogz(); |
//canvas1->SetLogz(); |
TH2F* h_threshold = (TH2F*) rootfile->Get("h_threshold"); |
h_threshold->SetTitleOffset(1.5,"y"); |
//h_threshold->GetZaxis()->SetRangeUser(0,1000); |
h_threshold->SetContour(20); |
h_threshold->Draw("colz"); |
TCanvas* canvas2 = new TCanvas("canvas2","",1600,800); |
canvas2->Divide(2); |
TH1D* h_projection1 = h_threshold->ProjectionY("Ch 37",38,38); |
canvas2->cd(1); |
//h_projection1->GetYaxis()->SetRangeUser(0,1000); |
//gPad->SetLogy(); |
h_projection1->Draw(); |
TH1D* h_projection2 = h_threshold->ProjectionY("Ch 38",39,39); |
canvas2->cd(2); |
//gPad->SetLogy(); |
h_projection2->Draw(); |
} |
if(strstr(plopt, "tdc") != NULL) { |
TCanvas *canvas21 = new TCanvas("canvas21","canvas21",1600,800); |
TH3F* h0 = (TH3F*) rootfile->Get("h_correctedTDC"); |
Int_t binsX = h0->GetXaxis()->GetNbins(); |
Int_t minX = h0->GetXaxis()->GetFirst(); |
Int_t maxX = h0->GetXaxis()->GetLast()+1; |
Int_t binsY = h0->GetYaxis()->GetNbins(); |
Int_t minY = h0->GetYaxis()->GetFirst(); |
Int_t maxY = h0->GetYaxis()->GetLast()+1; |
Int_t binsZ = h0->GetZaxis()->GetNbins(); |
Int_t minZ = h0->GetZaxis()->GetFirst(); |
Int_t maxZ = h0->GetZaxis()->GetLast()+1; |
Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
canvas21->Divide(2); |
//h0->Draw(); |
canvas21->cd(1); |
//h0->Project3D("xy")->Draw("colz"); |
//canvas21->cd(2); |
//h0->Project3D("yz")->Draw("colz"); |
//canvas21->cd(3); |
TH2D* h_correctedTDC = (TH2D*) h0->Project3D("xz"); |
//h0->GetZaxis()->SetRangeUser(-5,5); |
//h0->Project3D("xzo")->Draw("colz"); |
h_correctedTDC->SetTitle("; t [ns]; Channel"); |
h_correctedTDC->Draw("colz"); |
/* |
TH1F* tdc1 = new TH1F("tdc1",";TDC [ns];Events",binsZ, minZ, maxZ); |
//for(int j=minY; j<maxY; j++) { |
for(int k=minZ; k<maxZ; k++) { |
double signal = h0->GetBinContent(1,1,k); |
tdc1->Fill(k,signal); |
//} |
} |
tdc1->Draw(); |
*/ |
TH2D* h_timeWalk = (TH2D*) h0->Project3D("zy"); |
canvas21->cd(2); |
h_timeWalk->SetTitle(";Threshold [V]; t [ns]"); |
h_timeWalk->Draw("colz"); |
} |
if(strstr(plopt,"16") != NULL) |
{ |
TH2F* h_threshold = (TH2F*) rootfile->Get("h_threshold"); |
TCanvas* canvas = new TCanvas("canvas","",1600,1600); |
canvas->cd(); |
TVirtualPad *main = new TPad("main","main",0,0,1,1,10,1); |
main->Draw(); |
main->Divide(4,4); |
for(int i=chXstart; i<chXend; i++) { |
for(int j=chYstart; j<chYend; j++) { |
int channel = map[i][j]; |
TH1D* h_projection = h_threshold->ProjectionY("",channel+1,channel+1); |
int canvasPosition = i-chXend+4*(chYend-j)+1; |
printf(" %d ", canvasPosition); |
main->cd(canvasPosition); |
char name[128]; |
sprintf(name,"Channel %d",channel); |
h_projection->SetTitle(name); |
h_projection->DrawCopy(); |
} |
} |
printf("\n"); |
} |
if( strstr(plopt, "all") != NULL ) { |
TCanvas *canvas2 = new TCanvas("canvas2","Hits x;;",2000,2000); |
TCanvas *canvas3 = new TCanvas("canvas3","Hits y;;",2000,2000); |
canvas2->Divide(8,8); |
canvas3->Divide(8,8); |
TH1F* h_hitsx; |
TH1F* h_hitsy; |
for(int i=0; i<c_nChannels; i++) { |
canvas2->cd(i+1); |
char hname[128]; |
sprintf(hname, "hnhitsx%d", i); |
h_hitsx = (TH1F*)rootfile->Get(hname); |
h_hitsx->Draw(); |
canvas3->cd(i+1); |
sprintf(hname, "hnhitsy%d", i); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->Draw(); |
} |
} |
if( strstr(plopt, "x") != NULL ) { |
TCanvas *canvas10 = new TCanvas("canvas10","Ch x;;",500,500); |
TH1F* h_hitsx; |
canvas10->cd(); |
char hname[128]; |
sprintf(hname, "hnhitsx%d", chXstart); |
h_hitsx = (TH1F*)rootfile->Get(hname); |
h_hitsx->Draw("colz"); |
} |
if( strstr(plopt, "y") != NULL ) { |
TCanvas *canvas11 = new TCanvas("canvas11","Ch x;;",500,500); |
TH1F* h_hitsy; |
canvas11->cd(); |
char hname[128]; |
sprintf(hname, "hnhitsy%d", chXstart); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->Draw(); |
} |
if( strstr(plopt, "share") != NULL ) { |
/*TCanvas *canvas4 = new TCanvas("canvas1","canvas1",1000,1000); |
int nChannels = chYend-chYstart+1; |
int ncols = nChannels/2; |
printf("nch %d nch\\2 %d\n", nChannels, ncols); |
canvas4->Divide(2,ncols); |
TH1F* h_hitsy; |
for(int i=chYstart; i<=chYend; i++){ |
canvas4->cd(i-chYstart+1); |
char hname[128]; |
int chPosition = map[0][i]; |
sprintf(hname, "hnhitsy%d", chPosition); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->Draw(); |
}*/ |
TCanvas *canvas4 = new TCanvas("canvas4","canvas4",500,500); |
canvas4->cd(); |
for(int i=chXstart; i<=chXend; i++) { |
TH1F* h_hitsx; |
char hname[128]; |
int chPosition = map[i][chYstart]; |
sprintf(hname, "hnhitsx%d", chPosition); |
h_hitsx = (TH1F*)rootfile->Get(hname); |
h_hitsx->SetTitle("Scan X;x [mm]; Entries"); |
h_hitsx->GetYaxis()->SetTitleOffset(1.3); |
h_hitsx->SetStats(0); |
if (i == chXstart) |
h_hitsx->Draw(); |
else { |
h_hitsx->SetLineColor(i+1); |
h_hitsx->Draw("same"); |
} |
} |
//sprintf(fullname, "ps/%s_Yshare.eps", filename); |
//canvas4->SaveAs(fullname); |
TCanvas *canvas5 = new TCanvas("canvas5","canvas5",500,500); |
canvas5->cd(); |
for(int i=chYstart; i<=chYend; i++) { |
TH1F* h_hitsy; |
char hname[128]; |
int chPosition = map[chXstart][i]; |
sprintf(hname, "hnhitsy%d", chPosition); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->SetTitle("Scan Y;y [mm]; Entries"); |
h_hitsy->GetYaxis()->SetTitleOffset(1.3); |
h_hitsy->SetStats(0); |
if (i == chYstart) |
h_hitsy->Draw(); |
else { |
h_hitsy->SetLineColor(i+1); |
h_hitsy->Draw("same"); |
} |
} |
//sprintf(fullname, "ps/%s_Yshare.eps", filename); |
//canvas5->SaveAs(fullname); |
} |
/** Draws the signal from 8 channels in x-row |
* for one specific y bin, so the background and cross-talk |
* can be estimated. |
* Draws also a 2d scan of these channels. |
*/ |
if (strstr(plopt, "line") != NULL) { |
TCanvas* canvas6 = new TCanvas("canvas6","canvas6",500,500); |
canvas6->cd(0); |
gStyle->SetOptStat(0); |
TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
Int_t binsX = h0->GetXaxis()->GetNbins(); |
Int_t minX = h0->GetXaxis()->GetFirst(); |
Int_t maxX = h0->GetXaxis()->GetLast()+1; |
Int_t binsY = h0->GetYaxis()->GetNbins(); |
Int_t minY = h0->GetYaxis()->GetFirst(); |
Int_t maxY = h0->GetYaxis()->GetLast()+1; |
Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
//! 1-dimension position in x vs. hits |
TH2F* h[8]; |
TH1F* h_line[8]; |
for(int j=0; j<8; j++) { |
h_line[j] = new TH1F("h_line", "h_line", binsX, xLowUser, xUpUser); |
} |
for(int j=chXstart; j<=chXend; j++) { |
int chPosition = map[j][chYstart]; |
char hname[128]; |
sprintf(hname, "h2d%d", chPosition); |
int histogram = j; |
h[histogram] = (TH2F *) rootfile->Get(hname); |
int noise = getNoise(h[histogram], 1, 160); |
for(int k=minX; k<=maxX; k++) { |
int l=chYstart*20+12; |
//for(int l=12; l<=16; l++) { |
double signal = h[histogram]->GetBinContent(k,l); |
//signal -= noise; |
//signal /= 5*10000.0; |
double eta = -log(1 - signal); |
double x = k*(xUpUser-xLowUser)/double(binsX); |
//double y = l*(yUpUser-yLowUser)/double(binsY); |
h_line[j]->Fill(x-c_xOffset, signal); |
//} |
} |
if (j == chXstart) { |
h_line[j]->SetTitle("SiPM#2 w/o noise subtraction;x[mm];Hits"); |
//h_line[j]->GetYaxis()->SetRangeUser(-0.05, 0.3); |
//h_line[j]->GetYaxis()->SetRangeUser(-50, 2500); |
h_line[j]->Draw(""); |
} |
else { |
h_line[j]->SetLineColor(j+1); |
h_line[j]->Draw("same"); |
} |
} |
//! 2d scan |
TCanvas* canvas61 = new TCanvas("canvas61","canvas61",8*200,300); |
canvas61->cd(); |
TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
for(int i=chXstart; i<=chXend; i++) { |
//int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
//int canvasPosition = nX*(chYend-i)+chXstart+1; |
//if (debug) printf("canvas %d\n",canvasPosition); |
int chPosition = map[i][chYstart]; |
char hname[128]; |
sprintf(hname, "h2d%d", chPosition); |
int histogram = i; |
h[histogram] = (TH2F *) rootfile->Get(hname); |
int noise = getNoise(h[histogram], 1, 100); |
for(int k=minX; k<=maxX; k++) { |
for(int l=minY; l<=maxY; l++) { |
int signal = h[histogram]->GetBinContent(k,l); // detected |
//p /= 10000.; |
//double p0 = 1.0 - p; // events with zero photons |
//double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
//double eta = -log(p0); |
//printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
//double signal = ((p - noise) > 0.1) ? (p-noise) : 0.1; |
double p = signal - noise; |
p /= 10000.0; |
double eta = -log(1 - p); |
double x = k*(xUpUser-xLowUser)/double(binsX); |
double y = l*(yUpUser-yLowUser)/double(binsY); |
h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
} |
} |
} |
h_corrected->SetTitle("SiPM#2 n_pe = - ln(P0);x[mm];y[mm]"); |
h_corrected->GetZaxis()->SetRangeUser(-0.05,0.30); |
h_corrected->Draw("colz"); |
} |
/** Draws the sum of the channels |
* Each channel is a 2d plot |
* Intended for the study of 1 channel |
*/ |
if (strstr(plopt, "2d") != NULL) { |
int nX = chXend - chXstart + 1; |
int nY = chYend - chYstart + 1; |
TCanvas* canvas7 = new TCanvas("canvas7","canvas7", nX*400,nY*400); |
printf("nx %d ny %d\n",nX,nY); |
canvas7->Divide(nX,nY); |
for(int i=chYstart; i<=chYend; i++) { |
for(int j=chXstart; j<=chXend; j++) { |
//int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
int canvasPosition = nX*(chYend-i)+(j-chXstart)+1; |
if (debug) printf("canvas %d\n",canvasPosition); |
canvas7->cd(canvasPosition); |
char hname[128]; |
int chPosition = map[j][i]; |
sprintf(hname, "h2d%d", chPosition); |
TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
h_2d->Draw("colz"); |
} //x |
} |
// Number of photoelectrons - Poissonian correction |
TCanvas* canvas8 = new TCanvas("canvas8","canvas8", 1000,1000); |
canvas8->cd(); |
gStyle->SetOptStat(0); |
char hname[128]; |
int chPosition = map[1][2]; |
sprintf(hname, "h2d%d", chPosition); |
TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
Int_t binsX = h_2d->GetXaxis()->GetNbins(); |
Int_t minX = h_2d->GetXaxis()->GetFirst(); |
Int_t maxX = h_2d->GetXaxis()->GetLast()+1; |
Int_t binsY = h_2d->GetYaxis()->GetNbins(); |
Int_t minY = h_2d->GetYaxis()->GetFirst(); |
Int_t maxY = h_2d->GetYaxis()->GetLast()+1; |
Double_t xLowUser = h_2d->GetXaxis()->GetBinLowEdge(minX); |
Double_t xUpUser = h_2d->GetXaxis()->GetBinUpEdge(maxX); |
Double_t yLowUser = h_2d->GetYaxis()->GetBinLowEdge(minY); |
Double_t yUpUser = h_2d->GetYaxis()->GetBinUpEdge(maxY); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, xLowUser, xUpUser, binsY, yLowUser, yUpUser); |
double noise = getNoise(h_2d, 1, 89); |
if(debug) printf("Noise = %f\n", noise); |
for(int k=minX; k<=maxX; k++) { |
for(int j=minY; j<=maxY; j++) { |
double signal = h_2d->GetBinContent(k,j); // detected |
//double p = ((signal - noise) > 1) ? (signal-noise) : 1; |
double p = signal - noise; |
p /= 10000.; |
double p0 = 1.0 - p; // events with zero photons |
//double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
double eta = -log(p0); // constant of the poissonian statistics |
if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
double x = k*(xUpUser-xLowUser)/double(binsX); |
double y = j*(yUpUser-yLowUser)/double(binsY); |
h_corrected->Fill(x+3,y+8, eta); |
} |
} |
h_corrected->SetTitle("n_pe = - ln P(0);x[mm];y[mm]"); |
h_corrected->GetZaxis()->SetRangeUser(-0.05,0.3); |
h_corrected->Draw("colz"); |
// collection efficiency |
int nPoints =0; |
double efficiency=0; |
for (int i=18; i<=58; i++) { |
for (int j=19; j<=59; j++) { |
double signal = h_corrected->GetBinContent(i,j); |
if(debug) printf("signal %f\n",signal); |
efficiency += signal; |
nPoints++; |
} |
} |
printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
} |
/** Draws the sum of channel signals |
* Each channel is a 2d ('h2d') histogram |
* Suitable for 8x8 chs scan |
*/ |
if( strstr(plopt, "sum") != NULL ) { |
int nX = chXend - chXstart + 1; |
int nY = chYend - chYstart + 1; |
TCanvas* canvas12 = new TCanvas("canvas12","c2",8*200, 8*200); |
canvas12->cd(); |
gStyle->SetOptStat(0); |
// final histogram parameters |
TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
Int_t binsX = h0->GetXaxis()->GetNbins(); |
Int_t minX = h0->GetXaxis()->GetFirst(); |
Int_t maxX = h0->GetXaxis()->GetLast()+1; |
Int_t binsY = h0->GetYaxis()->GetNbins(); |
Int_t minY = h0->GetYaxis()->GetFirst(); |
Int_t maxY = h0->GetYaxis()->GetLast()+1; |
Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
TH2F* h[9]; |
// 2d histogram noise subtraction and poisson scaling |
for(int i=chYstart; i<=chYend; i++) { |
for(int j=chXstart; j<=chXend; j++) { |
int chPosition = map[j][i]; |
char hname[128]; |
sprintf(hname, "h2d%d", chPosition); |
int histogram = nX*(i-chYstart)+(j-chXstart); |
h[histogram] = (TH2F *) rootfile->Get(hname); |
int noise = getNoise(h[histogram], 1, 170); |
if (debug) printf("noise: %d\n",noise); |
for(int k=minX; k<=maxX; k++) { |
for(int l=minY; l<=maxY; l++) { |
int signal = h[histogram]->GetBinContent(k,l); // detected |
//double p = ((signal - noise) > 0.1) ? (signal-noise) : 0.1; |
double p = signal - noise; |
p /= 10000.; |
double p0 = 1.0 - p; // events with zero photons |
//double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
double eta = -log(p0); |
//printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
double x = k*(xUpUser-xLowUser)/double(binsX); |
double y = l*(yUpUser-yLowUser)/double(binsY); |
h_corrected->Fill(x-c_xOffset, y-c_yOffset, signal); |
} |
} |
} |
} |
h_corrected->SetTitle("SiPM#2 n_p.e.;x[mm];y[mm]"); |
//h_corrected->GetZaxis()->SetRangeUser(-0.05,.30); |
h_corrected->Draw("colz"); |
TCanvas* canvas13 = new TCanvas("canvas13","canvas13",600,300); |
canvas13->Divide(2); |
canvas13->cd(1); |
//h[16]->Draw("colz"); |
canvas13->cd(2); |
//h[8]->Draw("colz"); |
} |
return(0); |
} |
/** Function calculates the noise from one channel |
* it runs through the bins along x and returns the average value |
*/ |
double getNoise(TH2F* histogram, int yStart, int yEnd) |
{ |
double noise=0; |
int count=0; |
for(int j=yStart; j<yEnd; j++) { |
double value = histogram->GetBinContent(j,2); |
//if (noise < value) noise = value; |
noise += value; |
count++; |
} |
return (noise/double(count)); |
} |
/sipmScan/examples/threshold64.C |
---|
0,0 → 1,43 |
{ |
int map[8][8]={{32,34,53,55,40,42,61,63}, |
{48,50,37,39,56,58,45,47}, |
{33,35,52,54,41,43,60,62}, |
{49,51,36,38,57,59,44,46}, |
{17,19,4,6,25,27,12,14}, |
{1,3,20,22,9,11,28,30}, |
{16,18,5,7,24,26,13,15}, |
{0,2,21,23,8,10,29,31} |
}; |
TCanvas* canvas = new TCanvas("canvas","",1200,1200); |
canvas->cd(); |
//canvas->Divide(4,4); |
TVirtualPad *main = new TPad("main","main",0,0,1,1,10,1); |
main->Draw(); |
main->Divide(4,4); |
TH1F* h_thrProj[16]; |
for(int i=0; i<4; i++) { |
for(int j=0; j<4; j++) { |
int channel = map[i][j]; |
printf("%d ", channel); |
TH1D* histogram = h_threshold->ProjectionY("",channel+1,channel+1); |
int nBins = histogram->GetXaxis()->GetNbins(); |
char name[128]; |
int current = i*4+j; |
sprintf(name,"histo%d",current); |
//h_thrProj[i*4+j] = (TH1F*)gROOT->FindObject(name); |
if (h_thrProj[i*4+j]) delete h_thrProj[i*4+j]; |
h_thrProj[current] = new TH1F(name, name, nBins,0,nBins); |
for(int k=0; k< nBins; k++){ |
int signal = histogram->GetBinContent(k); |
h_thrProj[current]->Fill(k,signal); |
} |
int canvasPosition = i+4*(3-j)+1; |
main->cd(canvasPosition); |
//printf("%d \n", canvasPosition); |
h_thrProj[current]->Draw(); |
} |
} |
} |
/sipmScan/examples/sipm.C |
---|
0,0 → 1,538 |
#include "TROOT.h" |
#include "TFile.h" |
#include "TBenchmark.h" |
#include "TH1F.h" |
#include "TH2F.h" |
#include "TCanvas.h" |
#include "TStyle.h" |
#include "TPad.h" |
#include "TF1.h" |
#include "TGraph.h" |
#include "TSpectrum.h" |
#include "stdio.h" |
#include "include/RTUtil.h" |
double getNoise(TH2F*, int, int); |
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) |
{ |
const int c_nChannels = 64; |
const double c_xOffset = 0; // mm |
const double c_yOffset = 0; |
int map[8][8]={{32,34,53,55,40,42,61,63}, |
{48,50,37,39,56,58,45,47}, |
{33,35,52,54,41,43,60,62}, |
{49,51,36,38,57,59,44,46}, |
{17,19,4,6,25,27,12,14}, |
{1,3,20,22,9,11,28,30}, |
{16,18,5,7,24,26,13,15}, |
{0,2,21,23,8,10,29,31} |
}; |
char fnameroot[256]; |
TFile* rootfile; |
sprintf(fnameroot, "root/%s.root", filename); |
rootfile = (TFile *) gROOT->FindObject(filename); |
if(rootfile==NULL) rootfile = new TFile(fnameroot); |
if(rootfile==NULL) { |
printf("Cannot open root file %s!!!\n",fnameroot); |
return(0); |
} |
// set draw style |
gStyle->SetOptStat("ne"); |
gStyle->SetPalette(1, 0); |
gStyle->SetPaperSize(TStyle::kA4); |
gStyle->SetStatBorderSize(1); |
gStyle->SetFrameBorderMode(0); |
gStyle->SetFrameFillColor(0); |
gStyle->SetCanvasBorderMode(0); |
gStyle->SetPadBorderMode(0); |
gStyle->SetPadColor(0); |
gStyle->SetCanvasColor(0); |
gStyle->SetStatColor(0); |
gStyle->SetOptFit(11); |
gStyle->SetOptStat(); |
gStyle->SetPadRightMargin(0.15); |
gStyle->SetPadLeftMargin(0.12); |
//gStyle->SetTitleYOffset(1.4); |
if( strstr(plopt, "all") != NULL ) { |
TCanvas *canvas2 = new TCanvas("canvas2","Hits x;;",2000,2000); |
TCanvas *canvas3 = new TCanvas("canvas3","Hits y;;",2000,2000); |
canvas2->Divide(8,8); |
canvas3->Divide(8,8); |
TH1F* h_hitsx; |
TH1F* h_hitsy; |
for(int i=0; i<c_nChannels; i++) { |
canvas2->cd(i+1); |
char hname[128]; |
sprintf(hname, "hnhitsx%d", i); |
h_hitsx = (TH1F*)rootfile->Get(hname); |
h_hitsx->Draw(); |
canvas3->cd(i+1); |
sprintf(hname, "hnhitsy%d", i); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->Draw(); |
} |
} |
if( strstr(plopt, "x") != NULL ) { |
TCanvas *canvas10 = new TCanvas("canvas10","Ch x;;",500,500); |
TH1F* h_hitsx; |
canvas10->cd(); |
char hname[128]; |
sprintf(hname, "hnhitsx%d", chXstart); |
h_hitsx = (TH1F*)rootfile->Get(hname); |
h_hitsx->Draw(); |
} |
if( strstr(plopt, "y") != NULL ) { |
TCanvas *canvas11 = new TCanvas("canvas11","Ch x;;",500,500); |
TH1F* h_hitsy; |
canvas11->cd(); |
char hname[128]; |
sprintf(hname, "hnhitsy%d", chXstart); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->Draw(); |
} |
if( strstr(plopt, "share") != NULL ) { |
/*TCanvas *canvas4 = new TCanvas("canvas1","canvas1",1000,1000); |
int nChannels = chYend-chYstart+1; |
int ncols = nChannels/2; |
printf("nch %d nch\\2 %d\n", nChannels, ncols); |
canvas4->Divide(2,ncols); |
TH1F* h_hitsy; |
for(int i=chYstart; i<=chYend; i++){ |
canvas4->cd(i-chYstart+1); |
char hname[128]; |
int chPosition = map[0][i]; |
sprintf(hname, "hnhitsy%d", chPosition); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->Draw(); |
}*/ |
TCanvas *canvas4 = new TCanvas("canvas4","canvas4",500,500); |
canvas4->cd(); |
for(int i=chXstart; i<=chXend; i++) { |
TH1F* h_hitsx; |
char hname[128]; |
int chPosition = map[i][chYstart]; |
sprintf(hname, "hnhitsx%d", chPosition); |
h_hitsx = (TH1F*)rootfile->Get(hname); |
h_hitsx->SetTitle("Scan X;x [mm]; Entries"); |
h_hitsx->GetYaxis()->SetTitleOffset(1.3); |
h_hitsx->SetStats(0); |
if (i == chXstart) |
h_hitsx->Draw(); |
else { |
h_hitsx->SetLineColor(i+1); |
h_hitsx->Draw("same"); |
} |
} |
//sprintf(fullname, "ps/%s_Yshare.eps", filename); |
//canvas4->SaveAs(fullname); |
TCanvas *canvas5 = new TCanvas("canvas5","canvas5",500,500); |
canvas5->cd(); |
for(int i=chYstart; i<=chYend; i++) { |
TH1F* h_hitsy; |
char hname[128]; |
int chPosition = map[chXstart][i]; |
sprintf(hname, "hnhitsy%d", chPosition); |
h_hitsy = (TH1F*)rootfile->Get(hname); |
h_hitsy->SetTitle("Scan Y;y [mm]; Entries"); |
h_hitsy->GetYaxis()->SetTitleOffset(1.3); |
h_hitsy->SetStats(0); |
if (i == chYstart) |
h_hitsy->Draw(); |
else { |
h_hitsy->SetLineColor(i+1); |
h_hitsy->Draw("same"); |
} |
} |
//sprintf(fullname, "ps/%s_Yshare.eps", filename); |
//canvas5->SaveAs(fullname); |
} |
/** Draw one channel in x and y |
* Without and with the noise subtraction |
*/ |
if (strstr(plopt, "noise") != NULL) { |
TCanvas *canvas51 = new TCanvas("canvas51","canvas51",500,500); |
canvas51->cd(); |
// Get the filled histogram from scan |
TH1F* h_x; |
char name[32]; |
sprintf(name, "hnhitsx%d", chXstart); |
h_x = (TH1F*)rootfile->Get(name); |
h_x->DrawCopy(); |
TH1F* h_y; |
sprintf(name, "hnhitsy%d", chXstart); |
h_y = (TH1F*)rootfile->Get(name); |
// Create and fill corrected histogram |
Int_t binsX = h_x->GetXaxis()->GetNbins(); |
if(debug) printf("nBins: %d\n", binsX); |
//Double_t xLowUser = h_x->GetXaxis()->GetBinLowEdge(1); |
//Double_t xUpUser = h_x->GetXaxis()->GetBinUpEdge(binsX); |
Double_t xLowUser = h_x->GetXaxis()->GetXmin(); |
Double_t xUpUser = h_x->GetXaxis()->GetXmax(); |
TH1F* h_xCorrected = new TH1F("h_xCorrected",";;",binsX, xLowUser, xUpUser); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
double noise = h_x->GetBinContent(2); |
for (int jk = 0; jk < binsX; jk++) { |
double signal = h_x->GetBinContent(jk+1); |
signal -= noise; |
double x = h_x->GetXaxis()->GetBinCenter(jk+1); |
h_xCorrected->Fill(x, signal); |
} |
h_xCorrected->SetLineColor(2); |
h_xCorrected->DrawCopy("same"); |
//Create and fill corrected histogram |
Int_t binsY = h_y->GetXaxis()->GetNbins(); |
Double_t yLowUser = h_y->GetXaxis()->GetXmin(); |
Double_t yUpUser = h_y->GetXaxis()->GetXmax(); |
TH1F* h_yCorrected = new TH1F("h_yCorrected",";;",binsY, yLowUser, yUpUser); |
noise = h_y->GetBinContent(2); |
for (int jk = 0; jk < binsY; jk++) { |
double signal = h_y->GetBinContent(jk+1); |
signal -= noise; |
double x = h_y->GetXaxis()->GetBinCenter(jk+1); |
h_yCorrected->Fill(x, signal); |
} |
TCanvas *canvas52 = new TCanvas("canvas52","canvas52",500,500); |
canvas52->cd(); |
h_y->DrawCopy(); |
h_yCorrected->SetLineColor(2); |
h_yCorrected->DrawCopy("same"); |
} |
/** Draws the signal from 8 channels in x-row |
* for one specific y bin, so the background and cross-talk |
* can be estimated. |
* Draws also a 2d scan of these channels. |
*/ |
if (strstr(plopt, "line") != NULL) { |
TCanvas* canvas6 = new TCanvas("canvas6","canvas6",500,500); |
canvas6->cd(0); |
gStyle->SetOptStat(0); |
TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
Int_t binsX = h0->GetXaxis()->GetNbins(); |
Int_t minX = h0->GetXaxis()->GetFirst(); |
Int_t maxX = h0->GetXaxis()->GetLast()+1; |
Int_t binsY = h0->GetYaxis()->GetNbins(); |
Int_t minY = h0->GetYaxis()->GetFirst(); |
Int_t maxY = h0->GetYaxis()->GetLast()+1; |
Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
//! 1-dimension position in x vs. hits |
TH2F* h[8]; |
TH1F* h_line[8]; |
for(int j=0; j<8; j++) { |
h_line[j] = new TH1F("h_line", "h_line", binsX, xLowUser, xUpUser); |
} |
for(int j=chXstart; j<=chXend; j++) { |
int chPosition = map[j][chYstart]; |
char hname[128]; |
sprintf(hname, "h2d%d", chPosition); |
int histogram = j; |
h[histogram] = (TH2F *) rootfile->Get(hname); |
int noise = getNoise(h[histogram], 1, 160); |
for(int k=minX; k<=maxX; k++) { |
int l=chYstart*20+12; |
//for(int l=12; l<=16; l++) { |
double signal = h[histogram]->GetBinContent(k,l); |
//signal -= noise; |
//signal /= 5*10000.0; |
double eta = -log(1 - signal); |
double x = xLowUser + k*(xUpUser-xLowUser)/double(binsX); |
//double y = l*(yUpUser-yLowUser)/double(binsY); |
h_line[j]->Fill(x-c_xOffset, signal); |
//} |
} |
if (j == chXstart) { |
h_line[j]->SetTitle("SiPM#2 w/o noise subtraction;x[mm];Hits"); |
//h_line[j]->GetYaxis()->SetRangeUser(-0.05, 0.3); |
//h_line[j]->GetYaxis()->SetRangeUser(-50, 2500); |
h_line[j]->Draw(""); |
} |
else { |
h_line[j]->SetLineColor(j+1); |
h_line[j]->Draw("same"); |
} |
} |
//! 2d scan |
TCanvas* canvas61 = new TCanvas("canvas61","canvas61",8*200,300); |
canvas61->cd(); |
TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
for(int i=chXstart; i<=chXend; i++) { |
//int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
//int canvasPosition = nX*(chYend-i)+chXstart+1; |
//if (debug) printf("canvas %d\n",canvasPosition); |
int chPosition = map[i][chYstart]; |
char hname[128]; |
sprintf(hname, "h2d%d", chPosition); |
int histogram = i; |
h[histogram] = (TH2F *) rootfile->Get(hname); |
int noise = getNoise(h[histogram], 1, 100); |
for(int k=minX; k<=maxX; k++) { |
for(int l=minY; l<=maxY; l++) { |
int signal = h[histogram]->GetBinContent(k,l); // detected |
//p /= 10000.; |
//double p0 = 1.0 - p; // events with zero photons |
//double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
//double eta = -log(p0); |
//printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
//double signal = ((p - noise) > 0.1) ? (p-noise) : 0.1; |
double p = signal - noise; |
p /= 10000.0; |
double eta = -log(1 - p); |
//double x = k*(xUpUser-xLowUser)/double(binsX); |
//double y = l*(yUpUser-yLowUser)/double(binsY); |
double x = h[histogram]->GetXaxis()->GetBinCenter(k); |
double y = h[histogram]->GetYaxis()->GetBinCenter(l); |
h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
} |
} |
} |
h_corrected->SetTitle("SiPM#2 n_pe = - ln(P0);x[mm];y[mm]"); |
h_corrected->GetZaxis()->SetRangeUser(-0.05,0.30); |
h_corrected->Draw("colz"); |
} |
/** Draws the sum of the channels |
* Each channel is a 2d plot |
* Intended for the study of 1 channel |
* |
* Input: scanned channel coordinates (x,y) |
*/ |
if (strstr(plopt, "2d") != NULL) { |
int nX = 3; |
int nY = 3; |
TCanvas* canvas7 = new TCanvas("canvas7","canvas7", nX*400,nY*400); |
printf("nx %d ny %d\n",nX,nY); |
canvas7->Divide(nX,nY); |
for(int i=chXend-1; i<=chXend+1; i++) { |
for(int j=chXstart-1; j<=chXstart+1; j++) { |
//int canvasPosition = nX*(i-chYstart)+(j-chXstart)+1; |
int canvasPosition = nX*(chXend+1-i) + (j-chXstart+1) +1; |
if (debug) printf("canvas %d\n",canvasPosition); |
canvas7->cd(canvasPosition); |
char hname[128]; |
int chPosition = map[j][i]; |
sprintf(hname, "h2d%d", chPosition); |
TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
h_2d->Draw("colz"); |
} //x |
} |
// Number of photoelectrons - Poissonian correction |
TCanvas* canvas8 = new TCanvas("canvas8","canvas8", 1000,1000); |
canvas8->cd(); |
gStyle->SetOptStat(0); |
char hname[128]; |
int chPosition = map[chXstart][chXend]; |
sprintf(hname, "h2d%d", chPosition); |
TH2F* h_2d = (TH2F*)rootfile->Get(hname); |
Int_t binsX = h_2d->GetXaxis()->GetNbins(); |
Int_t minX = h_2d->GetXaxis()->GetFirst(); |
Int_t maxX = h_2d->GetXaxis()->GetLast()+1; |
Int_t binsY = h_2d->GetYaxis()->GetNbins(); |
Int_t minY = h_2d->GetYaxis()->GetFirst(); |
Int_t maxY = h_2d->GetYaxis()->GetLast()+1; |
Double_t xLowUser = h_2d->GetXaxis()->GetXmin(); |
Double_t xUpUser = h_2d->GetXaxis()->GetXmax(); |
Double_t yLowUser = h_2d->GetYaxis()->GetXmin(); |
Double_t yUpUser = h_2d->GetYaxis()->GetXmax(); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, xLowUser, xUpUser, binsY, yLowUser, yUpUser); |
//TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX, minX, maxX, binsY, minY, maxY); |
double noise = getNoise(h_2d, 1, 70); |
if(debug) printf("Noise = %f\n", noise); |
for(int k=minX; k<=maxX; k++) { |
for(int j=minY; j<=maxY; j++) { |
double signal = h_2d->GetBinContent(k,j); // detected |
//double p = ((signal - noise) > 1) ? (signal-noise) : 1; |
double p = signal - noise; |
p /= 10000.; |
double p0 = 1.0 - p; // events with zero photons |
//double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
double eta = -log(p0); // constant of the poissonian statistics |
if (debug) printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
//double x = xLowUser + k*(xUpUser - xLowUser) / double(binsX); |
double x = h_2d->GetXaxis()->GetBinCenter(k); |
//double y = yLowUser + j*(yUpUser-yLowUser)/double(binsY); |
double y = h_2d->GetYaxis()->GetBinCenter(j); |
if (debug) printf("x=%f y=%f\n", x, y); |
h_corrected->Fill(x + c_xOffset, y + c_yOffset, eta); |
} |
} |
//h_corrected->SetTitle("n_pe = - ln P(0);x[mm];y[mm]"); |
h_corrected->SetTitle(";x[mm];y[mm]"); |
h_corrected->GetZaxis()->SetRangeUser(-0.01,0.5); |
h_corrected->SetContour(50); |
h_corrected->Draw("colz"); |
// collection efficiency |
int nPoints =0; |
double efficiency=0; |
for (int i=18; i<=58; i++) { |
for (int j=19; j<=59; j++) { |
double signal = h_corrected->GetBinContent(i,j); |
if(debug) printf("signal %f\n",signal); |
efficiency += signal; |
nPoints++; |
} |
} |
printf("Signal sum = %f\n # of points = %d\n",efficiency,nPoints); |
} |
/** Draws the sum of channel signals |
* Each channel is a 2d ('h2d') histogram |
* Suitable for 8x8 chs scan |
*/ |
if( strstr(plopt, "sum") != NULL ) { |
int nX = chXend - chXstart + 1; |
int nY = chYend - chYstart + 1; |
TCanvas* canvas12 = new TCanvas("canvas12","c2",8*200, 8*200); |
canvas12->cd(); |
gStyle->SetOptStat(0); |
// final histogram parameters |
TH2F* h0 = (TH2F*) rootfile->Get("h2d0"); |
Int_t binsX = h0->GetXaxis()->GetNbins(); |
Int_t minX = h0->GetXaxis()->GetFirst(); |
Int_t maxX = h0->GetXaxis()->GetLast()+1; |
Int_t binsY = h0->GetYaxis()->GetNbins(); |
Int_t minY = h0->GetYaxis()->GetFirst(); |
Int_t maxY = h0->GetYaxis()->GetLast()+1; |
Double_t xLowUser = h0->GetXaxis()->GetBinLowEdge(minX); |
Double_t xUpUser = h0->GetXaxis()->GetBinUpEdge(maxX); |
Double_t yLowUser = h0->GetYaxis()->GetBinLowEdge(minY); |
Double_t yUpUser = h0->GetYaxis()->GetBinUpEdge(maxY); |
if (debug) printf("xLow %f xUp %f\n",xLowUser,xUpUser); |
TH2F* h_corrected = new TH2F("h_corrected","h_corrected",binsX,xLowUser,xUpUser, binsY,yLowUser,yUpUser); |
TH2F* h[nX*nY]; |
// 2d histogram noise subtraction and poisson scaling |
for(int i=chYstart; i<=chYend; i++) { |
for(int j=chXstart; j<=chXend; j++) { |
int chPosition = map[j][i]; |
char hname[128]; |
sprintf(hname, "h2d%d", chPosition); |
int histogram = nX*(i-chYstart)+(j-chXstart); |
h[histogram] = (TH2F *) rootfile->Get(hname); |
int noise = getNoise(h[histogram], 1, 170); |
if (debug) printf("noise: %d\n",noise); |
for(int k=minX; k<=maxX; k++) { |
for(int l=minY; l<=maxY; l++) { |
int signal = h[histogram]->GetBinContent(k,l); // detected |
//double p = ((signal - noise) > 0.1) ? (signal-noise) : 0.1; |
double p = signal - noise; |
p /= 10000.; |
double p0 = 1.0 - p; // events with zero photons |
//double eta = (-log(p0) * p0 ) / (1-p0-0.00001); |
double eta = -log(p0); |
//printf("p=%f p0=%f log(p0)=%f eta=%f\n",p,p0,log(p0),eta); |
double x = k*(xUpUser-xLowUser)/double(binsX); |
double y = l*(yUpUser-yLowUser)/double(binsY); |
h_corrected->Fill(x-c_xOffset, y-c_yOffset, eta); |
} |
} |
} |
} |
h_corrected->SetTitle("SiPM#2 n_p.e.;x[mm];y[mm]"); |
h_corrected->GetZaxis()->SetRangeUser(-0.05,.30); |
h_corrected->Draw("colz"); |
TCanvas* canvas13 = new TCanvas("canvas13","canvas13",600,300); |
canvas13->Divide(2); |
canvas13->cd(1); |
h[16]->Draw("colz"); |
canvas13->cd(2); |
h[8]->Draw("colz"); |
} |
/** Draws the beam profile and fits it with error function |
* on some background function |
*/ |
if (strstr(plopt, "beam") != NULL) { |
TCanvas* canvas9 = new TCanvas("canvas9","canvas9", 500,500); |
canvas9->cd(); |
char hname[128]; |
sprintf(hname, "hnhitsx%d", 36); |
TH1F* h_laser = (TH1F*)rootfile->Get(hname); |
h_laser->Draw(); |
h_laser->SetStats(1); |
TF1* err = new TF1("err","[0]+[1]*TMath::Erf((x-[2])/[3])",17.02,17.30); |
err->SetParameter(0,2500); |
err->SetParameter(1, h_laser->GetMaximum()); |
err->SetParameter(2, h_laser->GetBinCenter(h_laser->GetMaximumBin())); |
err->SetParameter(3, 0.001); |
h_laser->Fit(err,"qr"); |
h_laser->Fit(err,"lr"); |
double sigma = err->GetParameter(3); |
printf("sigma = %2.0f um, FWHM = %2.0f um\n", sigma*1000, 2.35*sigma*1000); |
} |
if (strstr(plopt, "map") != NULL) { |
TCanvas* canvas1 = new TCanvas("canvas1","canvas1",1000,1000); |
canvas1->cd(); |
TH2I* h_map = new TH2I("h_map","h_map",8,-0.5,7.5,8,-0.5,7.5); |
for (int i=7; i>=0; i--) { |
//for (int j=7; j>=0; j--) printf("(%d, %d) ", j,i); |
for (int j=7; j>=0; j--) { |
printf("%2d (%d %d)\n", map[j][i], j*10080+5040, i*10080+5040); |
h_map->Fill(j,i, map[j][i]); |
} |
printf("\n"); |
} |
gStyle->SetOptStat(0); |
h_map->Draw("text"); |
} |
return(0); |
} |
/** Function calculates the noise from one channel |
* it runs through the bins along x and returns the average value |
*/ |
double getNoise(TH2F* histogram, int yStart, int yEnd) |
{ |
double noise=0; |
int count=0; |
for(int j=yStart; j<yEnd; j++) { |
double value = histogram->GetBinContent(j,2); |
//if (noise < value) noise = value; |
noise += value; |
count++; |
} |
return (noise/double(count)); |
} |
/sipmScan/examples/tdc.C |
---|
0,0 → 1,109 |
#include "TROOT.h" |
#include "TFile.h" |
#include "TBenchmark.h" |
#include "TH1F.h" |
#include "TH2F.h" |
#include "TCanvas.h" |
#include "TStyle.h" |
#include "TPad.h" |
#include "TF1.h" |
#include "TGraph.h" |
#include "TSpectrum.h" |
#include "stdio.h" |
#include "include/RTUtil.h" |
int tdc(char filename[256] = "test", int chX=0, int chY=0, double rangeLeft=-16, double rangeRight=16, bool debug = false) |
{ |
//const int c_nChannels = 64; |
//const double c_xOffset = 2.2; // mm |
//const double c_yOffset = 2.3; |
int map[8][8]={{32,34,53,55,40,42,61,63}, |
{48,50,37,39,56,58,45,47}, |
{33,35,52,54,41,43,60,62}, |
{49,51,36,38,57,59,44,46}, |
{17,19,4,6,25,27,12,14}, |
{1,3,20,22,9,11,28,30}, |
{16,18,5,7,24,26,13,15}, |
{0,2,21,23,8,10,29,31} |
}; |
char fnameroot[256]; |
TFile* rootfile; |
sprintf(fnameroot, "root/%s.root", filename); |
rootfile = (TFile *) gROOT->FindObject(filename); |
if(rootfile==NULL) rootfile = new TFile(fnameroot); |
if(rootfile==NULL) { |
printf("Cannot open root file %s!!!\n",fnameroot); |
return(0); |
} |
// set draw style |
gStyle->SetPalette(1, 0); |
gStyle->SetPaperSize(TStyle::kA4); |
gStyle->SetStatBorderSize(1); |
gStyle->SetFrameBorderMode(0); |
gStyle->SetFrameFillColor(0); |
gStyle->SetCanvasBorderMode(0); |
gStyle->SetPadBorderMode(0); |
gStyle->SetPadColor(0); |
gStyle->SetCanvasColor(0); |
gStyle->SetStatColor(0); |
gStyle->SetOptFit(11); |
gStyle->SetOptStat(0); |
gStyle->SetPadRightMargin(0.15); |
gStyle->SetPadLeftMargin(0.12); |
//gStyle->SetTitleYOffset(1.4); |
TCanvas* canvas1 = new TCanvas("canvas1","canvas1",1000,1000); |
TH2F* htdc = (TH2F*) rootfile->Get("htdc"); |
canvas1->cd(); |
htdc->Draw("colz"); |
TH2F* h_correctedTDC = (TH2F*) rootfile->Get("h_correctedTDC"); |
TCanvas* canvas2 = new TCanvas("canvas2","canvas2",800,800); |
canvas2->cd(); |
h_correctedTDC->Draw("colz"); |
/* |
canvas2->cd(2); |
int binY = map[chX][chY]; |
TH1D* channelY = h_correctedTDC->ProjectionX("",binY+1,binY+1); |
//channelY->SetStats(0); |
char title[256]; |
sprintf(title,"Channel %d;t [ns];Events", binY); |
channelY->SetTitle(title); |
channelY->GetYaxis()->SetTitleOffset(1.7); |
//TAxis* xAxis = h_correctedTDC->GetXaxis(); |
//int range = xAxis->GetBinUpEdge(xAxis->GetLast()+1); |
//channelY->GetXaxis()->SetRangeUser(-range, range); |
channelY->Draw(); |
TF1* f_gaus1 = new TF1("f_gaus1","[0] + gaus(1)", rangeLeft,rangeRight); |
TF1* f_gaus2 = new TF1("f_gaus2","[0] + gaus(1) + gaus(4)",-8,8); |
f_gaus1->SetParNames("Linear","Norm","#mu","#sigma"); |
f_gaus2->SetParNames("Linear","Norm1","Mean1","Sigma1","Norm2","Mean2","Sigma2"); |
Int_t n = channelY->GetMaximum(); |
Float_t mean = channelY->GetBinCenter(channelY->GetMaximumBin()); |
f_gaus1->SetParameters(channelY->GetMinimum(), n, mean, 2.0); |
channelY->Fit(f_gaus1,"0"); |
channelY->Fit(f_gaus1,"r"); |
*/ |
//Fit in range +- 3 sigma |
//channelY->Fit(f_gaus1,"r","", f_gaus1->GetParameter(2)-3*f_gaus1->GetParameter(3), |
// f_gaus1->GetParameter(2)+3*f_gaus1->GetParameter(3)); |
//f_gaus2->SetParameters(-1, n, mean-0.1, 0.01, |
// n, mean+0.1, 0.1); |
//channelY->Fit(f_gaus2,"0"); |
//channelY->Fit(f_gaus2,"rl"); |
// set range to +- 3 sigma |
//channelY->Fit(f_gaus,"lr","",f_gaus->GetParameter(2)-3*f_gaus->GetParameter(3), |
// f_gaus->GetParameter(2)+3*f_gaus->GetParameter(3)); |
return (0); |
} |