From b8545bfea262e3fcd4e3796c8d61fdf3c6f9cb70 Mon Sep 17 00:00:00 2001
From: Marti Bolivar <mbolivar@leaflabs.com>
Date: Wed, 31 Aug 2011 11:26:53 -0400
Subject: test-session.cpp: Better analog noise readings.

Update measure_adc_noise() to actually use the Welford online
algorithm, instead of accumulating data in an array on the stack.
This allows us to increase the number of samples (to 1000).

Revised algorithm tested on host PC and compared (in Python) against
numpy with a list of 100 values in [0, 1) drawn using random.random().

Results (Python):

    >>> r = [random.random() for i in xrange(100)]
    >>> numpy.mean(r)
    0.50073064742634854
    >>> numpy.var(r)
    0.083726090293309297

Results (C++, x86 host PC):

     n: 100	mean: 0.500731	variance: 0.084572

So this algorithm for variance has some inaccuracies, but it appears
to be good to a couple of significant figures.
---
 examples/test-session.cpp | 17 +++++++++--------
 1 file changed, 9 insertions(+), 8 deletions(-)

(limited to 'examples')

diff --git a/examples/test-session.cpp b/examples/test-session.cpp
index 1a53c66..97d5754 100644
--- a/examples/test-session.cpp
+++ b/examples/test-session.cpp
@@ -759,7 +759,8 @@ void cmd_board_info(void) {     // TODO print more information
 // -- Helper functions --------------------------------------------------------
 
 void measure_adc_noise(uint8 pin) {
-    uint16 data[100];
+    const int N = 1000;
+    uint16 x;
     float mean = 0;
     float delta = 0;
     float M2 = 0;
@@ -767,21 +768,21 @@ void measure_adc_noise(uint8 pin) {
 
     // Variance algorithm from Welford, via Knuth, by way of Wikipedia:
     // http://en.wikipedia.org/wiki/Algorithms_for_calculating_variance#On-line_algorithm
-    for (int i = 0; i < 100; i++) {
-        data[i] = analogRead(pin);
-        delta = data[i] - mean;
-        mean = mean + delta / (i + 1);
-        M2 = M2 + delta * (data[i] - mean);
+    for (int i = 1; i <= N; i++) {
+        x = analogRead(pin);
+        delta = x - mean;
+        mean += delta / i;
+        M2 = M2 + delta * (x - mean);
     }
 
     SerialUSB.print("header: D");
     SerialUSB.print(pin, DEC);
     SerialUSB.print("\tn: ");
-    SerialUSB.print(100, DEC);
+    SerialUSB.print(N, DEC);
     SerialUSB.print("\tmean: ");
     SerialUSB.print(mean);
     SerialUSB.print("\tvariance: ");
-    SerialUSB.println(M2 / 99.0);
+    SerialUSB.println(M2 / (float)(N-1));
     pinMode(pin, OUTPUT);
 }
 
-- 
cgit v1.2.3