3 * Copyright 2002 Free Software Foundation, Inc.
5 * This file is part of GNU Radio
7 * GNU Radio is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 3, or (at your option)
12 * GNU Radio is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
17 * You should have received a copy of the GNU General Public License
18 * along with GNU Radio; see the file COPYING. If not, write to
19 * the Free Software Foundation, Inc., 51 Franklin Street,
20 * Boston, MA 02110-1301, USA.
27 #include <qa_gr_fir_fff.h>
28 #include <gr_fir_fff.h>
29 #include <gr_fir_util.h>
33 #include <cppunit/TestAssert.h>
41 typedef float tap_type;
42 typedef float acc_type;
44 #define ERR_DELTA (1e-6)
46 #define NELEM(x) (sizeof (x) / sizeof (x[0]))
50 // typedef for something logically "pointer to constructor".
51 // there may be a better way, please let me know...
53 typedef gr_fir_fff* (*fir_maker_t)(const std::vector<tap_type> &taps);
56 // ----------------------------------------------------------------
58 const static i_type input_1[] = {
59 234, -4, 23, -56, 45, 98, -23, -7
62 const static tap_type taps_1a[] = {
66 const static o_type expected_1a[] = {
67 -702, 12, -69, 168, -135, -294, 69, 21
70 const static tap_type taps_1b[] = {
74 const static o_type expected_1b[] = {
75 1186, -112, 339, -460, -167, 582, -87
78 // ----------------------------------------------------------------
81 test_known_io (fir_maker_t maker)
83 vector<tap_type> t1a (&taps_1a[0], &taps_1a[NELEM (taps_1a)]);
84 vector<tap_type> t1b (&taps_1b[0], &taps_1b[NELEM (taps_1b)]);
86 gr_fir_fff *f1 = maker (t1a); // create filter
87 CPPUNIT_ASSERT_EQUAL ((unsigned) 1, f1->ntaps ()); // check ntaps
89 // check filter output
90 int n = NELEM (input_1) - f1->ntaps () + 1;
91 for (int i = 0; i < n; i++)
92 CPPUNIT_ASSERT_DOUBLES_EQUAL (expected_1a[i], f1->filter (&input_1[i]), ERR_DELTA);
94 f1->set_taps (t1b); // set new taps
95 CPPUNIT_ASSERT_EQUAL ((unsigned) 2, f1->ntaps ()); // check ntaps
97 // check filter output
98 n = NELEM (input_1) - f1->ntaps () + 1;
99 for (int i = 0; i < n; i++)
100 CPPUNIT_ASSERT_DOUBLES_EQUAL (expected_1b[i], f1->filter (&input_1[i]), ERR_DELTA);
102 // test filterN interface
104 o_type output[NELEM (expected_1b)];
105 memset (output, 0, sizeof (output));
107 f1->filterN (output, input_1, n);
108 for (int i = 0; i < n; i++)
109 CPPUNIT_ASSERT_DOUBLES_EQUAL (expected_1b[i], output[i], ERR_DELTA);
115 // Test for ntaps in [0,9], and input lengths in [0,17].
116 // This ensures that we are building the shifted taps correctly,
117 // and exercises all corner cases on input alignment and length.
123 return 2.0 * ((float) random () / RANDOM_MAX - 0.5); // uniformly (-1, 1)
127 random_floats (float *buf, unsigned n)
129 for (unsigned i = 0; i < n; i++)
130 buf[i] = rint (uniform () * 32768);
134 ref_dotprod (const i_type input[], const tap_type taps[], int ntaps)
137 for (int i = 0; i < ntaps; i++)
138 sum += input[i] * taps[ntaps - i - 1];
144 test_random_io (fir_maker_t maker)
146 const int MAX_TAPS = 32;
147 const int OUTPUT_LEN = 17;
148 const int INPUT_LEN = MAX_TAPS + OUTPUT_LEN;
150 i_type input[INPUT_LEN];
151 o_type expected_output[OUTPUT_LEN];
152 o_type actual_output[OUTPUT_LEN];
153 tap_type taps[MAX_TAPS];
156 srandom (0); // we want reproducibility
158 for (int n = 0; n <= MAX_TAPS; n++){
159 for (int ol = 0; ol <= OUTPUT_LEN; ol++){
161 // cerr << "@@@ n:ol " << n << ":" << ol << endl;
163 // build random test case
164 random_floats (input, INPUT_LEN);
165 random_floats (taps, MAX_TAPS);
167 // compute expected output values
168 for (int o = 0; o < ol; o++){
169 expected_output[o] = ref_dotprod (&input[o], taps, n);
174 vector<tap_type> f1_taps (&taps[0], &taps[n]);
175 gr_fir_fff *f1 = maker (f1_taps);
177 // zero the output, then do the filtering
178 memset (actual_output, 0, sizeof (actual_output));
179 f1->filterN (actual_output, input, ol);
183 // we use a sloppy error margin because on the x86 architecture,
184 // our reference implementation is using 80 bit floating point
185 // arithmetic, while the SSE version is using 32 bit float point
188 for (int o = 0; o < ol; o++){
189 CPPUNIT_ASSERT_DOUBLES_EQUAL (expected_output[o], actual_output[o],
190 fabs (expected_output[o]) * 9e-3);
200 for_each (void (*f)(fir_maker_t))
202 std::vector<gr_fir_fff_info> info;
203 gr_fir_util::get_gr_fir_fff_info (&info); // get all known fff implementations
205 for (std::vector<gr_fir_fff_info>::iterator p = info.begin ();
209 std::cerr << " [" << p->name << "]";
213 std::cerr << std::endl;
219 for_each (test_known_io);
225 for_each (test_random_io);