Using Anderson-Darling and Kolmogorov-Smirnov goodness of fit tests 1 sample test is performed comparing data with a log-normal distribution and a 2 sample test is done comparing two gaussian data sets.
#include <cassert>
double landau(
double x) {
}
void goftest() {
TF1 *
f1 =
new TF1(
"logNormal",
"ROOT::Math::lognormal_pdf(x,[0],[1])",0,500);
TH1D* h1smp =
new TH1D(
"h1smp",
"LogNormal distribution histogram", 100, 0, 500);
for (
UInt_t i = 0; i < nEvents1; ++i) {
}
Double_t A2_1 = goftest_1-> AndersonDarlingTest(
"t");
assert(A2_1 == A2_2);
Double_t pvalueAD_1 = goftest_1-> AndersonDarlingTest();
assert(pvalueAD_1 == pvalueAD_2);
delete goftest_1;
Double_t Dn_1 = goftest_1-> KolmogorovSmirnovTest(
"t");
assert(Dn_1 == Dn_2);
Double_t pvalueKS_1 = goftest_1-> KolmogorovSmirnovTest();
assert(pvalueKS_1 == pvalueKS_2);
#ifdef TEST_ERROR_MESSAGE
assert(A2 == pvalueKS);
#endif
sprintf(str1, "p-value for A-D 1-smp test: %f", pvalueAD_1);
sprintf(str2, "p-value for K-S 1-smp test: %f", pvalueKS_1);
TH1D* h2smps_1 =
new TH1D(
"h2smps_1",
"Gaussian distribution histograms", 100, 0, 500);
TH1D* h2smps_2 =
new TH1D(
"h2smps_2",
"Gaussian distribution histograms", 100, 0, 500);
for (
UInt_t i = 0; i < nEvents1; ++i) {
}
h2smps_1->
Scale(1. / nEvents1,
"width");
for (
UInt_t i = 0; i < nEvents2; ++i) {
}
h2smps_2->
Scale(1. / nEvents2,
"width");
assert(A2_1 == A2_2);
pvalueAD_1 = goftest_2-> AndersonDarling2SamplesTest();
assert(pvalueAD_1 == pvalueAD_2);
Dn_1 = goftest_2-> KolmogorovSmirnov2SamplesTest("t");
assert(Dn_1 == Dn_2);
pvalueKS_1 = goftest_2-> KolmogorovSmirnov2SamplesTest();
assert(pvalueKS_1 == pvalueKS_2);
#ifdef TEST_ERROR_MESSAGE
assert(A2 == pvalueKS);
#endif
sprintf(str1, "p-value for A-D 2-smps test: %f", pvalueAD_1);
sprintf(str2, "p-value for K-S 2-smps test: %f", pvalueKS_1);
pt2-> AddText(str2);
for (
UInt_t i = 0; i < nEvents3; ++i) {
}
pvalueAD_1 = goftest_3a-> AndersonDarlingTest();
pvalueAD_2 = (*goftest_3b)();
std::cout << " \n\nTEST with LANDAU distribution:\t";
if (
TMath::Abs(pvalueAD_1 - pvalueAD_2) > 1.E-1 * pvalueAD_2) {
std::cout << "FAILED " << std::endl;
Error(
"goftest",
"Error in comparing testing using Landau and Landau CDF");
std::cerr << " pvalues are " << pvalueAD_1 << " " << pvalueAD_2 << std::endl;
}
else
std::cout << "OK ( pvalues = " << pvalueAD_2 << " )" << std::endl;
}
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
R__EXTERN TRandom * gRandom
Functor1D class for one-dimensional functions.
GoFTest class implementing the 1 sample and 2 sample goodness of fit tests for uni-variate distributi...
@ kLogNormal
Gaussian distribution with default mean=0, sigma=1.
@ kKS
Anderson-Darling 2-Samples Test.
@ kKS2s
Kolmogorov-Smirnov Test.
@ kAD2s
Anderson-Darling Test. Default value.
@ kPDF
Input distribution is a CDF : cumulative distribution function.
void SetDistribution(EDistribution dist, const std::vector< double > &distParams={})
Sets the distribution for the predefined distribution types and optionally its parameters for 1-sampl...
void AndersonDarling2SamplesTest(Double_t &pvalue, Double_t &testStat) const
Performs the Anderson-Darling 2-Sample Test.
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
virtual void SetLineColor(Color_t lcolor)
Set the line color.
virtual void SetTextColor(Color_t tcolor=1)
Set the text color.
virtual void SetTextFont(Font_t tfont=62)
Set the text font.
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
void Draw(Option_t *option="") override
Draw this function with its current attributes.
virtual void SetParameters(const Double_t *params)
1-D histogram with a double per channel (see TH1 documentation)
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
void Draw(Option_t *option="") override
Draw this histogram with options.
virtual void Scale(Double_t c1=1, Option_t *option="")
Multiply this histogram by a constant c1.
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
The most important graphics class in the ROOT system.
void SetLogy(Int_t value=1) override
Set Lin/Log scale for Y.
TVirtualPad * cd(Int_t subpadnumber=0) override
Set Current pad.
A Pave (see TPave) with text, lines or/and boxes inside.
virtual TText * AddText(Double_t x1, Double_t y1, const char *label)
Add a new Text line to this pavetext at given coordinates.
void Draw(Option_t *option="") override
Draw this pavetext with its current attributes.
Random number generator class based on M.
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
double landau_pdf(double x, double xi=1, double x0=0)
Probability density function of the Landau distribution:
double lognormal_cdf(double x, double m, double s, double x0=0)
Cumulative distribution function of the lognormal distribution (lower tail).
Double_t Exp(Double_t x)
Returns the base-e exponential function of x, which is e raised to the power x.
T MinElement(Long64_t n, const T *a)
Returns minimum of array a of length n.
T MaxElement(Long64_t n, const T *a)
Returns maximum of array a of length n.
Double_t LandauI(Double_t x)
Returns the cumulative (lower tail integral) of the Landau distribution function at point x.
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.