comparison spandsp-0.0.3/spandsp-0.0.3/tests/awgn_tests.c @ 5:f762bf195c4b

import spandsp-0.0.3
author Peter Meerwald <pmeerw@cosy.sbg.ac.at>
date Fri, 25 Jun 2010 16:00:21 +0200
parents
children
comparison
equal deleted inserted replaced
4:26cd8f1ef0b1 5:f762bf195c4b
1 /*
2 * SpanDSP - a series of DSP components for telephony
3 *
4 * awgn_tests.c
5 *
6 * Written by Steve Underwood <steveu@coppice.org>
7 *
8 * Copyright (C) 2001 Steve Underwood
9 *
10 * All rights reserved.
11 *
12 * This program is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU General Public License version 2, as
14 * published by the Free Software Foundation.
15 *
16 * This program is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 * GNU General Public License for more details.
20 *
21 * You should have received a copy of the GNU General Public License
22 * along with this program; if not, write to the Free Software
23 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
24 *
25 * $Id: awgn_tests.c,v 1.12 2006/11/19 14:07:26 steveu Exp $
26 */
27
28 /*! \page awgn_tests_page AWGN tests
29 \section awgn_tests_page_sec_1 What does it do?
30 */
31
32 #ifdef HAVE_CONFIG_H
33 #include "config.h"
34 #endif
35
36 #include <stdio.h>
37 #include <inttypes.h>
38 #include <stdlib.h>
39 #include <string.h>
40 #if defined(HAVE_TGMATH_H)
41 #include <tgmath.h>
42 #endif
43 #if defined(HAVE_MATH_H)
44 #include <math.h>
45 #endif
46 #include <tiffio.h>
47
48 #include "spandsp.h"
49
50 #if !defined(M_PI)
51 # define M_PI 3.14159265358979323846 /* pi */
52 #endif
53
54 #define OUT_FILE_NAME "awgn.wav"
55
56 /* Some simple sanity tests for the Gaussian noise generation routines */
57
58 int main (int argc, char *argv[])
59 {
60 int i;
61 int j;
62 int clip_high;
63 int clip_low;
64 int total_samples;
65 int idum = 1234567;
66 int16_t value;
67 double total;
68 double x;
69 double p;
70 double o;
71 double error;
72 int bins[65536];
73 awgn_state_t noise_source;
74
75 /* Generate noise at several RMS levels between -50dBm and 0dBm. Noise is
76 generated for a large number of samples (1,000,000), and the RMS value
77 of the noise is calculated along the way. If the resulting level is
78 close to the requested RMS level, at least the scaling of the noise
79 should be Ok. At high level some clipping may distort the result a
80 little. */
81 for (j = -50; j <= 0; j += 5)
82 {
83 clip_high = 0;
84 clip_low = 0;
85 total = 0.0;
86 awgn_init_dbm0(&noise_source, idum, (float) j);
87 total_samples = 1000000;
88 for (i = 0; i < total_samples; i++)
89 {
90 value = awgn(&noise_source);
91 if (value == 32767)
92 clip_high++;
93 else if (value == -32768)
94 clip_low++;
95 total += ((double) value)*((double) value);
96 }
97 error = 100.0*(1.0 - sqrt(total/total_samples)/noise_source.rms);
98 printf("RMS = %.3f (expected %d) %.2f%% error [clipped samples %d+%d]\n",
99 log10(sqrt(total/total_samples)/32768.0)*20.0 + DBM0_MAX_POWER,
100 j,
101 error,
102 clip_low,
103 clip_high);
104 /* We don't check the result at 0dBm0, as there will definitely be a lot of error due to clipping */
105 if (j < 0 && fabs(error) > 0.2)
106 {
107 printf("Test failed.\n");
108 exit(2);
109 }
110 }
111 /* Now look at the statistical spread of the results, by collecting data in
112 bins from a large number of samples. Use a fairly high noise level, but
113 low enough to avoid significant clipping. Use the Gaussian model to
114 predict the real probability, and present the results for graphing. */
115 memset(bins, 0, sizeof(bins));
116 clip_high = 0;
117 clip_low = 0;
118 awgn_init_dbm0(&noise_source, idum, -15);
119 total_samples = 10000000;
120 for (i = 0; i < total_samples; i++)
121 {
122 value = awgn(&noise_source);
123 if (value == 32767)
124 clip_high++;
125 else if (value == -32768)
126 clip_low++;
127 bins[value + 32768]++;
128 }
129 o = noise_source.rms;
130 for (i = 0; i < 65536 - 10; i++)
131 {
132 x = i - 32768;
133 /* Find the real probability for this bin */
134 p = (1.0/(o*sqrt(2.0*M_PI)))*exp(-(x*x)/(2.0*o*o));
135 /* Now do a little smoothing on the real data to get a reasonably
136 steady answer */
137 x = 0;
138 for (j = 0; j < 10; j++)
139 x += bins[i + j];
140 x /= 10.0;
141 x /= total_samples;
142 /* Now send it out for graphing. */
143 printf("%6d %.7f %.7f\n", i - 32768, x, p);
144 }
145
146 printf("Tests passed.\n");
147 return 0;
148 }
149 /*- End of function --------------------------------------------------------*/
150 /*- End of file ------------------------------------------------------------*/

Repositories maintained by Peter Meerwald, pmeerw@pmeerw.net.