Imported Upstream version 3.2.2
[debian/gnuradio] / gnuradio-core / src / lib / filter / qa_gr_fir_ccf.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 3, 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 gr_complex      i_type;
35 typedef gr_complex      o_type;
36 typedef float           tap_type;
37 typedef gr_complex      acc_type;
38
39
40 #include <qa_gr_fir_ccf.h>
41 #include <gr_fir_ccf.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 #include <malloc16.h>
50 #include <string.h>
51
52 using std::vector;
53
54 #define ERR_DELTA       (1e-5)
55
56 #define NELEM(x) (sizeof (x) / sizeof (x[0]))
57
58 //
59 // typedef for something logically "pointer to constructor".
60 // there may be a better way, please let me know...
61 //
62 typedef gr_fir_ccf* (*fir_maker_t)(const std::vector<tap_type> &taps);
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_floats (float *buf, unsigned n)
72 {
73   for (unsigned i = 0; i < n; i++)
74     buf[i] = (float) 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   // Our SIMD ccc kernel requires that the complex input be 64-bit (8-byte) aligned.
111   //i_type      input[INPUT_LEN];
112   i_type       *input = (i_type *)malloc16Align(INPUT_LEN * sizeof(i_type));
113   o_type        expected_output[OUTPUT_LEN];
114   o_type        actual_output[OUTPUT_LEN];
115   tap_type      taps[MAX_TAPS];
116
117
118   srandom (0);  // we want reproducibility
119
120   for (int n = 0; n <= MAX_TAPS; n++){
121     for (int ol = 0; ol <= OUTPUT_LEN; ol++){
122
123       // cerr << "@@@ n:ol " << n << ":" << ol << endl;
124
125       // build random test case
126       random_complex (input, INPUT_LEN);
127       random_floats (taps, MAX_TAPS);
128
129       // compute expected output values
130       for (int o = 0; o < ol; o++){
131         expected_output[o] = ref_dotprod (&input[o], taps, n);
132       }
133
134       // build filter
135
136       vector<tap_type> f1_taps (&taps[0], &taps[n]);
137       gr_fir_ccf *f1 = maker (f1_taps);
138
139       // zero the output, then do the filtering
140       memset (actual_output, 0, sizeof (actual_output));
141       f1->filterN (actual_output, input, ol);
142
143       // check results
144       //
145       // we use a sloppy error margin because on the x86 architecture,
146       // our reference implementation is using 80 bit floating point
147       // arithmetic, while the SSE version is using 32 bit float point
148       // arithmetic.
149       
150       for (int o = 0; o < ol; o++){
151         CPPUNIT_ASSERT_COMPLEXES_EQUAL(expected_output[o], actual_output[o],
152                                        abs (expected_output[o]) * ERR_DELTA);
153       }
154
155       delete f1;
156     }
157   }
158   free16Align(input);
159 }
160
161
162 static void
163 for_each (void (*f)(fir_maker_t))
164 {
165   std::vector<gr_fir_ccf_info>          info;
166   gr_fir_util::get_gr_fir_ccf_info (&info);     // get all known ccf implementations 
167
168   for (std::vector<gr_fir_ccf_info>::iterator p = info.begin ();
169        p != info.end ();
170        ++p){
171
172     std::cerr << " [" << p->name << "]";
173     f (p->create);
174   }
175
176   std::cerr << std::endl;
177 }
178
179
180 void
181 qa_gr_fir_ccf::t1 ()
182 {
183   for_each (test_random_io);
184 }