Imported Upstream version 3.0
[debian/gnuradio] / gnuradio-core / src / lib / filter / qa_gr_fir_fcc.cc
1 /* -*- c++ -*- */
2 /*
3  * Copyright 2002 Free Software Foundation, Inc.
4  * 
5  * This file is part of GNU Radio
6  * 
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 2, or (at your option)
10  * any later version.
11  * 
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.
16  * 
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.
21  */
22
23 /*
24  * FIXME.  This code is virtually identical to qa_gr_fir_?CC.cc
25  *   Kludge up some kind of macro to handle the minor differences.
26  */
27
28 #ifdef HAVE_CONFIG_H
29 #include <config.h>
30 #endif
31
32 #include <gr_types.h>
33
34 typedef float           i_type;
35 typedef gr_complex      o_type;
36 typedef gr_complex      tap_type;
37 typedef gr_complex      acc_type;
38
39
40 #include <qa_gr_fir_fcc.h>
41 #include <gr_fir_fcc.h>
42 #include <gr_fir_util.h>
43 #include <string.h>
44 #include <iostream>
45 #include <cmath>
46 #include <gr_types.h>
47 #include <cppunit/TestAssert.h>
48 #include <random.h>
49
50 using std::vector;
51
52 #define ERR_DELTA       (1e-5)
53
54 #define NELEM(x) (sizeof (x) / sizeof (x[0]))
55
56
57 //
58 // typedef for something logically "pointer to constructor".
59 // there may be a better way, please let me know...
60 //
61 typedef gr_fir_fcc* (*fir_maker_t)(const std::vector<tap_type> &taps);
62
63
64 static float
65 uniform ()
66 {
67   return 2.0 * ((float) random () / RANDOM_MAX - 0.5);  // uniformly (-1, 1)
68 }
69
70 static void
71 random_input (i_type *buf, unsigned n)
72 {
73   for (unsigned i = 0; i < n; i++)
74     buf[i] = (i_type) rint (uniform () * 32767);
75 }
76
77 static void
78 random_complex (gr_complex *buf, unsigned n)
79 {
80   for (unsigned i = 0; i < n; i++){
81     float re = rint (uniform () * 32767);
82     float im = rint (uniform () * 32767);
83     buf[i] = gr_complex (re, im);
84   }
85 }
86
87 static o_type
88 ref_dotprod (const i_type input[], const tap_type taps[], int ntaps)
89 {
90   acc_type      sum = 0;
91   for (int i = 0; i < ntaps; i++)
92     sum += input[i] * taps[ntaps - i - 1];
93
94   return sum;
95 }
96
97 //
98 // Test for ntaps in [0,9], and input lengths in [0,17].
99 // This ensures that we are building the shifted taps correctly,
100 // and exercises all corner cases on input alignment and length.
101 //
102
103 static void
104 test_random_io (fir_maker_t maker)  
105 {
106   const int     MAX_TAPS        = 9;
107   const int     OUTPUT_LEN      = 17;
108   const int     INPUT_LEN       = MAX_TAPS + OUTPUT_LEN;
109
110   i_type        input[INPUT_LEN];
111   o_type        expected_output[OUTPUT_LEN];
112   o_type        actual_output[OUTPUT_LEN];
113   tap_type      taps[MAX_TAPS];
114
115
116   srandom (0);  // we want reproducibility
117
118   for (int n = 0; n <= MAX_TAPS; n++){
119     for (int ol = 0; ol <= OUTPUT_LEN; ol++){
120
121       // cerr << "@@@ n:ol " << n << ":" << ol << endl;
122
123       // build random test case
124       random_input (input, INPUT_LEN);
125       random_complex (taps, MAX_TAPS);
126
127       // compute expected output values
128       for (int o = 0; o < ol; o++){
129         expected_output[o] = ref_dotprod (&input[o], taps, n);
130       }
131
132       // build filter
133
134       vector<tap_type> f1_taps (&taps[0], &taps[n]);
135       gr_fir_fcc *f1 = maker (f1_taps);
136
137       // zero the output, then do the filtering
138       memset (actual_output, 0, sizeof (actual_output));
139       f1->filterN (actual_output, input, ol);
140
141       // check results
142       //
143       // we use a sloppy error margin because on the x86 architecture,
144       // our reference implementation is using 80 bit floating point
145       // arithmetic, while the SSE version is using 32 bit float point
146       // arithmetic.
147       
148       for (int o = 0; o < ol; o++){
149         ASSERT_COMPLEXES_EQUAL (expected_output[o],
150                                 actual_output[o],
151                                 abs (expected_output[o]) * ERR_DELTA);
152       }
153
154       delete f1;
155     }
156   }
157 }
158
159 static void
160 for_each (void (*f)(fir_maker_t))
161 {
162   std::vector<gr_fir_fcc_info>          info;
163   gr_fir_util::get_gr_fir_fcc_info (&info);     // get all known fcc implementations 
164
165   for (std::vector<gr_fir_fcc_info>::iterator p = info.begin ();
166        p != info.end ();
167        ++p){
168
169     std::cerr << " [" << p->name << "]";
170     f (p->create);
171   }
172
173   std::cerr << std::endl;
174 }
175
176 void
177 qa_gr_fir_fcc::t1 ()
178 {
179   for_each (test_random_io);
180 }