59,6 → 59,7 |
else if(analTab->GetCurrent() == 1) // Relative PDE |
{ |
relPde->widgetChBox[0]->SetState(kButtonDown); |
relPde->widgetChBox[1]->SetState(kButtonUp); |
midPeak->widgetChBox[0]->SetState(kButtonUp); |
zeroAngle->widgetNE[0]->SetNumber(0.00); |
} |
96,6 → 97,9 |
else if( analtab == 3 ) |
analtype = 4; |
|
// Save analysis settings any time we run a new analysis |
SaveAnalSettings(); |
|
// Only integrate spectrum or make relative PDE |
if(type == 0) |
{ |
370,7 → 374,6 |
double range[4]; |
TGraph2D *gScan2D; |
gScan2D = new TGraph2D(); |
gScan2D->SetName("edgescan"); |
range[0] = TMath::MinElement(nrfiles, surfxy); |
range[1] = TMath::MaxElement(nrfiles, surfxy); |
range[2] = TMath::MinElement(nrfiles, surfz); |
679,6 → 682,7 |
int adcpedestal[2]; |
int zeromu = 0; |
int darkhist = -1; |
int nopeaks = -1; |
|
double pointest[12]; |
bool exclude = false; |
690,9 → 694,170 |
analysisProgress->widgetPB->SetPosition(progVal); |
gVirtualX->Update(1); |
|
// Check if the checkbox for no peaks is selected - TODO: Still need the situation when we do not have a peaked ADC spectrum |
if(relPde->widgetChBox[1]->IsDown()) |
{ |
printf("PhotonMu(): ADC spectrum has no peak structure.\n"); |
nopeaks = 1; |
|
// Error if there is no darkhist |
if(strcmp("", darkRun->widgetTE->GetText()) == 0) |
{ |
printf("PhotonMu(): Error! The no peak structure option needs a dark histogram.\n"); |
delete[] integralCount; |
delete[] integralPedestal; |
delete[] angle; |
delete[] pdeval; |
delete[] muval; |
return; |
} |
|
} |
|
// Start if we select at least one file |
if(nrfiles > 0) |
{ |
// Find the pedestal peak for the dark histogram, and use it for all if there are no peaks |
if(nopeaks != -1) |
{ |
// Replot the spectrum on analysisCanvas and do not close the input file |
DisplayHistogram( (char*)(darkRun->widgetTE->GetText()), 0, 1); |
analysisCanvas->GetCanvas()->Modified(); |
analysisCanvas->GetCanvas()->Update(); |
|
// Get the spectrum |
histtemp = (TH1F*)analysisCanvas->GetCanvas()->GetPrimitive(histname); |
npeaks = 15; |
double par[300]; |
spec = new TSpectrum(npeaks); |
// Find spectrum background |
histback = spec->Background(histtemp, (int)fitInter->widgetNE[0]->GetNumber(), "same"); |
// Clone histogram and subtract background from it if we select that option |
h2 = (TH1F*)histtemp->Clone("h2"); |
if(fitChecks->widgetChBox[0]->IsDown()) |
h2->Add(histback, -1); |
// Search for the peaks |
int found = spec->Search(h2, fitSigma->widgetNE[0]->GetNumber(), "goff", fitTresh->widgetNE[0]->GetNumber() ); |
printf("PhotonMu(): Found %d candidates to fit.\n",found); |
npeaks = found; |
|
// Set initial peak parameters |
xpeaks = spec->GetPositionX(); |
for(j = 0; j < found; j++) |
{ |
float xp = xpeaks[j]; |
int bin = h2->GetXaxis()->FindBin(xp); |
float yp = h2->GetBinContent(bin); |
par[3*j] = yp; |
par[3*j+1] = xp; |
par[3*j+2] = (double)fitSigma->widgetNE[0]->GetNumber(); |
} |
|
// Fit the histogram |
fit = new TF1("fit", FindPeaks, adcRange->widgetNE[0]->GetNumber(), adcRange->widgetNE[1]->GetNumber(), 3*npeaks); |
TVirtualFitter::Fitter(histtemp, 3*npeaks); |
fit->SetParameters(par); |
fit->SetNpx(300); |
h2->Fit("fit","Q"); |
// Get the fitted parameters |
fittingfunc = h2->GetFunction("fit"); |
fparam = fittingfunc->GetParameters(); |
fparamerr = fittingfunc->GetParErrors(); |
|
// Gather the parameters (mean peak value for now) |
int j = 1; |
int nrfit = 0; |
while(1) |
{ |
if( (fparam[j] < 1.E-30) || (nrfit > 8) ) |
break; |
else |
{ |
// Check if pedestal is above the lower limit and sigma is smaller than the mean |
if( (fparam[j] > pedesLow->widgetNE[0]->GetNumber()) && ((double)fparamerr[j]/fparam[j] < accError->widgetNE[0]->GetNumber()) ) |
{ |
// With the additional ADC offset, we can shift the mean values slightly, so they are not close to the X.5, but to the X.0 values |
meansel[nrfit] = fparam[j]+(adcOffset->widgetNE[0]->GetNumber()); |
sigmasel[nrfit] = fparam[j+1]; |
nrfit++; |
} |
} |
|
j+=3; |
} |
TMath::Sort(nrfit, meansel, sortindex, kFALSE); |
|
meanparam = meansel[sortindex[0]]; |
paramsigma = sigmasel[sortindex[0]]; |
|
for(j = 0; j < nrfit; j++) |
if(DBGSIG) |
printf("PhotonMu(): %d: peak mean = %lf\n", j, meansel[sortindex[j]]); |
|
j = 0; |
adcpedestal[0] = 0; |
adcpedestal[1] = -1; |
|
while(1) |
{ |
int bin = histtemp->GetXaxis()->FindBin((int)(j+meanparam+paramsigma)); |
int yp = histtemp->GetBinContent(bin); |
|
// Check where we get to first minimum after pedestal peak or where we get to the half maximum of the pedestal peak (in case there is only a pedestal peak) |
if(adcpedestal[1] == -1) |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else |
{ |
if( (npeaks > 1) && (adcpedestal[1] >= yp) ) |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else if( (npeaks == 1) && (adcpedestal[0] < meanparam+5*paramsigma) ) // TODO -> Determining the pedestal when only one peak |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else |
break; |
} |
|
j++; |
if(j > 50) break; |
} |
|
if(midPeak->widgetChBox[0]->IsDown()) |
{ |
if( (meanparam - (int)meanparam >= 0.) && (meanparam - (int)meanparam < 0.5) ) |
m = TMath::Floor(meanparam); |
else if( (meanparam - (int)meanparam >= 0.5) && (meanparam - (int)meanparam < 1.) ) |
m = TMath::Ceil(meanparam); |
int bin = histtemp->GetXaxis()->FindBin(m); |
adcpedestal[0] = m; |
printf("midpeak x = %d, ", adcpedestal[0]); |
adcpedestal[1] = histtemp->GetBinContent(bin); |
} |
|
// Option to show the fit |
fittingfunc->Draw("L SAME"); |
analysisCanvas->GetCanvas()->Modified(); |
analysisCanvas->GetCanvas()->Update(); |
|
printf("Pedestal ends = %d and nr. of counts = %d\n", adcpedestal[0], adcpedestal[1]); |
|
// Delete the opened histogram and spectrum |
delete spec; |
delete inroot; |
|
// return; |
} |
|
printf("PhotonMu(): Continuing with the rest of the spectra.\n"); |
|
// Check all histograms for pedestal peak values |
for(int i = 0; i < (int)nrfiles; i++) |
{ |
if( (nrfiles == 1) || (!multiSelect->widgetChBox[0]->IsDown()) ) |
788,145 → 953,148 |
darkhist = i; |
} |
|
// Replot the spectrum on analysisCanvas and do not close the input file |
DisplayHistogram( (char*)(files->At(i)->GetTitle()), 0, 1); |
analysisCanvas->GetCanvas()->Modified(); |
analysisCanvas->GetCanvas()->Update(); |
if(nopeaks == -1) |
{ |
// Replot the spectrum on analysisCanvas and do not close the input file |
DisplayHistogram( (char*)(files->At(i)->GetTitle()), 0, 1); |
analysisCanvas->GetCanvas()->Modified(); |
analysisCanvas->GetCanvas()->Update(); |
|
// Get the spectrum |
histtemp = (TH1F*)analysisCanvas->GetCanvas()->GetPrimitive(histname); |
npeaks = 15; |
double par[300]; |
spec = new TSpectrum(npeaks); |
// Find spectrum background |
histback = spec->Background(histtemp, (int)fitInter->widgetNE[0]->GetNumber(), "same"); |
// Clone histogram and subtract background from it if we select that option |
h2 = (TH1F*)histtemp->Clone("h2"); |
if(fitChecks->widgetChBox[0]->IsDown()) |
h2->Add(histback, -1); |
// Search for the peaks |
int found = spec->Search(h2, fitSigma->widgetNE[0]->GetNumber(), "goff", fitTresh->widgetNE[0]->GetNumber() ); |
printf("PhotonMu(): Found %d candidates to fit.\n",found); |
npeaks = found; |
// Get the spectrum |
histtemp = (TH1F*)analysisCanvas->GetCanvas()->GetPrimitive(histname); |
npeaks = 15; |
double par[300]; |
spec = new TSpectrum(npeaks); |
// Find spectrum background |
histback = spec->Background(histtemp, (int)fitInter->widgetNE[0]->GetNumber(), "same"); |
// Clone histogram and subtract background from it if we select that option |
h2 = (TH1F*)histtemp->Clone("h2"); |
if(fitChecks->widgetChBox[0]->IsDown()) |
h2->Add(histback, -1); |
// Search for the peaks |
int found = spec->Search(h2, fitSigma->widgetNE[0]->GetNumber(), "goff", fitTresh->widgetNE[0]->GetNumber() ); |
printf("PhotonMu(): Found %d candidates to fit.\n",found); |
npeaks = found; |
|
// Set initial peak parameters |
xpeaks = spec->GetPositionX(); |
for(j = 0; j < found; j++) |
{ |
float xp = xpeaks[j]; |
int bin = h2->GetXaxis()->FindBin(xp); |
float yp = h2->GetBinContent(bin); |
par[3*j] = yp; |
par[3*j+1] = xp; |
par[3*j+2] = (double)fitSigma->widgetNE[0]->GetNumber(); |
} |
// Set initial peak parameters |
xpeaks = spec->GetPositionX(); |
for(j = 0; j < found; j++) |
{ |
float xp = xpeaks[j]; |
int bin = h2->GetXaxis()->FindBin(xp); |
float yp = h2->GetBinContent(bin); |
par[3*j] = yp; |
par[3*j+1] = xp; |
par[3*j+2] = (double)fitSigma->widgetNE[0]->GetNumber(); |
} |
|
// Fit the histogram |
fit = new TF1("fit", FindPeaks, adcRange->widgetNE[0]->GetNumber(), adcRange->widgetNE[1]->GetNumber(), 3*npeaks); |
TVirtualFitter::Fitter(histtemp, 3*npeaks); |
fit->SetParameters(par); |
fit->SetNpx(300); |
h2->Fit("fit","Q"); |
// Get the fitted parameters |
fittingfunc = h2->GetFunction("fit"); |
fparam = fittingfunc->GetParameters(); |
fparamerr = fittingfunc->GetParErrors(); |
// Fit the histogram |
fit = new TF1("fit", FindPeaks, adcRange->widgetNE[0]->GetNumber(), adcRange->widgetNE[1]->GetNumber(), 3*npeaks); |
TVirtualFitter::Fitter(histtemp, 3*npeaks); |
fit->SetParameters(par); |
fit->SetNpx(300); |
h2->Fit("fit","Q"); |
// Get the fitted parameters |
fittingfunc = h2->GetFunction("fit"); |
fparam = fittingfunc->GetParameters(); |
fparamerr = fittingfunc->GetParErrors(); |
|
// Gather the parameters (mean peak value for now) |
int j = 1; |
int nrfit = 0; |
while(1) |
{ |
if( (fparam[j] < 1.E-30) || (nrfit > 8) ) |
break; |
else |
{ |
// Check if pedestal is above the lower limit and sigma is smaller than the mean |
if( (fparam[j] > pedesLow->widgetNE[0]->GetNumber()) && ((double)fparamerr[j]/fparam[j] < accError->widgetNE[0]->GetNumber()) ) |
// Gather the parameters (mean peak value for now) |
int j = 1; |
int nrfit = 0; |
while(1) |
{ |
if( (fparam[j] < 1.E-30) || (nrfit > 8) ) |
break; |
else |
{ |
// With the additional ADC offset, we can shift the mean values slightly, so they are not close to the X.5, but to the X.0 values |
meansel[nrfit] = fparam[j]+(adcOffset->widgetNE[0]->GetNumber()); |
sigmasel[nrfit] = fparam[j+1]; |
nrfit++; |
// Check if pedestal is above the lower limit and sigma is smaller than the mean |
if( (fparam[j] > pedesLow->widgetNE[0]->GetNumber()) && ((double)fparamerr[j]/fparam[j] < accError->widgetNE[0]->GetNumber()) ) |
{ |
// With the additional ADC offset, we can shift the mean values slightly, so they are not close to the X.5, but to the X.0 values |
meansel[nrfit] = fparam[j]+(adcOffset->widgetNE[0]->GetNumber()); |
sigmasel[nrfit] = fparam[j+1]; |
nrfit++; |
} |
} |
} |
|
j+=3; |
} |
TMath::Sort(nrfit, meansel, sortindex, kFALSE); |
j+=3; |
} |
TMath::Sort(nrfit, meansel, sortindex, kFALSE); |
|
meanparam = meansel[sortindex[0]]; |
paramsigma = sigmasel[sortindex[0]]; |
meanparam = meansel[sortindex[0]]; |
paramsigma = sigmasel[sortindex[0]]; |
|
for(j = 0; j < nrfit; j++) |
if(DBGSIG) |
printf("PhotonMu(): %d: peak mean = %lf\n", j, meansel[sortindex[j]]); |
for(j = 0; j < nrfit; j++) |
if(DBGSIG) |
printf("PhotonMu(): %d: peak mean = %lf\n", j, meansel[sortindex[j]]); |
|
j = 0; |
adcpedestal[0] = 0; |
adcpedestal[1] = -1; |
j = 0; |
adcpedestal[0] = 0; |
adcpedestal[1] = -1; |
|
while(1) |
{ |
int bin = histtemp->GetXaxis()->FindBin((int)(j+meanparam+paramsigma)); |
int yp = histtemp->GetBinContent(bin); |
while(1) |
{ |
int bin = histtemp->GetXaxis()->FindBin((int)(j+meanparam+paramsigma)); |
int yp = histtemp->GetBinContent(bin); |
|
// Check where we get to first minimum after pedestal peak or where we get to the half maximum of the pedestal peak (in case there is only a pedestal peak) |
if(adcpedestal[1] == -1) |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else |
{ |
if( (npeaks > 1) && (adcpedestal[1] >= yp) ) |
// Check where we get to first minimum after pedestal peak or where we get to the half maximum of the pedestal peak (in case there is only a pedestal peak) |
if(adcpedestal[1] == -1) |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else if( (npeaks == 1) && (adcpedestal[0] < meanparam+5*paramsigma) ) // TODO -> Determining the pedestal when only one peak |
else |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
if( (npeaks > 1) && (adcpedestal[1] >= yp) ) |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else if( (npeaks == 1) && (adcpedestal[0] < meanparam+5*paramsigma) ) // TODO -> Determining the pedestal when only one peak |
{ |
adcpedestal[0] = j+meanparam+paramsigma; |
adcpedestal[1] = yp; |
} |
else |
break; |
} |
else |
break; |
|
j++; |
if(j > 50) break; |
} |
|
j++; |
if(j > 50) break; |
} |
|
if( (npeaks > 1) && (nrfit > 1) ) |
{ |
int bin = histtemp->GetXaxis()->FindBin((int)(meanparam+meansel[sortindex[1]])/2); |
adcpedestal[0] = (meanparam+meansel[sortindex[1]])/2; |
printf("PhotonMu(): multipeak x = %d, ", adcpedestal[0]); |
adcpedestal[1] = histtemp->GetBinContent(bin); |
} |
if(npeaks > 1) |
{ |
int bin = histtemp->GetXaxis()->FindBin((int)(meanparam+meansel[sortindex[1]])/2); |
adcpedestal[0] = (meanparam+meansel[sortindex[1]])/2; |
printf("PhotonMu(): multipeak x = %d, ", adcpedestal[0]); |
adcpedestal[1] = histtemp->GetBinContent(bin); |
} |
|
if(midPeak->widgetChBox[0]->IsDown()) |
{ |
if( (meanparam - (int)meanparam >= 0.) && (meanparam - (int)meanparam < 0.5) ) |
m = TMath::Floor(meanparam); |
else if( (meanparam - (int)meanparam >= 0.5) && (meanparam - (int)meanparam < 1.) ) |
m = TMath::Ceil(meanparam); |
int bin = histtemp->GetXaxis()->FindBin(m); |
adcpedestal[0] = m; |
printf("midpeak x = %d, ", adcpedestal[0]); |
adcpedestal[1] = histtemp->GetBinContent(bin); |
} |
if(midPeak->widgetChBox[0]->IsDown()) |
{ |
if( (meanparam - (int)meanparam >= 0.) && (meanparam - (int)meanparam < 0.5) ) |
m = TMath::Floor(meanparam); |
else if( (meanparam - (int)meanparam >= 0.5) && (meanparam - (int)meanparam < 1.) ) |
m = TMath::Ceil(meanparam); |
int bin = histtemp->GetXaxis()->FindBin(m); |
adcpedestal[0] = m; |
printf("midpeak x = %d, ", adcpedestal[0]); |
adcpedestal[1] = histtemp->GetBinContent(bin); |
} |
|
/* // Option to show the fit |
fittingfunc->Draw("L SAME"); |
analysisCanvas->GetCanvas()->Modified(); |
analysisCanvas->GetCanvas()->Update();*/ |
/* // Option to show the fit |
fittingfunc->Draw("L SAME"); |
analysisCanvas->GetCanvas()->Modified(); |
analysisCanvas->GetCanvas()->Update();*/ |
|
printf("Pedestal ends = %d and nr. of counts = %d\n", adcpedestal[0], adcpedestal[1]); |
printf("Pedestal ends = %d and nr. of counts = %d\n", adcpedestal[0], adcpedestal[1]); |
|
// Delete the opened histogram and spectrum |
delete spec; |
delete inroot; |
// Delete the opened histogram and spectrum |
delete spec; |
delete inroot; |
} |
|
// Open the input file and read header, ADC and TDC values |
sprintf(ctemp, "%s", files->At(i)->GetTitle()); |
994,14 → 1162,7 |
if( (angle[i] == zeroAngle->widgetNE[0]->GetNumber()) && (darkhist != i) ) |
zeromu = i; |
|
// Checking for errors when fitting a histogram |
if(k2 == 0) |
{ |
printf("PhotonMu(): No pedestal entries found. Check the fitting results.\n"); |
muval[i] = -1; |
} |
else |
muval[i] = -TMath::Log((double)k2/(double)k); |
muval[i] = -TMath::Log((double)k2/(double)k); |
printf("PhotonMu(): %lf: muval = %lf\n", angle[i], muval[i]); |
|
inroot->Close(); |
1018,7 → 1179,7 |
|
printf("PhotonMu(): angle\tmu\trelative PDE\n"); |
m = 0; |
|
|
// Set the 0 degree muval, reuse meansel[1] |
meansel[1] = muval[zeromu]; |
printf("Zero value (id=%d, angle=%lf) = %lf\n", zeromu, angle[zeromu], meansel[1]); |
1048,10 → 1209,6 |
if(i == darkhist) |
exclude = false; |
|
// Wrong fit |
if(muval[i] == -1) |
exclude = true; |
|
// If nothing excluded, pass the points in pointest variable like in a FIFO |
if(!exclude) |
{ |
1312,6 → 1469,7 |
{ |
// With the additional ADC offset, we can shift the mean values slightly, so they are not close to the X.5, but to the X.0 values |
meansel[nrfit] = fparam[j]+(adcOffset->widgetNE[0]->GetNumber()); |
meanselerr[nrfit] = fparamerr[j]; |
sigmasel[nrfit] = fparam[j+1]; |
nrfit++; |
} |
1488,6 → 1646,8 |
exclude = false; |
} |
|
printf("BreakdownVolt(): Calculated separation between peaks %d and %d: sep = %lf, seperr = %lf\n", (int)peakSepCalc->widgetNE[0]->GetNumber()-1, (int)peakSepCalc->widgetNE[0]->GetNumber(), sep[(int)peakSepCalc->widgetNE[0]->GetNumber()-1][p], seperr[(int)peakSepCalc->widgetNE[0]->GetNumber()-1][p]); |
|
// Write out parameters to a file |
fp = fopen(paramname, "a"); |
if(exclude) |
1540,7 → 1700,7 |
bdplot = new TGraphErrors(p, volt, sep[k-1], volterr, seperr[k-1]); |
else |
{ |
printf("BreakdownVold(): Unsupported peak separation selected (%d).\n", k); |
printf("BreakdownVolt(): Unsupported peak separation selected (%d).\n", k); |
return; |
} |
|
1577,7 → 1737,7 |
sprintf(ctemp, "#splitline{#Delta_{p}(U) = (%.2lf #pm %.2lf)#timesU + (%.2lf #pm %.3lf)}{U_{0} = %.2lf #pm %.3lf}", meansel[0], meanselerr[0], meansel[1], meanselerr[1], meansel[2], meansel[2]*(TMath::Abs(meanselerr[0]/meansel[0]) + TMath::Abs(meanselerr[1]/meansel[1])) ); |
latex = new TLatex(); |
latex->SetTextSize(0.039); |
latex->DrawLatex(volt[0], 0.97*sep[0][sortindex[p-1]], ctemp); |
latex->DrawLatex(volt[0], 0.97*sep[0][sortindex[1]], ctemp); |
printf("#Delta_{p}(U) = (%.6lf #pm %.8lf)#timesU + (%.6lf #pm %.8lf)}{U_{0} = %.6lf #pm %.8lf\n", meansel[0], meanselerr[0], meansel[1], meanselerr[1], meansel[2], meansel[2]*(TMath::Abs(meanselerr[0]/meansel[0]) + TMath::Abs(meanselerr[1]/meansel[1])) ); |
} |
else |
1752,7 → 1912,6 |
double range[4]; |
TGraph2D *gScan2D; |
gScan2D = new TGraph2D(); |
gScan2D->SetName("surfscan"); |
range[0] = TMath::MinElement(nrfiles, surfx); |
range[1] = TMath::MaxElement(nrfiles, surfx); |
range[2] = TMath::MinElement(nrfiles, surfy); |
1855,8 → 2014,6 |
{ |
gCanvas->Modified(); |
gCanvas->Update(); |
|
UpdateIntegrateSurface(-1); |
} |
} |
} |