+/* -*- c++ -*- */
+/*
+ * Copyright 2010 Free Software Foundation, Inc.
+ *
+ * This file is part of GNU Radio
+ *
+ * GNU Radio is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 3, or (at your option)
+ * any later version.
+ *
+ * GNU Radio is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with GNU Radio; see the file COPYING. If not, write to
+ * the Free Software Foundation, Inc., 51 Franklin Street,
+ * Boston, MA 02110-1301, USA.
+ */
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <gr_types.h>
+#include <qa_gri_fir_filter_with_buffer_scc.h>
+#include <gri_fir_filter_with_buffer_scc.h>
+#include <string.h>
+#include <iostream>
+#include <cmath>
+#include <cppunit/TestAssert.h>
+#include <random.h>
+#include <malloc16.h>
+#include <string.h>
+
+typedef short i_type;
+typedef gr_complex o_type;
+typedef gr_complex tap_type;
+typedef gr_complex acc_type;
+
+using std::vector;
+
+#define ERR_DELTA (1e-5)
+
+#define NELEM(x) (sizeof (x) / sizeof (x[0]))
+
+static float
+uniform ()
+{
+ return 2.0 * ((float) random () / RANDOM_MAX - 0.5); // uniformly (-1, 1)
+}
+
+static void
+random_shorts (short *buf, unsigned n)
+{
+ for (unsigned i = 0; i < n; i++)
+ buf[i] = (short) rint (uniform () * 16384);
+}
+
+static void
+random_complex (gr_complex *buf, unsigned n)
+{
+ for (unsigned i = 0; i < n; i++){
+ float re = rint (uniform () * 32767);
+ float im = rint (uniform () * 32767);
+ buf[i] = gr_complex (re, im);
+ }
+}
+
+static o_type
+ref_dotprod (const i_type input[], const tap_type taps[], int ntaps)
+{
+ acc_type sum = 0;
+ for (int i = 0; i < ntaps; i++) {
+ sum += (float)input[i] * taps[i];
+ }
+
+ return sum;
+}
+
+void
+qa_gri_fir_filter_with_buffer_scc::t1 ()
+{
+ test_decimate(1);
+}
+
+void
+qa_gri_fir_filter_with_buffer_scc::t2 ()
+{
+ test_decimate(2);
+}
+
+void
+qa_gri_fir_filter_with_buffer_scc::t3 ()
+{
+ test_decimate(5);
+}
+
+//
+// Test for ntaps in [0,9], and input lengths in [0,17].
+// This ensures that we are building the shifted taps correctly,
+// and exercises all corner cases on input alignment and length.
+//
+void
+qa_gri_fir_filter_with_buffer_scc::test_decimate (unsigned int decimate)
+{
+ const int MAX_TAPS = 9;
+ const int OUTPUT_LEN = 17;
+ const int INPUT_LEN = MAX_TAPS + OUTPUT_LEN;
+
+ // Mem aligned buffer not really necessary, but why not?
+ i_type *input = (i_type *)malloc16Align(INPUT_LEN * sizeof(i_type));
+ i_type *dline = (i_type*)malloc16Align(INPUT_LEN * sizeof(i_type));
+ o_type expected_output[OUTPUT_LEN];
+ o_type actual_output[OUTPUT_LEN];
+ tap_type taps[MAX_TAPS];
+
+ srandom (0); // we want reproducibility
+ memset(dline, 0, INPUT_LEN*sizeof(i_type));
+
+ for (int n = 0; n <= MAX_TAPS; n++){
+ for (int ol = 0; ol <= OUTPUT_LEN; ol++){
+
+ // cerr << "@@@ n:ol " << n << ":" << ol << endl;
+
+ // build random test case
+ random_shorts (input, INPUT_LEN);
+ random_complex (taps, MAX_TAPS);
+
+ // compute expected output values
+ memset(dline, 0, INPUT_LEN*sizeof(i_type));
+ for (int o = 0; o < (int)(ol/decimate); o++){
+ // use an actual delay line for this test
+ for(int dd = 0; dd < (int)decimate; dd++) {
+ for(int oo = INPUT_LEN-1; oo > 0; oo--)
+ dline[oo] = dline[oo-1];
+ dline[0] = input[decimate*o+dd];
+ }
+ expected_output[o] = ref_dotprod (dline, taps, n);
+ }
+
+ // build filter
+ vector<tap_type> f1_taps(&taps[0], &taps[n]);
+ gri_fir_filter_with_buffer_scc *f1 = new gri_fir_filter_with_buffer_scc(f1_taps);
+
+ // zero the output, then do the filtering
+ memset (actual_output, 0, sizeof (actual_output));
+ f1->filterNdec (actual_output, input, ol/decimate, decimate);
+
+ // check results
+ //
+ // we use a sloppy error margin because on the x86 architecture,
+ // our reference implementation is using 80 bit floating point
+ // arithmetic, while the SSE version is using 32 bit float point
+ // arithmetic.
+
+ for (int o = 0; o < (int)(ol/decimate); o++){
+ CPPUNIT_ASSERT_COMPLEXES_EQUAL(expected_output[o], actual_output[o],
+ abs (expected_output[o]) * ERR_DELTA);
+ }
+ delete f1;
+ }
+ }
+ free16Align(input);
+}