From 1206251231696359270a260508551e044f3af33a Mon Sep 17 00:00:00 2001
From: Stefan <stefan.wunsch@student.kit.edu>
Date: Tue, 1 Sep 2015 12:46:09 +0200
Subject: include random.h in swig; add qa_random testcase
---
gnuradio-runtime/python/gnuradio/gr/qa_random.py | 99 ++++++++++++++++++++++++
1 file changed, 99 insertions(+)
create mode 100644 gnuradio-runtime/python/gnuradio/gr/qa_random.py
(limited to 'gnuradio-runtime/python')
diff --git a/gnuradio-runtime/python/gnuradio/gr/qa_random.py b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
new file mode 100644
index 0000000000..39d75f3afa
--- /dev/null
+++ b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
@@ -0,0 +1,99 @@
+#!/usr/bin/env python
+#
+# Copyright 2006,2007,2010 Free Software Foundation, Inc.
+#
+# This file is part of GNU Radio
+#
+# GNU Radio is free software; you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation; either version 3, or (at your option)
+# any later version.
+#
+# GNU Radio is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with GNU Radio; see the file COPYING. If not, write to
+# the Free Software Foundation, Inc., 51 Franklin Street,
+# Boston, MA 02110-1301, USA.
+#
+
+from gnuradio import gr, gr_unittest
+import numpy as np
+from scipy.stats import norm, laplace
+
+class test_random(gr_unittest.TestCase):
+
+ # Disclaimer
+ def test_0(self):
+ print 'NOTE: Following tests are not statistically significant! Check out fulltest_random.py for full testing.'
+ self.assertEqual(1,1)
+
+ # Check for range [0,1) of uniform distributed random numbers and print minimal and maximal value
+ def test_1(self):
+ print '# TEST 1'
+ print 'Uniform distributed numbers: Range'
+ num_tests = 10000
+ values = np.zeros(num_tests)
+ rndm = gr.random()
+ for k in range(num_tests):
+ values[k] = rndm.ran1()
+ for value in values:
+ self.assertLess(value, 1)
+ self.assertGreaterEqual(value, 0)
+ print 'Uniform random numbers (num/min/max):', num_tests, min(values), max(values)
+
+ # Check uniformly distributed random numbers on uniformity (without assert, only printing)
+ def test_2(self):
+ print '# TEST 2'
+ print 'Uniform random numbers: Distribution'
+ num_tests = 10000
+ num_bins = 11
+ values = np.zeros(num_tests)
+ rndm = gr.random()
+ for k in range(num_tests):
+ values[k] = rndm.ran1()
+ bins = np.linspace(0,1,num_bins) # These are the bin edges!
+ hist = np.histogram(values,bins)
+ print 'Lower edge bin / upper edge bin / count / expected'
+ for k in range(len(hist[0])):
+ print hist[1][k], hist[1][k+1], hist[0][k], float(num_tests)/(num_bins-1)
+
+ # Check distribution of normally (gaussian, mean=0, variance=1) distributed random numbers (no assert)
+ def test_3(self):
+ print '# TEST 3'
+ print 'Normal random numbers: Distribution'
+ num_tests = 10000
+ num_bins = 11
+ hist_range = [-5,5]
+ values = np.zeros(num_tests)
+ rndm = gr.random()
+ for k in range(num_tests):
+ values[k] = rndm.gasdev()
+ bins = np.linspace(hist_range[0],hist_range[1],num_bins)
+ hist = np.histogram(values,bins)
+ print 'Lower edge bin / upper edge bin / count / expected'
+ for k in range(len(hist[0])):
+ print hist[1][k], hist[1][k+1], hist[0][k], float(norm.cdf(hist[1][k+1])-norm.cdf(hist[1][k]))*num_tests
+
+ # Check distribution of laplacian (mean=0, variance=1) distributed random numbers (no assert)
+ def test_4(self):
+ print '# TEST 4'
+ print 'Laplacian random numbers: Distribution'
+ num_tests = 100000
+ num_bins = 11
+ hist_range = [-5,5]
+ values = np.zeros(num_tests)
+ rndm = gr.random()
+ for k in range(num_tests):
+ values[k] = rndm.laplacian()
+ bins = np.linspace(hist_range[0],hist_range[1],num_bins)
+ hist = np.histogram(values,bins)
+ print 'Lower edge bin / upper edge bin / count / expected'
+ for k in range(len(hist[0])):
+ print hist[1][k], hist[1][k+1], hist[0][k], float(laplace.cdf(hist[1][k+1])-laplace.cdf(hist[1][k]))*num_tests
+
+if __name__ == '__main__':
+ gr_unittest.run(test_random, "test_random.xml")
--
cgit v1.2.3
From ebad2162b40b2b144449ff2927dbe89c683c4972 Mon Sep 17 00:00:00 2001
From: Stefan <stefan.wunsch@student.kit.edu>
Date: Tue, 1 Sep 2015 13:21:51 +0200
Subject: fix wrong laplacian random numbers and add testcase
---
gnuradio-runtime/include/gnuradio/random.h | 21 ++++++++--
gnuradio-runtime/lib/math/random.cc | 23 ++---------
gnuradio-runtime/python/gnuradio/gr/qa_random.py | 51 +++++++++++++++---------
3 files changed, 54 insertions(+), 41 deletions(-)
(limited to 'gnuradio-runtime/python')
diff --git a/gnuradio-runtime/include/gnuradio/random.h b/gnuradio-runtime/include/gnuradio/random.h
index e01fcb7be7..454c88ee6c 100644
--- a/gnuradio-runtime/include/gnuradio/random.h
+++ b/gnuradio-runtime/include/gnuradio/random.h
@@ -60,18 +60,33 @@ namespace gr {
void reseed(long seed);
/*!
- * \brief uniform random deviate in the range [0.0, 1.0)
+ * \brief Uniform random numbers in the range [0.0, 1.0)
*/
float ran1();
/*!
- * \brief normally distributed deviate with zero mean and variance 1
+ * \brief Normally distributed random numbers (Gaussian distribution with zero mean and variance 1)
*/
float gasdev();
+ /*!
+ * \brief Laplacian distributed random numbers with zero mean and variance 1
+ */
float laplacian();
- float impulse(float factor);
+
+ /*!
+ * \brief Rayleigh distributed random numbers (zero mean and variance 1 for the underlying Gaussian distributions)
+ */
float rayleigh();
+
+ /*!
+ * \brief FIXME: add description
+ */
+ float impulse(float factor);
+
+ /*!
+ * \brief Normally distributed random numbers with zero mean and variance 1 on real and imaginary part. This results in a Rayleigh distribution for the amplitude and an uniform distribution for the phase.
+ */
gr_complex rayleigh_complex();
};
diff --git a/gnuradio-runtime/lib/math/random.cc b/gnuradio-runtime/lib/math/random.cc
index 7170f27edf..ecf70e0f9a 100644
--- a/gnuradio-runtime/lib/math/random.cc
+++ b/gnuradio-runtime/lib/math/random.cc
@@ -130,18 +130,12 @@ namespace gr {
return d_gset;
}
- /*
- * Copied from The KC7WW / OH2BNS Channel Simulator
- * FIXME Need to check how good this is at some point
- */
float
random::laplacian()
{
- float z = ran1();
- if(z < 0.5)
- return log(2.0 * z) / M_SQRT2;
- else
- return -log(2.0 * (1.0 - z)) / M_SQRT2;
+ float z = ran1()-0.5;
+ if(z>0) return -log(1-2*z);
+ else return log(1+2*z);
}
/*
@@ -161,10 +155,6 @@ namespace gr {
/*
* Complex rayleigh is really gaussian I and gaussian Q
- * It can also be generated by real rayleigh magnitude and
- * uniform random angle
- * Adapted from The KC7WW / OH2BNS Channel Simulator
- * FIXME Need to check how good this is at some point
*/
gr_complex
random::rayleigh_complex()
@@ -172,13 +162,6 @@ namespace gr {
return gr_complex(gasdev(),gasdev());
}
- /* Other option
- mag = rayleigh();
- ang = 2.0 * M_PI * RNG();
- *Rx = rxx * cos(z);
- *Iy = rxx * sin(z);
- */
-
float
random::rayleigh()
{
diff --git a/gnuradio-runtime/python/gnuradio/gr/qa_random.py b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
index 39d75f3afa..ee4018327b 100644
--- a/gnuradio-runtime/python/gnuradio/gr/qa_random.py
+++ b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
@@ -22,78 +22,93 @@
from gnuradio import gr, gr_unittest
import numpy as np
-from scipy.stats import norm, laplace
+from scipy.stats import norm, laplace, rayleigh
class test_random(gr_unittest.TestCase):
+ num_tests = 10000
+
# Disclaimer
def test_0(self):
- print 'NOTE: Following tests are not statistically significant! Check out fulltest_random.py for full testing.'
+ print 'NOTE: Following tests are not statistically significant!'
+ print 'Realisations per test:',self.num_tests
self.assertEqual(1,1)
# Check for range [0,1) of uniform distributed random numbers and print minimal and maximal value
def test_1(self):
print '# TEST 1'
print 'Uniform distributed numbers: Range'
- num_tests = 10000
- values = np.zeros(num_tests)
+ values = np.zeros(self.num_tests)
rndm = gr.random()
- for k in range(num_tests):
+ for k in range(self.num_tests):
values[k] = rndm.ran1()
for value in values:
self.assertLess(value, 1)
self.assertGreaterEqual(value, 0)
- print 'Uniform random numbers (num/min/max):', num_tests, min(values), max(values)
+ print 'Uniform random numbers (num/min/max):', self.num_tests, min(values), max(values)
# Check uniformly distributed random numbers on uniformity (without assert, only printing)
def test_2(self):
print '# TEST 2'
print 'Uniform random numbers: Distribution'
- num_tests = 10000
num_bins = 11
- values = np.zeros(num_tests)
+ values = np.zeros(self.num_tests)
rndm = gr.random()
- for k in range(num_tests):
+ for k in range(self.num_tests):
values[k] = rndm.ran1()
bins = np.linspace(0,1,num_bins) # These are the bin edges!
hist = np.histogram(values,bins)
print 'Lower edge bin / upper edge bin / count / expected'
for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(num_tests)/(num_bins-1)
+ print hist[1][k], hist[1][k+1], hist[0][k], float(self.num_tests)/(num_bins-1)
# Check distribution of normally (gaussian, mean=0, variance=1) distributed random numbers (no assert)
def test_3(self):
print '# TEST 3'
print 'Normal random numbers: Distribution'
- num_tests = 10000
num_bins = 11
hist_range = [-5,5]
- values = np.zeros(num_tests)
+ values = np.zeros(self.num_tests)
rndm = gr.random()
- for k in range(num_tests):
+ for k in range(self.num_tests):
values[k] = rndm.gasdev()
bins = np.linspace(hist_range[0],hist_range[1],num_bins)
hist = np.histogram(values,bins)
print 'Lower edge bin / upper edge bin / count / expected'
for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(norm.cdf(hist[1][k+1])-norm.cdf(hist[1][k]))*num_tests
+ print hist[1][k], hist[1][k+1], hist[0][k], float(norm.cdf(hist[1][k+1])-norm.cdf(hist[1][k]))*self.num_tests
# Check distribution of laplacian (mean=0, variance=1) distributed random numbers (no assert)
def test_4(self):
print '# TEST 4'
print 'Laplacian random numbers: Distribution'
- num_tests = 100000
num_bins = 11
hist_range = [-5,5]
- values = np.zeros(num_tests)
+ values = np.zeros(self.num_tests)
rndm = gr.random()
- for k in range(num_tests):
+ for k in range(self.num_tests):
values[k] = rndm.laplacian()
bins = np.linspace(hist_range[0],hist_range[1],num_bins)
hist = np.histogram(values,bins)
print 'Lower edge bin / upper edge bin / count / expected'
for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(laplace.cdf(hist[1][k+1])-laplace.cdf(hist[1][k]))*num_tests
+ print hist[1][k], hist[1][k+1], hist[0][k], float(laplace.cdf(hist[1][k+1])-laplace.cdf(hist[1][k]))*self.num_tests
+
+ # Check distribution of laplacian (mean=0, variance=1) distributed random numbers (no assert)
+ def test_5(self):
+ print '# TEST 5'
+ print 'Rayleigh random numbers: Distribution'
+ num_bins = 11
+ hist_range = [0,10]
+ values = np.zeros(self.num_tests)
+ rndm = gr.random()
+ for k in range(self.num_tests):
+ values[k] = rndm.rayleigh()
+ bins = np.linspace(hist_range[0],hist_range[1],num_bins)
+ hist = np.histogram(values,bins)
+ print 'Lower edge bin / upper edge bin / count / expected'
+ for k in range(len(hist[0])):
+ print hist[1][k], hist[1][k+1], hist[0][k], float(rayleigh.cdf(hist[1][k+1])-rayleigh.cdf(hist[1][k]))*self.num_tests
if __name__ == '__main__':
gr_unittest.run(test_random, "test_random.xml")
--
cgit v1.2.3
From e172d55e2bfe63e5b99e7d432f9386c37bff8c84 Mon Sep 17 00:00:00 2001
From: Stefan <stefan.wunsch@student.kit.edu>
Date: Tue, 1 Sep 2015 15:09:17 +0200
Subject: add boost.random as random number generator
---
gnuradio-runtime/include/gnuradio/random.h | 21 +++--
gnuradio-runtime/lib/math/random.cc | 106 +++++++++--------------
gnuradio-runtime/python/gnuradio/gr/qa_random.py | 19 ++++
3 files changed, 73 insertions(+), 73 deletions(-)
(limited to 'gnuradio-runtime/python')
diff --git a/gnuradio-runtime/include/gnuradio/random.h b/gnuradio-runtime/include/gnuradio/random.h
index 454c88ee6c..f15b12a492 100644
--- a/gnuradio-runtime/include/gnuradio/random.h
+++ b/gnuradio-runtime/include/gnuradio/random.h
@@ -37,6 +37,8 @@ static const int RANDOM_MAX = 2147483647; // 2^31-1
#endif /* RANDOM_MAX */
#include <stdlib.h>
+#include <boost/random.hpp>
+#include <ctime>
namespace gr {
@@ -47,17 +49,22 @@ namespace gr {
class GR_RUNTIME_API random
{
protected:
- static const int NTAB = 32;
long d_seed;
- long d_iy;
- long d_iv[NTAB];
- int d_iset;
- float d_gset;
+ bool d_gauss_stored;
+ float d_gauss_value;
+
+ boost::mt19937 *d_rng; // mersenne twister as random number generator
+ boost::uniform_real<float> *d_uniform; // choose uniform distribution, default is [0,1)
+ boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > *d_generator;
public:
- random(long seed=3021);
+ random(unsigned int seed=0);
+ ~random();
- void reseed(long seed);
+ /*!
+ * \brief Change the seed for the initialized number generator
+ */
+ void reseed(unsigned int seed);
/*!
* \brief Uniform random numbers in the range [0.0, 1.0)
diff --git a/gnuradio-runtime/lib/math/random.cc b/gnuradio-runtime/lib/math/random.cc
index ecf70e0f9a..76de3bdde5 100644
--- a/gnuradio-runtime/lib/math/random.cc
+++ b/gnuradio-runtime/lib/math/random.cc
@@ -44,90 +44,67 @@
namespace gr {
-#define IA 16807
-#define IM 2147483647
-#define AM (1.0/IM)
-#define IQ 127773
-#define IR 2836
-#define NDIV (1+(IM-1)/NTAB)
-#define EPS 1.2e-7
-#define RNMX (1.0-EPS)
+ random::random(unsigned int seed)
+ {
+ d_gauss_stored = false; // set gasdev (gauss distributed numbers) on calculation state
+
+ // Setup random number generator
+ d_rng = new boost::mt19937;
+ reseed(seed); // set seed for random number generator
+ d_uniform = new boost::uniform_real<float>;
+ d_generator = new boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > (*d_rng,*d_uniform); // create number generator in [0,1) from boost.random
+ }
- random::random(long seed)
+ random::~random()
{
- reseed(seed);
+ delete d_rng;
+ delete d_uniform;
+ delete d_generator;
}
+ /*
+ * Seed is initialized with time if the given seed is 0. Otherwise the seed is taken directly. Sets the seed for the random number generator.
+ */
void
- random::reseed(long seed)
+ random::reseed(unsigned int seed)
{
- d_seed = seed;
- d_iy = 0;
- for(int i = 0; i < NTAB; i++)
- d_iv[i] = 0;
- d_iset = 0;
- d_gset = 0;
+ if(seed==0) d_seed = static_cast<unsigned int>(std::time(0)); // FIXME: add seed method correctly
+ else d_seed = seed;
+ d_rng->seed(d_seed);
}
/*
- * This looks like it returns a uniform random deviate between 0.0 and 1.0
- * It looks similar to code from "Numerical Recipes in C".
+ * Returns uniformly distributed numbers in [0,1) taken from boost.random using a Mersenne twister
*/
float
random::ran1()
{
- int j;
- long k;
- float temp;
-
- if(d_seed <= 0 || !d_iy) {
- if(-d_seed < 1)
- d_seed=1;
- else
- d_seed = -d_seed;
- for(j=NTAB+7;j>=0;j--) {
- k=d_seed/IQ;
- d_seed=IA*(d_seed-k*IQ)-IR*k;
- if(d_seed < 0)
- d_seed += IM;
- if(j < NTAB)
- d_iv[j] = d_seed;
- }
- d_iy=d_iv[0];
- }
- k=(d_seed)/IQ;
- d_seed=IA*(d_seed-k*IQ)-IR*k;
- if(d_seed < 0)
- d_seed += IM;
- j=d_iy/NDIV;
- d_iy=d_iv[j];
- d_iv[j] = d_seed;
- temp=AM * d_iy;
- if(temp > RNMX)
- temp = RNMX;
- return temp;
+ return (*d_generator)();
}
/*
* Returns a normally distributed deviate with zero mean and variance 1.
- * Also looks like it's from "Numerical Recipes in C".
+ * Used is the Marsaglia polar method.
+ * Every second call a number is stored because the transformation works only in pairs. Otherwise half calculation is thrown away.
*/
float
random::gasdev()
{
- float fac,rsq,v1,v2;
- d_iset = 1 - d_iset;
- if(d_iset) {
- do {
- v1=2.0*ran1()-1.0;
- v2=2.0*ran1()-1.0;
- rsq=v1*v1+v2*v2;
- } while(rsq >= 1.0 || rsq == 0.0);
- fac= sqrt(-2.0*log(rsq)/rsq);
- d_gset=v1*fac;
- return v2*fac;
- }
- return d_gset;
+ if(d_gauss_stored){ // just return the stored value if available
+ d_gauss_stored = false;
+ return d_gauss_value;
+ }
+ else{ // generate a pair of gaussian distributed numbers
+ float x,y,s;
+ do{
+ x = 2.0*ran1()-1.0;
+ y = 2.0*ran1()-1.0;
+ s = x*x+y*y;
+ }while(not(s<1.0));
+ d_gauss_stored = true;
+ d_gauss_value = x*sqrt(-2.0*log(s)/s);
+ return y*sqrt(-2.0*log(s)/s);
+ }
}
float
@@ -153,9 +130,6 @@ namespace gr {
return z;
}
- /*
- * Complex rayleigh is really gaussian I and gaussian Q
- */
gr_complex
random::rayleigh_complex()
{
diff --git a/gnuradio-runtime/python/gnuradio/gr/qa_random.py b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
index ee4018327b..ef85337ea6 100644
--- a/gnuradio-runtime/python/gnuradio/gr/qa_random.py
+++ b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
@@ -110,5 +110,24 @@ class test_random(gr_unittest.TestCase):
for k in range(len(hist[0])):
print hist[1][k], hist[1][k+1], hist[0][k], float(rayleigh.cdf(hist[1][k+1])-rayleigh.cdf(hist[1][k]))*self.num_tests
+ # Check seeds (init with time and seed as fix number, no assert)
+ def test_6(self):
+ print '# TEST 6'
+ rndm0 = gr.random(0); # init with time
+ rndm1 = gr.random(42); # init with fix seed
+ num = 5
+
+ print 'Some random numbers in [0,1), should change every run:'
+ for k in range(num):
+ print rndm0.ran1(),
+ print ' '
+
+ print 'Some random numbers in [0,1), should be the same every run:'
+ for k in range(num):
+ print rndm1.ran1(),
+ print '== '
+ print '0.374540120363 0.796543002129 0.950714290142 0.183434784412 0.731993913651'
+
+
if __name__ == '__main__':
gr_unittest.run(test_random, "test_random.xml")
--
cgit v1.2.3
From 30965aab6c6b28fabb1403c635d016fc532caa0d Mon Sep 17 00:00:00 2001
From: Stefan <stefan.wunsch@student.kit.edu>
Date: Tue, 1 Sep 2015 17:01:29 +0200
Subject: add test-case for reseed feature
---
gnuradio-runtime/include/gnuradio/random.h | 4 ++-
gnuradio-runtime/python/gnuradio/gr/qa_random.py | 38 +++++++++++++++++++-----
2 files changed, 33 insertions(+), 9 deletions(-)
(limited to 'gnuradio-runtime/python')
diff --git a/gnuradio-runtime/include/gnuradio/random.h b/gnuradio-runtime/include/gnuradio/random.h
index f15b12a492..3ecdf6a2b7 100644
--- a/gnuradio-runtime/include/gnuradio/random.h
+++ b/gnuradio-runtime/include/gnuradio/random.h
@@ -36,6 +36,8 @@
static const int RANDOM_MAX = 2147483647; // 2^31-1
#endif /* RANDOM_MAX */
+// FIXME: RANDOM_MAX still necessary? Some test-cases (use grep, e.g. in gr-filter or gr-atsc) use this constant, but combine it with the random() function from the std namespace. That should not the way to use this.
+
#include <stdlib.h>
#include <boost/random.hpp>
#include <ctime>
@@ -62,7 +64,7 @@ namespace gr {
~random();
/*!
- * \brief Change the seed for the initialized number generator
+ * \brief Change the seed for the initialized number generator. seed = 0 initializes the random number generator with the system time. Note that a fast initialization of various instances can result in the same seed.
*/
void reseed(unsigned int seed);
diff --git a/gnuradio-runtime/python/gnuradio/gr/qa_random.py b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
index ef85337ea6..c0d9a7f34c 100644
--- a/gnuradio-runtime/python/gnuradio/gr/qa_random.py
+++ b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
@@ -23,6 +23,7 @@
from gnuradio import gr, gr_unittest
import numpy as np
from scipy.stats import norm, laplace, rayleigh
+#from time import sleep
class test_random(gr_unittest.TestCase):
@@ -110,24 +111,45 @@ class test_random(gr_unittest.TestCase):
for k in range(len(hist[0])):
print hist[1][k], hist[1][k+1], hist[0][k], float(rayleigh.cdf(hist[1][k+1])-rayleigh.cdf(hist[1][k]))*self.num_tests
- # Check seeds (init with time and seed as fix number, no assert)
+ # Check seeds (init with time and seed as fix number)
def test_6(self):
print '# TEST 6'
- rndm0 = gr.random(0); # init with time
- rndm1 = gr.random(42); # init with fix seed
num = 5
print 'Some random numbers in [0,1), should change every run:'
+ rndm0 = gr.random(0); # init with time
+ # NOTE: the sleep increases the executiont time massively, remove assert for convenience
+ #sleep(1)
+ #rndm1 = gr.random(0); # init with fix seed
for k in range(num):
- print rndm0.ran1(),
+ x = rndm0.ran1();
+ print x,
+ # y = rndm1.ran1();
+ # print x, '!=', y
+ # self.assertNotEqual(x,y)
print ' '
- print 'Some random numbers in [0,1), should be the same every run:'
+ print 'Some random numbers in [0,1) (seed two instances), should be the same every run:'
+ rndm0 = gr.random(42); # init with time
+ rndm1 = gr.random(42); # init with fix seed
for k in range(num):
- print rndm1.ran1(),
- print '== '
- print '0.374540120363 0.796543002129 0.950714290142 0.183434784412 0.731993913651'
+ x = rndm0.ran1();
+ y = rndm1.ran1();
+ print x, '=', y
+ self.assertEqual(x,y)
+ print 'Some random numbers in [0,1) (reseed one instance), should be the same every run:'
+ x = np.zeros(num)
+ y = np.zeros(num)
+ rndm0 = gr.random(42); # init with time
+ for k in range(num):
+ x[k] = rndm0.ran1();
+ rndm1.reseed(43); # init with fix seed
+ for k in range(num):
+ y[k] = rndm0.ran1();
+ for k in range(num):
+ print x[k], '!=', y[k]
+ self.assertNotEqual(x[k],y[k])
if __name__ == '__main__':
gr_unittest.run(test_random, "test_random.xml")
--
cgit v1.2.3
From 44fb1cb0482fa778c8e652164551711818db5476 Mon Sep 17 00:00:00 2001
From: Stefan <stefan.wunsch@student.kit.edu>
Date: Fri, 4 Sep 2015 11:22:13 +0200
Subject: redo qa_random without print statements and scipy; add stand-alone
evaluation script in gnuradio-runtime/apps
---
gnuradio-runtime/apps/evaluation_random_numbers.py | 139 +++++++++++++++++++++
gnuradio-runtime/python/gnuradio/gr/qa_random.py | 109 ++--------------
2 files changed, 148 insertions(+), 100 deletions(-)
create mode 100644 gnuradio-runtime/apps/evaluation_random_numbers.py
(limited to 'gnuradio-runtime/python')
diff --git a/gnuradio-runtime/apps/evaluation_random_numbers.py b/gnuradio-runtime/apps/evaluation_random_numbers.py
new file mode 100644
index 0000000000..069493c73e
--- /dev/null
+++ b/gnuradio-runtime/apps/evaluation_random_numbers.py
@@ -0,0 +1,139 @@
+#!/usr/bin/env python
+#
+# Copyright 2015 Free Software Foundation, Inc.
+#
+# This file is part of GNU Radio
+#
+# GNU Radio is free software; you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation; either version 3, or (at your option)
+# any later version.
+#
+# GNU Radio is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with GNU Radio; see the file COPYING. If not, write to
+# the Free Software Foundation, Inc., 51 Franklin Street,
+# Boston, MA 02110-1301, USA.
+#
+
+from gnuradio import gr
+import numpy as np
+from scipy.stats import norm, laplace, rayleigh
+from matplotlib import pyplot as plt
+
+# NOTE: scipy and matplotlib are optional packages and not included in the default gnuradio dependencies
+
+#*** SETUP ***#
+
+# Number of realisations per histogram
+num_tests = 1000000
+
+# Set number of bins in histograms
+uniform_num_bins = 31
+gauss_num_bins = 31
+rayleigh_num_bins = 31
+laplace_num_bins = 31
+
+rndm = gr.random() # instance of gnuradio random class (gr::random)
+
+print 'All histograms contain',num_tests,'realisations.'
+
+#*** GENERATE DATA ***#
+
+uniform_values = np.zeros(num_tests)
+gauss_values = np.zeros(num_tests)
+rayleigh_values = np.zeros(num_tests)
+laplace_values = np.zeros(num_tests)
+
+for k in range(num_tests):
+ uniform_values[k] = rndm.ran1()
+ gauss_values[k] = rndm.gasdev()
+ rayleigh_values[k] = rndm.rayleigh()
+ laplace_values[k] = rndm.laplacian()
+
+#*** HISTOGRAM DATA AND CALCULATE EXPECTED COUNTS ***#
+
+uniform_bins = np.linspace(0,1,uniform_num_bins)
+gauss_bins = np.linspace(-8,8,gauss_num_bins)
+laplace_bins = np.linspace(-8,8,laplace_num_bins)
+rayleigh_bins = np.linspace(0,10,rayleigh_num_bins)
+
+uniform_hist = np.histogram(uniform_values,uniform_bins)
+gauss_hist = np.histogram(gauss_values,gauss_bins)
+rayleigh_hist = np.histogram(rayleigh_values,rayleigh_bins)
+laplace_hist = np.histogram(laplace_values,laplace_bins)
+
+uniform_expected = np.zeros(uniform_num_bins-1)
+gauss_expected = np.zeros(gauss_num_bins-1)
+rayleigh_expected = np.zeros(rayleigh_num_bins-1)
+laplace_expected = np.zeros(laplace_num_bins-1)
+
+for k in range(len(uniform_hist[0])):
+ uniform_expected[k] = num_tests/float(uniform_num_bins-1)
+
+for k in range(len(gauss_hist[0])):
+ gauss_expected[k] = float(norm.cdf(gauss_hist[1][k+1])-norm.cdf(gauss_hist[1][k]))*num_tests
+
+for k in range(len(rayleigh_hist[0])):
+ rayleigh_expected[k] = float(rayleigh.cdf(rayleigh_hist[1][k+1])-rayleigh.cdf(rayleigh_hist[1][k]))*num_tests
+
+for k in range(len(laplace_hist[0])):
+ laplace_expected[k] = float(laplace.cdf(laplace_hist[1][k+1])-laplace.cdf(laplace_hist[1][k]))*num_tests
+
+#*** PLOT HISTOGRAMS AND EXPECTATIONS TAKEN FROM SCIPY ***#
+
+uniform_bins_center = uniform_bins[0:-1]+(uniform_bins[1]-uniform_bins[0])/2.0
+gauss_bins_center = gauss_bins[0:-1]+(gauss_bins[1]-gauss_bins[0])/2.0
+rayleigh_bins_center = rayleigh_bins[0:-1]+(rayleigh_bins[1]-rayleigh_bins[0])/2.0
+laplace_bins_center = laplace_bins[0:-1]+(laplace_bins[1]-laplace_bins[0])/2.0
+
+plt.figure(1)
+
+plt.subplot(2,1,1)
+plt.plot(uniform_bins_center,uniform_hist[0],'s--',uniform_bins_center,uniform_expected,'o:')
+plt.xlabel('Bins'), plt.ylabel('Count'), plt.title('Uniform: Distribution')
+plt.legend(['histogram gr::random','calculation scipy'],loc=1)
+
+plt.subplot(2,1,2)
+plt.plot(uniform_bins_center,uniform_hist[0]/uniform_expected,'rs--')
+plt.xlabel('Bins'), plt.ylabel('Relative deviation'), plt.title('Uniform: Relative deviation to scipy')
+
+plt.figure(2)
+
+plt.subplot(2,1,1)
+plt.plot(gauss_bins_center,gauss_hist[0],'s--',gauss_bins_center,gauss_expected,'o:')
+plt.xlabel('Bins'), plt.ylabel('Count'), plt.title('Gauss: Distribution')
+plt.legend(['histogram gr::random','calculation scipy'],loc=1)
+
+plt.subplot(2,1,2)
+plt.plot(gauss_bins_center,gauss_hist[0]/gauss_expected,'rs--')
+plt.xlabel('Bins'), plt.ylabel('Relative deviation'), plt.title('Gauss: Relative deviation to scipy')
+
+plt.figure(3)
+
+plt.subplot(2,1,1)
+plt.plot(rayleigh_bins_center,rayleigh_hist[0],'s--',rayleigh_bins_center,rayleigh_expected,'o:')
+plt.xlabel('Bins'), plt.ylabel('Count'), plt.title('Rayleigh: Distribution')
+plt.legend(['histogram gr::random','calculation scipy'],loc=1)
+
+
+plt.subplot(2,1,2)
+plt.plot(rayleigh_bins_center,rayleigh_hist[0]/rayleigh_expected,'rs--')
+plt.xlabel('Bins'), plt.ylabel('Relative deviation'), plt.title('Rayleigh: Relative deviation to scipy')
+
+plt.figure(4)
+
+plt.subplot(2,1,1)
+plt.plot(laplace_bins_center,laplace_hist[0],'s--',laplace_bins_center,laplace_expected,'o:')
+plt.xlabel('Bins'), plt.ylabel('Count'), plt.title('Laplace: Distribution')
+plt.legend(['histogram gr::random','calculation scipy'],loc=1)
+
+plt.subplot(2,1,2)
+plt.plot(laplace_bins_center,laplace_hist[0]/laplace_expected,'rs--')
+plt.xlabel('Bins'), plt.ylabel('Relative deviation'), plt.title('Laplace: Relative deviation to scipy')
+
+plt.show()
diff --git a/gnuradio-runtime/python/gnuradio/gr/qa_random.py b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
index c0d9a7f34c..83fee56181 100644
--- a/gnuradio-runtime/python/gnuradio/gr/qa_random.py
+++ b/gnuradio-runtime/python/gnuradio/gr/qa_random.py
@@ -1,6 +1,6 @@
#!/usr/bin/env python
#
-# Copyright 2006,2007,2010 Free Software Foundation, Inc.
+# Copyright 2006,2007,2010,2015 Free Software Foundation, Inc.
#
# This file is part of GNU Radio
#
@@ -22,133 +22,42 @@
from gnuradio import gr, gr_unittest
import numpy as np
-from scipy.stats import norm, laplace, rayleigh
-#from time import sleep
class test_random(gr_unittest.TestCase):
- num_tests = 10000
+ # NOTE: For tests on the output distribution of the random numbers, see gnuradio-runtime/apps/evaluation_random_numbers.py.
- # Disclaimer
- def test_0(self):
- print 'NOTE: Following tests are not statistically significant!'
- print 'Realisations per test:',self.num_tests
- self.assertEqual(1,1)
-
- # Check for range [0,1) of uniform distributed random numbers and print minimal and maximal value
+ # Check for range [0,1) of uniform distributed random numbers
def test_1(self):
- print '# TEST 1'
- print 'Uniform distributed numbers: Range'
- values = np.zeros(self.num_tests)
+ num_tests = 10000
+ values = np.zeros(num_tests)
rndm = gr.random()
- for k in range(self.num_tests):
+ for k in range(num_tests):
values[k] = rndm.ran1()
for value in values:
self.assertLess(value, 1)
self.assertGreaterEqual(value, 0)
- print 'Uniform random numbers (num/min/max):', self.num_tests, min(values), max(values)
- # Check uniformly distributed random numbers on uniformity (without assert, only printing)
+ # Check reseed method (init with time and seed as fix number)
def test_2(self):
- print '# TEST 2'
- print 'Uniform random numbers: Distribution'
- num_bins = 11
- values = np.zeros(self.num_tests)
- rndm = gr.random()
- for k in range(self.num_tests):
- values[k] = rndm.ran1()
- bins = np.linspace(0,1,num_bins) # These are the bin edges!
- hist = np.histogram(values,bins)
- print 'Lower edge bin / upper edge bin / count / expected'
- for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(self.num_tests)/(num_bins-1)
-
- # Check distribution of normally (gaussian, mean=0, variance=1) distributed random numbers (no assert)
- def test_3(self):
- print '# TEST 3'
- print 'Normal random numbers: Distribution'
- num_bins = 11
- hist_range = [-5,5]
- values = np.zeros(self.num_tests)
- rndm = gr.random()
- for k in range(self.num_tests):
- values[k] = rndm.gasdev()
- bins = np.linspace(hist_range[0],hist_range[1],num_bins)
- hist = np.histogram(values,bins)
- print 'Lower edge bin / upper edge bin / count / expected'
- for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(norm.cdf(hist[1][k+1])-norm.cdf(hist[1][k]))*self.num_tests
-
- # Check distribution of laplacian (mean=0, variance=1) distributed random numbers (no assert)
- def test_4(self):
- print '# TEST 4'
- print 'Laplacian random numbers: Distribution'
- num_bins = 11
- hist_range = [-5,5]
- values = np.zeros(self.num_tests)
- rndm = gr.random()
- for k in range(self.num_tests):
- values[k] = rndm.laplacian()
- bins = np.linspace(hist_range[0],hist_range[1],num_bins)
- hist = np.histogram(values,bins)
- print 'Lower edge bin / upper edge bin / count / expected'
- for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(laplace.cdf(hist[1][k+1])-laplace.cdf(hist[1][k]))*self.num_tests
-
- # Check distribution of laplacian (mean=0, variance=1) distributed random numbers (no assert)
- def test_5(self):
- print '# TEST 5'
- print 'Rayleigh random numbers: Distribution'
- num_bins = 11
- hist_range = [0,10]
- values = np.zeros(self.num_tests)
- rndm = gr.random()
- for k in range(self.num_tests):
- values[k] = rndm.rayleigh()
- bins = np.linspace(hist_range[0],hist_range[1],num_bins)
- hist = np.histogram(values,bins)
- print 'Lower edge bin / upper edge bin / count / expected'
- for k in range(len(hist[0])):
- print hist[1][k], hist[1][k+1], hist[0][k], float(rayleigh.cdf(hist[1][k+1])-rayleigh.cdf(hist[1][k]))*self.num_tests
-
- # Check seeds (init with time and seed as fix number)
- def test_6(self):
- print '# TEST 6'
num = 5
- print 'Some random numbers in [0,1), should change every run:'
- rndm0 = gr.random(0); # init with time
- # NOTE: the sleep increases the executiont time massively, remove assert for convenience
- #sleep(1)
- #rndm1 = gr.random(0); # init with fix seed
- for k in range(num):
- x = rndm0.ran1();
- print x,
- # y = rndm1.ran1();
- # print x, '!=', y
- # self.assertNotEqual(x,y)
- print ' '
-
- print 'Some random numbers in [0,1) (seed two instances), should be the same every run:'
rndm0 = gr.random(42); # init with time
rndm1 = gr.random(42); # init with fix seed
for k in range(num):
x = rndm0.ran1();
y = rndm1.ran1();
- print x, '=', y
self.assertEqual(x,y)
- print 'Some random numbers in [0,1) (reseed one instance), should be the same every run:'
x = np.zeros(num)
y = np.zeros(num)
- rndm0 = gr.random(42); # init with time
+ rndm0 = gr.random(42); # init with fix seed 1
for k in range(num):
x[k] = rndm0.ran1();
- rndm1.reseed(43); # init with fix seed
+ rndm1.reseed(43); # init with fix seed 2
for k in range(num):
y[k] = rndm0.ran1();
for k in range(num):
- print x[k], '!=', y[k]
self.assertNotEqual(x[k],y[k])
if __name__ == '__main__':
--
cgit v1.2.3