// runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramwoRTS->pixeltimefiredsorted, "", 0);
// Integralberechnung vom Veto Peak
- runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogram->Seed, "GaussTail", true);
- runs[runi]->plot1DHistogram(runs[runi]->histogram->normalized, runs[runi]->histogram->normalized->Seed, "GaussTail", true);
+// runs[runi]->plot1DHistogram(runs[runi]->histogramwoRTS, runs[runi]->histogramwoRTS->Seed, "GaussTail", true);
+// runs[runi]->plot1DHistogram(runs[runi]->histogramwoRTS, runs[runi]->histogramwoRTS->Sum, "GaussTail", true);
- compareHistogramVector.push_back(runs[runi]->histogram->normalized->Seed);
- compareHistogramVector.push_back(runs[runi]->histogram->Seed);
+ // compareHistogramVector.push_back(runs[runi]->histogram->normalized->Seed);
+// compareHistogramVector.push_back(runs[runi]->histogram->normalized->Seed);
+// compareHistogramVector.push_back(runs[runi]->histogramthreshold->normalized->Seed);
+// compareHistogramVector.push_back(runs[runi]->histogramwoRTS->normalized->Seed);
+// compareHistogramVector.push_back(runs[runi]->histogramwoRTS->normalized->Sum);
+// compareHistogramVector.push_back(runs[runi]->histogramwoRTSAggresive->normalized->Seed);
+
+ runs[runi]->plot1DHistogram(runs[runi]->histogramwoRTS, runs[runi]->histogramwoRTS->Seed, "GaussTail", true, true, false, 500);
+ runs[runi]->plot1DHistogram(runs[runi]->histogramwoRTS, runs[runi]->histogramwoRTS->Sum, "GaussTail", true, true, false, 500);
+
+// compareHistogramVector2.push_back(runs[runi]->histogramwoRTS->normalized->Seed);
+// compareHistogramVector3.push_back(runs[runi]->histogramwoRTS->normalized->Sum);
+
+// compareHistogramVector3.push_back(runs[runi]->histogram->pixeltimefired);
+// compareHistogramVector3.push_back(runs[runi]->histogramwoRTS->pixeltimefired);
+// compareHistogramVector3.push_back(runs[runi]->histogramwoRTSAggresive->pixeltimefired);
+ //
+ // compareHistogramVector5.push_back(runs[runi]->histogram->pixeltimefiredsorted);
+
+// compareHistogramVector2.push_back(runs[runi]->histogramwoRTS->Seed);
+// compareHistogramVector3.push_back(runs[runi]->histogramwoRTS->Sum);
// runs[runi]->plot1DHistogram(runs[runi]->histogramwoRTS->normalized->calibrated, runs[runi]->histogramwoRTS->normalized->calibrated->Sum, "gaus", true);
//
// // compareHistogramClassVector.push_back(runs[runi]->histogram);
// // compareHistogramVector2.push_back(runs[runi]->histogram->LeakageCurrentInPixel);
-// compareHistogramVector3.push_back(runs[runi]->histogram->LeakageCurrentInPixelSorted);
-// compareHistogramVector4.push_back(runs[runi]->histogramwoRTS->LeakageCurrentInPixelSorted);
//
// compareHistogramVector5.push_back(runs[runi]->histogram->pixeltimefiredsorted);
// compareHistogramVector6.push_back(runs[runi]->histogramwoRTS->pixeltimefiredsorted);
curhistogramclone->SetAxisRange(0,lastbin1*1.1,"X");
// curhistogramclone->GetYaxis()->SetRangeUser(5,heighestval1*4);
- // gPad->SetLogy(1);
+ gPad->SetLogy(1);
owntitle->Clear();
owntitle->AddText(trimRunnumberAtBegin(curhistogramclone->GetName()));
owntitle->Draw("SAME");
calibrated->posSum = posSum * gain;
calibrated->posVeto = posVeto * gain;
calibrated->integralSeed = integralSeed;
+ calibrated->integralSeedErr = integralSeedErr;
calibrated->integralVeto = integralVeto;
calibrated->integralSum = integralVeto;
+ calibrated->integralSumErr = integralSumErr;
calibrated->posSeedPerc = posSeedPerc;
calibrated->sigmaSeedPerc = sigmaSeedPerc;
calibrated->avgNoise = avgNoise * gain;
calibrated->RTSpixel = RTSpixel;
calibrated->percentageofRTSpixel = percentageofRTSpixel;
calibrated->avgLeakageCurrentInChip = avgLeakageCurrentInChip * gain;
- calibrated->medianLeakageCurrent = medianLeakageCurrent * gain / 3.7 / pow10(-3) * (1.6*pow10(-18)) * pow10(15);
- calibrated->medianLeakageCurrentMinus = medianLeakageCurrentMinus * gain / 3.7 / pow10(-3) * (1.6*pow10(-18)) * pow10(15);
- calibrated->medianLeakageCurrentPlus = medianLeakageCurrentPlus * gain / 3.7 / pow10(-3) * (1.6*pow10(-18)) * pow10(15);
+ calibrated->medianLeakageCurrent = medianLeakageCurrent * gain / 3.7 / pow10(-3) * (1.6*pow10(-19)) * pow10(15);
+ calibrated->medianLeakageCurrentMinus = medianLeakageCurrentMinus * gain / 3.7 / pow10(-3) * (1.6*pow10(-19)) * pow10(15);
+ calibrated->medianLeakageCurrentPlus = medianLeakageCurrentPlus * gain / 3.7 / pow10(-3) * (1.6*pow10(-19)) * pow10(15);
calibrated->iscalibrated = true;
normalized->posSum = posSum;
normalized->posVeto = posVeto;
normalized->integralSeed = integralSeed/frames_found*pow10(6);
+ normalized->integralSeedErr = integralSeedErr/frames_found*pow10(6);
normalized->integralVeto= integralVeto/frames_found*pow10(6);
normalized->integralVeto= integralVeto/frames_found*pow10(6);
normalized->integralSum = integralSum/frames_found*pow10(6);
+ normalized->integralSumErr = integralSumErr/frames_found*pow10(6);
normalized->posSeedPerc = posSeedPerc;
normalized->sigmaSeedPerc = sigmaSeedPerc;
normalized->noisethresholdborder = noisethresholdborder;
}
-Double_t* HistogramType::FitPerform(TH1F* histogrampointer, TString fitFuncType, Bool_t verbose, Float_t fitstart) {
+//Double_t* HistogramType::FitPerform(TH1F* histogrampointer, TString fitFuncType, TFitResult* fit_result, Bool_t verbose, Float_t fitstart) {
+Double_t* HistogramType::FitPerform(TH1F* histogrampointer, TString fitFuncType, TFitResult* fit_result_giveback, Bool_t verbose, Float_t fitstart) {
+ TFitResultPtr fit_result_ptr;
Double_t* parameters = 0;
Float_t posMax = 0;
Float_t integralPeak = 0;
+ Double_t integralPeakError = 0;
Float_t posMax2 = 0;
Float_t posMaxValHist = histogrampointer->GetXaxis()->GetXmax();
// empirical value for USB system is 50, should be overwritten in 99% of cases
if (fitFuncType.EqualTo("gaus"))
{
+ parameters = (Double_t *)calloc(10, sizeof(Double_t));
+
histogrampointer->GetXaxis()->SetRange(histogrampointer->GetXaxis()->FindBin(noiseborder),histogrampointer->GetXaxis()->FindBin(posMaxValHist)); // look only for maxima with x greater than noiseborder, cut away noise
// Int_t xValMax = histogrampointer->GetBinCenter(histogrampointer->GetMaximumBin());
TF1* fitFunc = new TF1("fitFunc",fitFuncType,noiseborder,posMaxValHist);
int fittries = 0;
do {
// cout << fittries << ": New range: " << histogrampointer->GetXaxis()->GetBinCenter(min+fittries*min/20) << " to " << posMaxValHist << endl;
- histogrampointer->Fit(fitFunc, "N,Q,W", "", histogrampointer->GetXaxis()->GetBinCenter(min+fittries*min/20), posMaxValHist);
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "NQWS", "", histogrampointer->GetXaxis()->GetBinCenter(min+fittries*min/20), posMaxValHist);
// cout << colorred << gMinuit->fCstatu.Data() << endlr;
} while (gMinuit->fCstatu.Contains("FAILED") && fittries++ < 8);
posMax = histogrampointer->GetXaxis()->GetBinCenter(histogrampointer->GetMaximumBin()); // Methode 1
- integralPeak = histogrampointer->Integral(min+fittries*min/20, histogrampointer->GetXaxis()->FindBin(posMaxValHist), "width");
+ integralPeak = histogrampointer->IntegralAndError(min+fittries*min/20, histogrampointer->GetXaxis()->FindBin(posMaxValHist), integralPeakError, "width");
// if (verbose)
// {
// cout << coloryellow << "min " << histogrampointer->GetXaxis()->GetBinCenter( min ) << endlr;
if (verbose)
fitFunc->DrawCopy("same");
} else {
- histogrampointer->Fit(fitFunc, "N,Q,W", "", noiseborder, posMaxValHist);
- noiseborder = FindBorderToPeak(histogrampointer, noiseborder,fitFunc->GetParameter(1), verbose); // starting point of histogram integration
- integralPeak = histogrampointer->Integral(histogrampointer->GetXaxis()->FindBin(noiseborder), histogrampointer->GetXaxis()->FindBin(posMaxValHist), "width");
+ // modified the gaussian approuch to be more powerful, escpecially the
+ // integral is calculated in a +- 2 sigma region.
+
+ const Double_t def_amplitude=15;
+ const Double_t def_peakcenter=histogrampointer->GetBinCenter(histogrampointer->GetMaximumBin());
+ const Double_t def_gausssig=21;
+ // set start values for some parameters
+ fitFunc->SetParName(0,"amplitude of peak");
+ fitFunc->SetParameter(0,def_amplitude);
+ fitFunc->SetParName(1,"peak centroid");
+ fitFunc->SetParameter(1,def_peakcenter);
+ fitFunc->SetParLimits(1,def_peakcenter*0.8,def_peakcenter*1.2);
+ fitFunc->SetParName(2,"Gaussian sigma");
+ fitFunc->SetParameter(2,def_gausssig);
+ fitFunc->SetParLimits(2,0,100);
+
+ fitFunc->SetLineWidth(4);
+ fitFunc->SetLineColor(kGreen);
+
+ int fittries = 0;
+ Bool_t failed = false;
+ do {
+ if (failed)
+ {
+ fitFunc->SetParameter(0,def_amplitude*(1.0-0.1*++fittries));
+ fitFunc->SetParameter(1,def_peakcenter);
+ fitFunc->SetParameter(2,def_gausssig);
+ } else fittries = 100;
+ failed = false;
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "NMWQS", "", noiseborder, posMaxValHist);
+ if (gMinuit == nullptr) {
+ failed = true;
+ } else
+ {
+ if (gMinuit->fCstatu == "NOT POSDEF" || gMinuit->fCstatu.Contains("FAILED")) {
+ failed = true;
+ }
+ }
+ } while (fittries < 10);
+ // get parameters
+ for (Int_t pari=0; pari<3; pari++) {
+ //cout << colorcyan << fitFunc->GetParameter(pari) << endlr;
+ parameters[pari]=fitFunc->GetParameter(pari);
+ }
+ if (failed || ( parameters[2]/parameters[0] > 100 ) )
+ {
+ fitFunc->SetParameter(0,def_amplitude);
+ fitFunc->SetParameter(1,def_peakcenter);
+ fitFunc->SetParameter(2,def_gausssig);
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "NMWQS", "", def_peakcenter*0.7, def_peakcenter*1.3);
+ }
+ if (failed)
+ {
+ cout << colorred << " Could not find " << histogrampointer->GetName() << " peak" << endlr;
+ parameters = (Double_t *)calloc(8, sizeof(Double_t));
+ return parameters;
+ }
+
+ // get parameters
+ for (Int_t pari=0; pari<3; pari++) {
+ //cout << colorcyan << fitFunc->GetParameter(pari) << endlr;
+ parameters[pari]=fitFunc->GetParameter(pari);
+ }
+ // https://suchideas.com/articles/maths/applied/histogram-errors/#The_Sensible_Way
+ // make an error estimate, in case of rare events one can use the poisson distribution
+ // error bars become +- sqrt(n) for each bin, where n is the bin content.
+
+ parameters[7] = parameters[1] - 2*parameters[2] ; // starting point of histogram integration
+ parameters[8] = parameters[1] + 2*parameters[2] ; // end point of histogram integration
+ integralPeak = histogrampointer->IntegralAndError(histogrampointer->GetXaxis()->FindBin( parameters[7] ), histogrampointer->GetXaxis()->FindBin( parameters[8]), integralPeakError, "width"); // integral value of histogram (NOT fit), normaliezed with bin size!
posMax = fitFunc->GetMaximumX(); // Methode 2
fitFunc->SetLineStyle(1); // normal for the following fits
if (verbose)
fitFunc->Draw("SAME");
}
- parameters = (Double_t *)calloc(5, sizeof(Double_t));
+
+
+ // parameters are
+ // parameters[0]: amplitude of peak
+ // parameters[1]: x position of peak
+ // parameters[2]: sigma of gaussian
+ // parameters[3]: integral of peak
+ // parameters[6]: integral of peak
+ // parameters[7]: start of integration
+ // parameters[8]: end of integratio
for (Int_t pari=0; pari<3; pari++) {
//cout << colorcyan << fitFunc->GetParameter(pari) << endlr;
parameters[pari]=fitFunc->GetParameter(pari);
}
- parameters[3]=integralPeak;
- parameters[4] = noiseborder;
+ parameters[6]=integralPeak;
+ parameters[9]=integralPeakError;
Float_t sigma = fitFunc->GetParameter(2);
posMax2 = fitFunc->GetMaximumX(); // Methode 2
// if (verbose)
cout << coloryellow << " Please check " << histogrampointer->GetName() << " peak position: " << colorreset << colorwhite << posMax << endlr;
}
// }
- }
+ }
+ if (verbose) {
+ cout << colorcyan << "Depletion voltage: " << labbook->depletionV << endlr;
+ cout << colorcyan << "Integral from bin : " << histogrampointer->FindBin(parameters[7]) << " to " << histogrampointer->GetXaxis()->FindBin(parameters[8]) << endlr;
+ cout << colorcyan << "Integral from val : " << parameters[7] << " to " << parameters[8] << endlr;
+ cout << colorcyan << "Integral value without bg: " << parameters[6] << endlr;
+ cout << colorcyan << "Integral error: " << parameters[9] << endlr;
+ cout << colorcyan << "Gauss Sigma: " << parameters[2] << endlr;
+ cout << colorcyan << "Number of events : " << frames_found << endlr;
+
+ (fit_result_ptr.Get())->Print("V");
+
+// cout << fitFunc->GetParError(0) << endlr; // retrieve the error for the parameter 0
+ }
// TString legendEntry = TString(Form("%s","Gaus fit"));
// TLegend* leg = new TLegend(0.5,0.5,0.89,0.89);//(0.6,0.7,0.89,0.89);
Float_t fitMax3 = 1000;
Float_t minFitMax = 1000;
Float_t maxFitMax = 1000;
- histogrampointer->Fit(fitFunc, "N,Q,W", "", noiseborder, posMaxValHist);
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "N,Q,W,S", "", noiseborder, posMaxValHist);
for (Int_t pari=0; pari<3; pari++)
parameters[pari]=fitFunc->GetParameter(pari);
// posMax = fitFunc->GetParameter(1); // new method, just use MPV of fit, as noise border is well known. Compare with old method anyway and warn if something suspecious.
// cout << colorred << "noiseborder: " << noiseborder << endlr;
fitMax1 = fitFunc->GetMaximumX();
// fitFunc->DrawCopy("same");
- histogrampointer->Fit(fitFunc, "N,Q,W", "", noiseborder, fitMax1*1.1);
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "N,Q,W,S", "", noiseborder, fitMax1*1.1);
fitMax2 = fitFunc->GetMaximumX();
fitFunc->SetLineColor(kBlue);
fitFunc->SetLineStyle(2); // dashed
// fitFunc->DrawCopy("same");
- histogrampointer->Fit(fitFunc, "N,Q,W", "", fitMax1*0.9, posMaxValHist);
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "N,Q,W,S", "", fitMax1*0.9, posMaxValHist);
// histogrampointer->Fit(fitFunc, "N,Q,W", "", fitMax1, histogrampointer->GetBinCenter(bini));
fitMax3 = fitFunc->GetMaximumX();
fitFunc->SetLineColor(kGreen);
// fitFcn->FixParameter(3,-22.43016284);
//histogrampointer->Fit("fitFcn","V+","ep");
- histogrampointer->Fit("fitFcn","Q,W","ep");
+ fit_result_ptr = histogrampointer->Fit("fitFcn","Q,W,S","ep");
parameters = (Double_t *)calloc(4, sizeof(Double_t));
for (Int_t pari=0; pari<4; pari++)
parameters[pari]=fitFcn->GetParameter(pari);
// histogrampointer->GetXaxis()->SetRangeUser(noiseborder,posMaxValHist);
TF1* fitFunc = new TF1("fitFunc",GaussTail,noiseborder,posMaxValHist,6);
- parameters = (Double_t *)calloc(9, sizeof(Double_t));
+ parameters = (Double_t *)calloc(10, sizeof(Double_t));
const Double_t def_amplitude=459.951;
const Double_t def_peakcenter=histogrampointer->GetBinCenter(histogrampointer->GetMaximumBin());
// cout << colorcyan << "def_peakcenter: " << def_peakcenter << endlr;
fitFunc->SetParLimits(1,def_peakcenter*0.8,def_peakcenter*1.2);
fitFunc->SetParName(2,"Gaussian sigma");
fitFunc->SetParameter(2,def_gausssig);
+ fitFunc->SetParLimits(2,2,50);
fitFunc->SetParLimits(2,0,10000);
fitFunc->SetParName(3,"Distance from Gauss");
fitFunc->SetParameter(3,def_distgauss);
fitFunc->SetParameter(5,def_bgoffs);
fitFunc->SetParLimits(5,0,def_bgoffs*4);
- parameters = (Double_t *)calloc(9, sizeof(Double_t));
+ // TODO: This fix disables the background
+ fitFunc->FixParameter(4,0);
+ // fitFunc->FixParameter(5,0);
+
int fittries = 0;
Bool_t failed = false;
do {
failed = false;
// cout << fittries << ": New range: " << histogrampointer->GetXaxis()->GetBinCenter(min+fittries*min/20) << " to " << posMaxValHist << endl;
- histogrampointer->Fit(fitFunc, "N,M,W,Q", "", noiseborder, posMaxValHist);
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "N,M,W,Q,S", "", noiseborder, posMaxValHist);
// cout << colorcyan << " AFTER fit " << endlr;
if (gMinuit == nullptr) {
failed = true;
fitFunc->SetParameter(2,def_gausssig);
fitFunc->SetParameter(4,def_bgslope);
fitFunc->SetParameter(5,def_bgoffs);
+
+ // TODO: This fix disables the background
+ fitFunc->FixParameter(4,0);
+// fitFunc->FixParameter(5,0);
}
else
fittries = 100;
fitFunc->SetParameter(3,def_distgauss);
fitFunc->SetParameter(4,def_bgslope);
fitFunc->SetParameter(5,def_bgoffs);
- histogrampointer->Fit(fitFunc, "N,M,W,Q", "", noiseborder, posMaxValHist);
+
+ // TODO: This fix disables the background
+ fitFunc->FixParameter(4,0);
+// fitFunc->FixParameter(5,0);
+
+ fit_result_ptr = histogrampointer->Fit(fitFunc, "N,M,W,Q,S", "", noiseborder, posMaxValHist);
if (gMinuit == nullptr) {
failed = true;
} else
if (failed)
{
cout << colorred << " Could not find " << histogrampointer->GetName() << " peak" << endlr;
- parameters = (Double_t *)calloc(8, sizeof(Double_t));
+ parameters = (Double_t *)calloc(10, sizeof(Double_t));
return parameters;
}
// fitFunc->FixParameter(1,parameters[1]+histogrampointer->GetBinWidth(0));
// histogrampointer->Fit(fitFunc, "N,M,W,Q", "", noiseborder, posMaxValHist);
+
+ // parameters are
+ // parameters[0]: amplitude of peak
+ // parameters[1]: x position of peak
+ // parameters[2]: sigma of gaussian
+ // parameters[6]: integral of peak
+ // parameters[7]: start of integral
+ // parameters[8]: end of integratio
+ // parameters[9]: error of integral
for (Int_t pari=0; pari<6; pari++)
parameters[pari]=fitFunc->GetParameter(pari);
//parameters[1] = parameters[1] - histogrampointer->GetBinWidth(0)*2;
//parameters[7] = FindBorderToPeak(histogrampointer, noiseborder,def_peakcenter, verbose); // starting point of histogram integration
- parameters[7] = parameters[1] - 2*parameters[2] ; // starting point of histogram integration
- parameters[8] = parameters[1] + 2*parameters[2] ; // end point of histogram integration
- parameters[6] = histogrampointer->Integral(histogrampointer->GetXaxis()->FindBin( parameters[1] - 2*parameters[2] ), histogrampointer->GetXaxis()->FindBin( parameters[1] + 2*parameters[2]), "width"); // integral value of histogram (NOT fit)
+ parameters[7] = parameters[1] - 3*parameters[2] ; // starting point of histogram integration
+ parameters[8] = parameters[1] + 3*parameters[2] ; // end point of histogram integration
+ integralPeak = histogrampointer->IntegralAndError(histogrampointer->GetXaxis()->FindBin( parameters[7] ), histogrampointer->GetXaxis()->FindBin( parameters[8]), integralPeakError, "width"); // integral value of histogram (NOT fit),
+
+ parameters[6]=integralPeak;
+ parameters[9]=integralPeakError;
TF1 *bgfct = new TF1("f1","[0] +[1]*x",0,posMaxValHist);
bgfct->SetParameters(parameters[5],parameters[4]);
cout << colorcyan << "Integral value without bg: " << parameters[6] << endlr;
cout << colorcyan << "Integral value with bg: " << parameters[6] + integralbg << endlr;
cout << colorcyan << "Gauss Sigma: " << parameters[2] << endlr;
- cout << colorcyan << "Number of events : " << frames_found << endlr;
+ cout << colorcyan << "Number of events : " << frames_found << endlr;
+
+ (fit_result_ptr.Get())->Print("V");
+
}
// DEBUG
fitFunc->Draw("SAME");
}
+
+
+ fit_result_giveback = fit_result_ptr.Get();
+// cout << "Correct address:" << endl;
+// cout << fit_result_ptr.Get() << endl;
+//
+// fit_result_ptr_giveback = new TFitResultPtr();
+//
+// fit_result_ptr_giveback = &fit_result_ptr;
+// cout << "Gave pointer address:" << endl;
+// cout << fit_result_ptr_giveback->Get() << endl;
+// TFitResult resultnew = new TFitResult(fit_result_ptr.Get());
+// cout << "saved address:" << endl;
+// cout << &resultnew << endl;
return parameters;
}
#include <TMath.h>
#include <TStyle.h>
#include <TMinuit.h>
+#include <TFitResult.h>
#include "help.h"
/// fitted position of the most probable value of the seed spectrum
Float_t posSeed=0;
- /// fintegral over the Seed peak, 4 sigma area
+ /// fintegral over the Seed peak, 2 sigma area
Double_t integralSeed=0;
+ /// fintegral error over the Seed peak, poisson dist. assumed -> sqrt(N) for each bin
+ Double_t integralSeedErr=0;
/// fitted position of the most probable value in the over cluster summed spectrum
Float_t posSum=0;
/// fintegral over the sum peak, normalized to number of events
Double_t integralSum=0;
+ /// fintegral error over the Seed peak, poisson dist. assumed -> sqrt(N) for each bin
+ Double_t integralSumErr=0;
/// fitted position of the calibration peak of Fe55-beta-photons in the seed spectrum, from a run best suited to the current
Float_t posVeto=0;
/// fintegral over the veto peak, normalized to number of events
/**
* @brief rescales one specific histogram bin content from ADU to electrons */
void calibrate2DHistogramCounts(TH2F* &histogrampointernew, TH2F* &histogrampointerold );
+
+
+
/**
* @brief intern function to calculate and plot fit curve to given histogram
*
* @return peak position of the fit
*
*/
- Double_t* FitPerform(TH1F* histogrampointer, TString fitFuncType, Bool_t verbose = 0, Float_t fitstart = -1);
+ Double_t* FitPerform(TH1F* histogrampointer, TString fitFuncType, TFitResult* r = 0, Bool_t verbose = 0, Float_t fitstart = -1);
+ //Double_t* FitPerform(TH1F* histogrampointer, TString fitFuncType, TFitResultPtr* r = 0, Bool_t verbose = 0, Float_t fitstart = -1);
/** @brief calculates Charge Collection Efficiency if Fe55 run */
Bool_t calculteCCE(Bool_t verbose=false);
/**
};
-#endif
\ No newline at end of file
+#endif
// cout << endl;
// if cluster charge > clusternoise * const
- if(fOrderCode=="FSBB" || fOrderCode=="Mi19")
+ //if(fOrderCode=="FSBB" || fOrderCode=="Mi19") // if you need CCE_4, uncomment this and comment the line below for Mi19
+ if(fOrderCode=="FSBB")
{
if (1.0*chargesumincluster4 > noisesumincluster4*6.0)
fFrameInfo.pixelthreshold[fHits] = Hitlist[hit];
/// Size of cluster. Don't change without code modification
const Int_t clustersize = 5;
};
-#endif
\ No newline at end of file
+#endif
} else {
storepathRAWLinux=Form("/d/garlic/%s/%d",labbook.chipGen.Data(),runnumber); // default empty path
}
+
storepathRootLinux = Form("/d/garlic/%s/rootFiles",labbook.chipGen.Data()); // default empty path
labbook.posVetoDB = (rowsql->GetField(10) != NULL)?atof(rowsql->GetField(10)):-1;
labbook.posSeedDB = (rowsql->GetField(11) != NULL)?atof(rowsql->GetField(11)):-1;
HistogramClassVector.push_back(histogramwoRTSthreshold);
- // histogramwoRTSAggresive = new HistogramType(" more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, this, rootcolors[plotStyle], rootlinestyle[plotStyle]);
-// histogramwoRTSAggresive->maskRTSpixel = true;
-// histogramwoRTSAggresive->RTSthreshold = 1.5;
-// HistogramClassVector.push_back(histogramwoRTSAggresive);
+ histogramwoRTSAggresive = new HistogramType(" more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, this, rootcolors[plotStyle], rootlinestyle[plotStyle]);
+ histogramwoRTSAggresive->maskRTSpixel = true;
+ histogramwoRTSAggresive->RTSthreshold = 1.5;
+ HistogramClassVector.push_back(histogramwoRTSAggresive);
// histogramwoRTSAggresivethreshold = new HistogramType(" Threshold, more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, this, rootcolors[plotStyle], rootlinestyle[plotStyle] );
// histogramwoRTSAggresivethreshold->maskRTSpixel = true;
// histogramwoRTSAggresivethreshold->RTSthreshold = 1.0;
systemparam systemparamPXI (800*16,800,75,150,150);
systemparam systemparamFSBB (2800,2800/4,25,10,100);
systemparam systemparamUSBMi19 (400/*maxbin*/,400/1/*nbins*/, 25/*vetothreshold*/, 10/*maxbinnoise*/, 100/*nbinsnoise*/);
- systemparam systemparamPipper2 (1000,500,15,10,100);
+ systemparam systemparamPipper2 (1000,200,15,10,100);
if (labbook.system.EqualTo("USB") && labbook.chipGen.EqualTo("Mi34") )
cursystemparam = systemparamUSB;
else if (labbook.system.EqualTo("USB") && labbook.chipGen.EqualTo("FSBB") )
}
}
+
// calculate the mean firing times and bin fireing times in a propability histogram, first approach
for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > 0) {
// HistogramTypepointer->pixeltimefiredsorted->Fill(2000);
}
}
-// TH1F* pixeltimefiredWithOverflow = ShowOverflow(HistogramTypepointer->pixeltimefiredsorted);
+ meanpixeltimesfired /= numberofactivepixel; // Very rough estimate of a mean firing time
+
+ // vrey rough estimate on standard deviation
+ for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
+ if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > 0)
+ stdeviation += pow(HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired,2);
+ }
+ stdeviation /= numberofactivepixel;
+ stdeviation = sqrt(stdeviation);
+ cout << colorcyan << "Mean firing rate first approximation: " << meanpixeltimesfired << endlr;
+ cout << colorcyan << "Deviation value first approximation: " << stdeviation << endlr;
+ cout << colorcyan << "number of considered pixel: " << numberofactivepixel << endlr;
+ cout << colorcyan << endlr;
+
+ // better estimate on average firing time
+ Double_t meanpixeltimesfired2, stdeviation2;
+ uint numberofactivepixel2;
+ Bool_t RTSpixel;
+ u_int numberofconsideredpixel;
+ u_int poscounter;
+
+
+ for (UInt_t i=0; i < 3; ++i) // 3 times remove most fired pixel
+ {
+ meanpixeltimesfired2 = 0;
+ numberofactivepixel2=0;
+ HistogramTypepointer->RTSpixel.clear();
+ cout << "i " << i << endl;
+ cout << "HistogramTypepointer->pixeltimefired->GetBinContent(0) " << HistogramTypepointer->pixeltimefired->GetBinContent(0) << endl;
+ cout << "meanpixeltimesfired " << meanpixeltimesfired << endl;
+ cout << "HistogramTypepointer->RTSthreshold " << HistogramTypepointer->RTSthreshold << endl;
+ cout << "stdeviation " << stdeviation << endl;
+ for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
+ if (abs(HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired) > HistogramTypepointer->RTSthreshold * stdeviation) {
+ HistogramTypepointer->RTSpixel.push_back(pixeli);
+ RTSpixelHits+=HistogramTypepointer->pixeltimefired->GetBinContent(pixeli);
+ }
+ else
+ {
+ meanpixeltimesfired2 += HistogramTypepointer->pixeltimefired->GetBinContent(pixeli);
+ numberofactivepixel2++;
+ }
+
+ }
+ if (numberofactivepixel2 > 0) {
+ meanpixeltimesfired2 /= numberofactivepixel2;
+ meanpixeltimesfired = meanpixeltimesfired2;
+ }
// Double_t const probabilities[] = {1-erf(HistogramTypepointer->RTSthreshold/sqrt(2)), 0.5, erf(HistogramTypepointer->RTSthreshold/sqrt(2))}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17};
// HistogramTypepointer->pixeltimefiredsorted->GetQuantiles( 3, leakagequantiles, probabilities);
+
// HistogramTypepointer->medianLeakageCurrent = leakagequantiles[1];
// HistogramTypepointer->medianLeakageCurrentPlus = leakagequantiles[2] - leakagequantiles[1];
// HistogramTypepointer->medianLeakageCurrentMinus = leakagequantiles[1] - leakagequantiles[0];
// cout << "last value above 1 in overflow TH1F: " << pixeltimefiredWithOverflow->GetBinCenter(pixeltimefiredWithOverflow->FindLastBinAbove(4)) << endl;
// cout << colorcyan << "Deviation value after one more loop: " << stdeviation << endlr;
- HistogramTypepointer->RTSpixel.clear();
- for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
- if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > leakagequantiles[2]) {
-// cout << coloryellow << numberToString(5) << " mal drüber, Böser pixel #" << pixeli << " Wert: " << HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) << endlr;
- HistogramTypepointer->RTSpixel.push_back(pixeli);
- }
- }
-// cout << colorred << " number of RTS pixel in quantile run: " << HistogramTypepointer->RTSpixel.size() << endlr;
+// HistogramTypepointer->RTSpixel.clear();
+// for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
+// if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > leakagequantiles[2]) {
+// // cout << coloryellow << numberToString(5) << " mal drüber, Böser pixel #" << pixeli << " Wert: " << HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) << endlr;
+// HistogramTypepointer->RTSpixel.push_back(pixeli);
+// }
+// }
+// cout << colorred << " number of RTS pixel in quantile run: " << HistogramTypepointer->RTSpixel.size() << endlr;
- meanpixeltimesfired /= numberofactivepixel; // Very rough estimate of a mean firing time
+// meanpixeltimesfired /= HistogramTypepointer->pixeltimefired->GetNbinsX(); // Very rough estimate of a mean firing time
// meanpixeltimesfired = leakagequantiles[1];
// cout << colorcyan << "Mean value own: " << meanpixeltimesfired << endlr;
- // vrey rough estimate on standard deviation
- for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
- if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > 0)
- stdeviation += pow(HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired,2);
- }
- stdeviation /= numberofactivepixel;
- stdeviation = sqrt(stdeviation);
-// cout << colorcyan << "Deviation value own 1: " << stdeviation << endlr;
-
- // better estimate on average firing time
- Double_t meanpixeltimesfired2 = 0;
- uint numberofactivepixel2=0;
- HistogramTypepointer->RTSpixel.clear();
- for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
- if (abs(HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired) > 2 * stdeviation) {
- HistogramTypepointer->RTSpixel.push_back(pixeli);
- }
- else
- {
- meanpixeltimesfired2 += HistogramTypepointer->pixeltimefired->GetBinContent(pixeli);
- numberofactivepixel2++;
- }
- }
- meanpixeltimesfired2 /= numberofactivepixel2;
-// cout << colorred << " number of RTS pixel in first run: " << HistogramTypepointer->RTSpixel.size() << endlr;
-// cout << colorcyan << " mean value new: " << meanpixeltimesfired2 << endl;
-
-// Estimate new standard deviation
- Bool_t RTSpixel =false;
- u_int numberofconsideredpixel=0;
- u_int poscounter = 0;
- stdeviation = 0;
- for (u_int pixeli=0; (int)pixeli<cursensorinfo.columns*cursensorinfo.rows ; pixeli++) {
- if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > 0) {
- RTSpixel = false;
- for (u_int RTSpixeli=poscounter; RTSpixeli < HistogramTypepointer->RTSpixel.size(); RTSpixeli++) {
- if (pixeli == HistogramTypepointer->RTSpixel[RTSpixeli]) {
- poscounter = RTSpixeli;
- RTSpixel = true;
+ // Estimate new standard deviation
+ RTSpixel =false;
+ numberofconsideredpixel=0;
+ poscounter = 0;
+ stdeviation2 = 0;
+ for (UInt_t pixeli=0; pixeli < (UInt_t)HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
+ if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > 0) {
+ RTSpixel = false;
+ for (u_int RTSpixeli=0; RTSpixeli < (u_int) HistogramTypepointer->RTSpixel.size(); RTSpixeli++) {
+ if (pixeli == HistogramTypepointer->RTSpixel[RTSpixeli]) {
+ poscounter = RTSpixeli;
+ RTSpixel = true;
+ }
+ }
+ if (!RTSpixel) {
+ stdeviation2 += pow(HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired2,2);
+ numberofconsideredpixel++;
}
}
- if (!RTSpixel) {
- stdeviation += pow(HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired2,2);
- numberofconsideredpixel++;
- }
- }
- }
- stdeviation /= numberofconsideredpixel;
- stdeviation = sqrt(stdeviation);
-// cout << colorcyan << "Deviation value own 2: " << stdeviation << endlr;
-
-
- // better estimate on RTS pixel
- meanpixeltimesfired = 0;
- numberofconsideredpixel = 0;
- HistogramTypepointer->RTSpixel.clear();
- for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
- totalHits+=HistogramTypepointer->pixeltimefired->GetBinContent(pixeli);
- if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli)-meanpixeltimesfired2 > HistogramTypepointer->RTSthreshold * stdeviation) {
- HistogramTypepointer->RTSpixel.push_back(pixeli);
- RTSpixelHits+=HistogramTypepointer->pixeltimefired->GetBinContent(pixeli);
-// cout << coloryellow << numberToString(HistogramTypepointer->RTSthreshold) << " mal drüber, Böser pixel #" << pixeli << " Wert: " << HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) << endlr;
}
- else
- {
- meanpixeltimesfired += HistogramTypepointer->pixeltimefired->GetBinContent(pixeli);
- numberofconsideredpixel++;
+ if (numberofconsideredpixel > 0) {
+ stdeviation2 /= numberofconsideredpixel;
+ stdeviation2 = sqrt(stdeviation);
+ stdeviation = stdeviation2;
}
+ cout << colorred << " number of RTS pixel in " << i << " run: " << HistogramTypepointer->RTSpixel.size() << endlr;
+ cout << colorred << " number of considered pixel: " << numberofconsideredpixel << endlr;
+ cout << colorcyan << " mean value: " << meanpixeltimesfired << endl;
+ cout << colorcyan << " Deviation value: " << stdeviation << endlr;
+
+
+
+ // Double_t const probabilities[] = {1-erf(HistogramTypepointer->RTSthreshold/sqrt(2)), 0.5, erf(HistogramTypepointer->RTSthreshold/sqrt(2))}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17};
+ // HistogramTypepointer->pixeltimefiredsorted->GetQuantiles( 3, leakagequantiles, probabilities);
+
+ // HistogramTypepointer->medianLeakageCurrent = leakagequantiles[1];
+ // HistogramTypepointer->medianLeakageCurrentPlus = leakagequantiles[2] - leakagequantiles[1];
+ // HistogramTypepointer->medianLeakageCurrentMinus = leakagequantiles[1] - leakagequantiles[0];
+ // cout << "Median is at: " << HistogramTypepointer->medianLeakageCurrent << endl;
+ // cout << "Upper limit is at: " << leakagequantiles[2] << endl;
+ // cout << "Lower limit is at: " << leakagequantiles[0] << endl;
+ // cout << "Plus: " << HistogramTypepointer->medianLeakageCurrentPlus << endl;
+ // cout << "minus: " << HistogramTypepointer->medianLeakageCurrentMinus << endl;
+
+ // cout << "last value above 1: " << HistogramTypepointer->pixeltimefiredsorted->GetBinCenter(HistogramTypepointer->pixeltimefiredsorted->FindLastBinAbove(4)) << endl;
+ // cout << "last value above 1 in overflow TH1F: " << pixeltimefiredWithOverflow->GetBinCenter(pixeltimefiredWithOverflow->FindLastBinAbove(4)) << endl;
+ // cout << colorcyan << "Deviation value after one more loop: " << stdeviation << endlr;
+
+ // HistogramTypepointer->RTSpixel.clear();
+ // for (Int_t pixeli=0; pixeli < HistogramTypepointer->pixeltimefired->GetNbinsX(); pixeli++) {
+ // if (HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) > leakagequantiles[2]) {
+ // // cout << coloryellow << numberToString(5) << " mal drüber, Böser pixel #" << pixeli << " Wert: " << HistogramTypepointer->pixeltimefired->GetBinContent(pixeli) << endlr;
+ // HistogramTypepointer->RTSpixel.push_back(pixeli);
+ // }
+ // }
+ // cout << colorred << " number of RTS pixel in quantile run: " << HistogramTypepointer->RTSpixel.size() << endlr;
+
+
+ // meanpixeltimesfired /= HistogramTypepointer->pixeltimefired->GetNbinsX(); // Very rough estimate of a mean firing time
+ // meanpixeltimesfired = leakagequantiles[1];
+ // cout << colorcyan << "Mean value own: " << meanpixeltimesfired << endlr;
}
- meanpixeltimesfired /= numberofconsideredpixel;
-// cout << colorred << " number of RTS pixel in second run: " << HistogramTypepointer->RTSpixel.size() << endlr;
-// cout << colorcyan << " mean value new new: " << meanpixeltimesfired << endl;
+
+
HistogramTypepointer->percentageofRTSpixel = (HistogramTypepointer->RTSpixel.size()*1.0)/(numberofactivepixel*1.0)*100.0;
cout << " cutted away evertyhing with more then " << std::right<< Form("%.1f",HistogramTypepointer->RTSthreshold * stdeviation+meanpixeltimesfired2) << " hits" << endlr;
} else {
histogramclassToUseForDB = histogramwoRTSthreshold;
}
-
-
-
+
string sqlupdatequery = "";
constructUpdateString(&sqlupdatequery, "Gain", histogramclassToUseForDB->normalized->gain);
constructUpdateString(&sqlupdatequery, "SumPeak", histogramclassToUseForDB->normalized->posSum, 4);
constructUpdateString(&sqlupdatequery, "AvgF0", labbook.averageF0, 6);
constructUpdateString(&sqlupdatequery, "SigmaF0", labbook.sigmaF0, 6);
constructUpdateString(&sqlupdatequery, "SeedIntegralMinusBG2Sigma", histogramclassToUseForDB->normalized->integralSeed, 10);
+ constructUpdateString(&sqlupdatequery, "SeedIntegralMinusBG2SigmaErr", histogramclassToUseForDB->normalized->integralSeedErr, 5);
constructUpdateString(&sqlupdatequery, "VetoPeak", histogramclassToUseForDB->normalized->posVeto, 4);
constructUpdateString(&sqlupdatequery, "VetoIntegral", histogramclassToUseForDB->normalized->integralVeto, 10);
constructUpdateString(&sqlupdatequery, "SumIntegral", histogramclassToUseForDB->normalized->integralSum, 10);
+ constructUpdateString(&sqlupdatequery, "SumIntegralErr", histogramclassToUseForDB->normalized->integralSumErr, 5);
if (histogramclassToUseForDB->normalized->calibrated != 0)
// if (histogramthreshold->normalized->calibrated->avgNoise < 100)
constructUpdateString(&sqlupdatequery, "Avg.Noise", histogramclassToUseForDB->normalized->calibrated->avgNoise);
Bool_t Run::binLeakageCurrent(HistogramType* oneHistogramClass)
{
Double_t const probabilities[] = {0.158655254, 0.5, 1-0.158655254}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17};
-// Double_t const probabilities[] = {0.25, 0.5, 0.75}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17};
+// Double_t const probabilities[] = {0.25, 0.5, 0.75}; // sigma/2 from gaus to the left and to the right, an area of 68% must be within the interesting range //{0.17, 0.5, 1-0.17};
Double_t pedestals [cursensorinfo.columns*cursensorinfo.rows];
for (Int_t pixeli=0; pixeli<cursensorinfo.columns*cursensorinfo.rows ; pixeli++) { // loop over all pixel
pedestals[pixeli]=0;
// sum histogram
pixelSum = 0;
notSeedSum = 0;
- if(labbook.chipGen=="FSBB" || labbook.chipGen=="Mi19")
+ //if(labbook.chipGen=="FSBB" || labbook.chipGen=="Mi19")
+ if(labbook.chipGen=="FSBB")
{
Float_t clusterArray[processed->clustersize*processed->clustersize];// temp variable clusterArray necessary, because Sort only accepts 1-dim arrays
Int_t index[processed->clustersize*processed->clustersize];
if (processed->fFrameInfo.pixel[hiti] == histogramwoRTS->RTSpixel[RTSpixeli])
{
RTSpixel = true;
-// cout << "not binned! RTS pixel #" << processed->fFrameInfo.pixel[hiti] << " with charge: " << processed->fFrameInfo.p[12][hiti] << endl;
+// cout << "not binned! RTS pixel #" << processed->fFrameInfo.pixel[hiti] << " with charge: " << processed->fFrameInfo.p[12][hiti] << endl;
}
}
if (!RTSpixel) {
if (processed->fFrameInfo.pixel[hiti] == histogramwoRTSAggresive->RTSpixel[RTSpixeli])
{
RTSpixel = true;
+// cout << "not binned! RTS pixel #" << processed->fFrameInfo.pixel[hiti] << " with charge: " << processed->fFrameInfo.p[12][hiti] << endl;
// cout << "not binned! RTS pixel #" << processed->fFrameInfo.pixel[hiti] << " with charge: " << processed->fFrameInfo.p[12][hiti] << endl;
}
}
if (!RTSpixel) {
+// cout << "binned! pixel #" << processed->fFrameInfo.pixel[hiti] << " with charge: " << processed->fFrameInfo.p[12][hiti] << endl;
histogramwoRTSAggresive->Seed->Fill(processed->fFrameInfo.p[12][hiti]);
histogramwoRTSAggresive->Sum->Fill(pixelSum);
histogramwoRTSAggresive->SeedPerc->Fill(processed->fFrameInfo.p[12][hiti]/pixelSum*100);
} // end loop over hits in frame
}
} // end loop over all frames
+ cout << colorred << " histogramwoRTSAggresive->Seed size: " << histogramwoRTSAggresive->Seed->GetEntries() << endlr;
// gROOT->SetBatch(kTRUE);
for (vector<HistogramType*>::iterator curHistogramClass = HistogramClassVector.begin(); curHistogramClass != HistogramClassVector.end(); curHistogramClass++) {
- Double_t* parameters = (Double_t *)calloc(10, sizeof(Double_t)); // allocate 10 parameters for safety, maximum 9 used at the moment
+ Double_t* parameters = (Double_t *)calloc(11, sizeof(Double_t)); // allocate 11 parameters for safety, maximum 10 used at the moment
(*curHistogramClass)->FindNoisethresholdborder((*curHistogramClass)->Seed, false, false);
if (labbook.runnumber == 343056)
(*curHistogramClass)->noisethresholdborder = 34;
(*curHistogramClass)->integralVeto = parameters[6];
}
if (labbook.chipGen.EqualTo("Pipper2"))
- parameters = (*curHistogramClass)->FitPerform((*curHistogramClass)->Seed, "GaussTail");
+ parameters = (*curHistogramClass)->FitPerform((*curHistogramClass)->Seed, "GaussTail", 0, false, 500);
else
parameters = (*curHistogramClass)->FitPerform((*curHistogramClass)->Seed, "landau");
(*curHistogramClass)->integralSeed = parameters[6];
(*curHistogramClass)->posSeed = parameters[1];
- parameters = (*curHistogramClass)->FitPerform((*curHistogramClass)->Sum, "gaus");
+ (*curHistogramClass)->integralSeedErr = parameters[9];
+ parameters = (*curHistogramClass)->FitPerform((*curHistogramClass)->Sum, "GaussTail", 0, false, 500); //TODO change back to gauss
(*curHistogramClass)->posSum = parameters[1];
- (*curHistogramClass)->integralSum = parameters[3];
+ (*curHistogramClass)->integralSum = parameters[6];
+ (*curHistogramClass)->integralSumErr = parameters[9];
parameters =(*curHistogramClass)->FitPerform((*curHistogramClass)->SeedPerc, "landau");
(*curHistogramClass)->posSeedPerc = parameters[1];
// create lorentz fit of a slice
TH1F* xslicetroughcluster = new TH1F(Form("%s_xslice", HistogramTypepointer->histAvgCluster->GetTitle()),"X slice",HistogramTypepointer->histAvgCluster->GetXaxis()->GetNbins(),HistogramTypepointer->histAvgCluster->GetXaxis()->GetBinLowEdge(1),HistogramTypepointer->histAvgCluster->GetXaxis()->GetBinUpEdge(HistogramTypepointer->histAvgCluster->GetXaxis()->GetNbins()));
- Double_t* xslicetroughclusterpar = (Double_t *)calloc(4, sizeof(Double_t)); // allocate 6 parameters for safety, maximum 4 used at the moment
+ Double_t* xslicetroughclusterpar = (Double_t *)calloc(11, sizeof(Double_t)); // allocate 11 parameters for safety, maximum 10 used at the moment
Int_t middlebin = (int)(HistogramTypepointer->histAvgCluster->GetXaxis()->GetNbins()/2.0 + 0.5);
for (Int_t bini=1; bini <= HistogramTypepointer->histAvgCluster->GetXaxis()->GetNbins(); bini++)
xslicetroughcluster->SetBinContent(bini,HistogramTypepointer->histAvgCluster->GetBinContent(bini,middlebin));
TH1F* yslicetroughcluster = new TH1F(Form("%s_yslice", HistogramTypepointer->histAvgCluster->GetTitle()),"Y slice",HistogramTypepointer->histAvgCluster->GetYaxis()->GetNbins(),HistogramTypepointer->histAvgCluster->GetYaxis()->GetBinLowEdge(1),HistogramTypepointer->histAvgCluster->GetYaxis()->GetBinUpEdge(HistogramTypepointer->histAvgCluster->GetYaxis()->GetNbins()));
- Double_t* yslicetroughclusterpar = (Double_t *)calloc(4, sizeof(Double_t)); // allocate 6 parameters for safety, maximum 4 used at the moment
+ Double_t* yslicetroughclusterpar = (Double_t *)calloc(11, sizeof(Double_t)); // allocate 11 parameters for safety, maximum 10 used at the moment
middlebin = (int)(HistogramTypepointer->histAvgCluster->GetYaxis()->GetNbins()/2.0 + 0.5);
for (Int_t bini=1; bini <= HistogramTypepointer->histAvgCluster->GetYaxis()->GetNbins(); bini++)
yslicetroughcluster->SetBinContent(bini,HistogramTypepointer->histAvgCluster->GetBinContent(middlebin,bini));
if (onehistogram != nullptr)
{
Int_t random = random1->Rndm()*100000;
+ TFitResult* fit_result;
+// TFitResultPtr* fit_result_ptr;
+ fit_result = new TFitResult();
+// fit_result_ptr = new TFitResultPtr();
TString canvastitle = Form("%s %s", onehistogram->GetName(), runcode.Data());
TString canvasname = Form("%s %s %d", onehistogram->GetName(), runcode.Data(), random);
TCanvas* canvas = new TCanvas(canvasname, canvastitle, 900, 700);
onehistogram=ShowOverflow(onehistogram);
onehistogram->Draw();
- Double_t* parameters = (Double_t *)calloc(10, sizeof(Double_t)); // allocate 10 parameters for safety, maximum 9 used at the moment
+ Double_t* parameters = (Double_t *)calloc(11, sizeof(Double_t)); // allocate 10 parameters for safety, maximum 9 used at the moment
Float_t maxValue = 0;
if (fitFuncType!="")
{
- parameters = HistogramTypepointer->FitPerform(onehistogram, fitFuncType, verbose, fitstart);
+ parameters = HistogramTypepointer->FitPerform(onehistogram, fitFuncType, fit_result, verbose, fitstart);
maxValue = parameters[1];
plotVerticalLine(onehistogram, maxValue);
}
// canvas->Draw();
// canvas->Update();
- if (fitFuncType=="GaussTail")
+ if (fitFuncType=="GaussTail" || fitFuncType=="gaus")
{
Float_t integralstart = parameters[7];
Float_t integralend = parameters[8];
cout << "Distance under : " << parameters[3] << endl;
cout << "background slope : " << parameters[4] << endl;
cout << "background offset: " << parameters[5] << endl;
+ cout << "FWHM : " << parameters[2]*sqrt(8*log(2)) << endl;
+// cout << "got address2:" << endl;
+// cout << fit_result << endl;
+// fit_result->Print("V");
+
//TF1 *f1 = new TF1("f1","[0] +[1]*x +gaus(2)",0,5);
TF1 *f1 = new TF1("f1","[0] +[1]*x",0,1000);
f1->SetParameters(parameters[5],parameters[4]);
f1->Draw("SAME");
}
}
- if (fitFuncType=="gaus")
- {
- Float_t integralstart = parameters[4];
- Float_t integralval = parameters[3];
- TString integrallbl = Form("Integral: %.0f", integralval);
-// cout << colorcyan << " " << HistogramTypepointer->histogramdescription << endlr;
-// cout << colorcyan << ": " << parameters[4] << endlr;
- plotVerticalLineHeigher(onehistogram, integralstart);
- plotTextAtVal(onehistogram, maxValue, integrallbl);
- }
+// if (fitFuncType=="gaus")
+// {
+// Float_t integralstart = parameters[4];
+// Float_t integralval = parameters[3];
+// TString integrallbl = Form("Integral: %.0f", integralval);
+// // cout << colorcyan << " " << HistogramTypepointer->histogramdescription << endlr;
+// // cout << colorcyan << ": " << parameters[4] << endlr;
+// plotVerticalLineHeigher(onehistogram, integralstart);
+// plotTextAtVal(onehistogram, maxValue, integrallbl);
+// }
// canvas -> SaveAs( savepathresults + "/" + runcode + " " + onehistogram->GetName() + ".eps");
#include <TH2F.h>
#include <TF1.h>
#include <TImageDump.h>
+#include <TFitResult.h>
#include <TMath.h>
#include <TF1.h>
}
//####################################################################
-#endif
\ No newline at end of file
+#endif
jspc29 / radhard.git / commitdiff