git.gag.com Git - debian/gnuradio/blob - gnuradio-core/src/lib/filter/qa_complex_dotprod_x86.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 #ifdef HAVE_CONFIG_H
  24 #include <config.h>
  25 #endif
  26
  27 #include <cppunit/TestAssert.h>
  28 #include <qa_complex_dotprod_x86.h>
  29 #include <complex_dotprod_x86.h>
  30 #include <string.h>
  31 #include <iostream>
  32 #include <malloc16.h>
  33 #include <sse_debug.h>
  34 #include <cmath>
  35 #include <gr_cpu.h>
  36 #include <random.h>
  37
  38 using std::cerr;
  39
  40 /// Macro for primitive value comparisons
  41 #define assertcomplexEqual(expected0,expected1,actual,delta)    \
  42   CPPUNIT_ASSERT_DOUBLES_EQUAL (expected0, actual[0], delta);   \
  43   CPPUNIT_ASSERT_DOUBLES_EQUAL (expected1, actual[1], delta);
  44
  45
  46 #define MAX_BLKS        10
  47 #define FLOATS_PER_BLK   4
  48 #define SHORTS_PER_BLK   2
  49
  50 #define ERR_DELTA       (1e-6)
  51
  52 static float
  53 uniform ()
  54 {
  55   return 2.0 * ((float) random () / RANDOM_MAX - 0.5);  // uniformly (-1, 1)
  56 }
  57
  58 static void
  59 random_floats (float *buf, unsigned n)
  60 {
  61   for (unsigned i = 0; i < n; i++)
  62     buf[i] = rint (uniform () * 32767);
  63 }
  64
  65 static void
  66 zero_floats (float *buf, unsigned n)
  67 {
  68   for (unsigned i = 0; i < n; i++)
  69     buf[i] = 0.0;
  70 }
  71
  72 static void
  73 random_shorts (short *buf, unsigned n)
  74 {
  75   for (unsigned i = 0; i < n; i++)
  76     buf[i] = (short) rint (uniform () * 32767);
  77 }
  78
  79 static void
  80 zero_shorts (short *buf, unsigned n)
  81 {
  82   for (unsigned i = 0; i < n; i++)
  83     buf[i] = 0;
  84 }
  85
  86 void
  87 ref_complex_dotprod (const short *input,
  88                    const float *taps, unsigned n_2_complex_blocks,
  89                    float *result)
  90 {
  91   float sum0[2] = {0,0};
  92   float sum1[2] = {0,0};
  93
  94   do {
  95
  96     sum0[0] += input[0] * taps[0];
  97     sum0[1] += input[0] * taps[1];
  98     sum1[0] += input[1] * taps[2];
  99     sum1[1] += input[1] * taps[3];
 100
 101     input += 2;
 102     taps += 4;
 103
 104   } while (--n_2_complex_blocks != 0);
 105
 106
 107   result[0] = sum0[0] + sum1[0];
 108   result[1] = sum0[1] + sum1[1];
 109 }
 110
 111 void
 112 qa_complex_dotprod_x86::setUp ()
 113 {
 114   taps = (float *) calloc16Align (MAX_BLKS,
 115                                   sizeof (float) * FLOATS_PER_BLK);
 116
 117   input = (short *) calloc16Align (MAX_BLKS,
 118                                    sizeof (short) * SHORTS_PER_BLK);
 119
 120   if (taps == 0 || input == 0)
 121     abort ();
 122 }
 123
 124 void
 125 qa_complex_dotprod_x86::tearDown ()
 126 {
 127   free16Align (taps);
 128   free16Align (input);
 129   taps = 0;
 130   input = 0;
 131 }
 132
 133
 134 void
 135 qa_complex_dotprod_x86::zb ()   // "zero both"
 136 {
 137   zero_floats (taps, MAX_BLKS * FLOATS_PER_BLK);
 138   zero_shorts (input, MAX_BLKS * SHORTS_PER_BLK);
 139 }
 140
 141 //
 142 // t1
 143 //
 144
 145 void
 146 qa_complex_dotprod_x86::t1_base (complex_dotprod_t complex_dotprod)
 147 {
 148   float result[2];
 149
 150   // cerr << "Testing dump_xmm_regs\n";
 151   // dump_xmm_regs ();
 152
 153   // test basic cases, 1 block
 154
 155   zb ();
 156   complex_dotprod (input, taps, 1, result);
 157   assertcomplexEqual (0.0, 0.0, result, ERR_DELTA);
 158
 159   // vary each input
 160
 161   zb ();
 162   input[0] = 1; taps[0] = 1.0; taps[1] = -1.0;
 163   complex_dotprod (input, taps, 1, result);
 164   //cerr << result[0] << " " << result[1] << "\n";
 165   assertcomplexEqual (1.0, -1.0, result, ERR_DELTA);
 166
 167   zb ();
 168   input[1] = 2; taps[2] = 1.0; taps[3] = -1.0;
 169   complex_dotprod (input, taps, 1, result);
 170   assertcomplexEqual (2.0, -2.0, result, ERR_DELTA);
 171
 172   zb ();
 173   input[2] = 3; taps[4] = 1.0; taps[5] = -1.0;
 174   complex_dotprod (input, taps, 2, result);
 175   assertcomplexEqual (3.0, -3.0, result, ERR_DELTA);
 176
 177   zb ();
 178   input[3] = 4; taps[6] = 1.0; taps[7] = -1.0;
 179   complex_dotprod (input, taps, 2, result);
 180   assertcomplexEqual (4.0, -4.0, result, ERR_DELTA);
 181
 182   // vary each tap
 183
 184   zb ();
 185   input[0] = 1; taps[0] = 0.5; taps[1] = -0.5;
 186   complex_dotprod (input, taps, 1, result);
 187   assertcomplexEqual (0.5, -0.5, result, ERR_DELTA);
 188
 189   zb ();
 190   input[0] = 1; taps[0] = 2.0; taps[1] = -2.0;
 191   complex_dotprod (input, taps, 1, result);
 192   assertcomplexEqual (2.0, -2.0, result, ERR_DELTA);
 193
 194   zb ();
 195   input[0] = 1; taps[0] = 3.0; taps[1] = -3.0;
 196   complex_dotprod (input, taps, 1, result);
 197   assertcomplexEqual (3.0, -3.0, result, ERR_DELTA);
 198
 199   zb ();
 200   input[0] = 1; taps[0] = 4.0; taps[1] = -4.0;
 201   complex_dotprod (input, taps, 1, result);
 202   assertcomplexEqual (4.0, -4.0, result, ERR_DELTA);
 203 }
 204
 205 //
 206 // t2
 207 //
 208 void
 209 qa_complex_dotprod_x86::t2_base (complex_dotprod_t complex_dotprod)
 210 {
 211   float result[2];
 212
 213   zb ();
 214   input[0] =  1;        taps[0] =  2.0; taps[1] =  -2.0;
 215   input[1] =  3;        taps[2] =  5.0; taps[3] =  -5.0;
 216   input[2] =  7;        taps[4] = 11.0; taps[5] =  -11.0;
 217   input[3] = 13;        taps[6] = 17.0; taps[7] =  -17.0;
 218
 219   complex_dotprod (input, taps, 2, result);
 220   assertcomplexEqual (315.0, -315.0, result, ERR_DELTA);
 221
 222   input[4] = 19;        taps[8] = 23.0; taps[9] = -23.0;
 223   complex_dotprod (input, taps, 3, result);
 224   assertcomplexEqual (752.0, -752.0, result, ERR_DELTA);
 225
 226 }
 227
 228 //
 229 // t3
 230 //
 231 void
 232 qa_complex_dotprod_x86::t3_base (complex_dotprod_t complex_dotprod)
 233 {
 234   srandom (0);  // we want reproducibility
 235
 236   for (unsigned int i = 0; i < 10; i++){
 237     random_shorts (input, MAX_BLKS * SHORTS_PER_BLK);
 238     random_floats (taps, MAX_BLKS * FLOATS_PER_BLK);
 239
 240     // we use a sloppy error margin because on the x86 architecture,
 241     // our reference implementation is using 80 bit floating point
 242     // arithmetic, while the SSE version is using 32 bit float point
 243     // arithmetic.
 244
 245     float ref[2];
 246     ref_complex_dotprod (input, taps, MAX_BLKS, ref);
 247     float calc[2];
 248     complex_dotprod (input, taps, MAX_BLKS, calc);
 249     CPPUNIT_ASSERT_DOUBLES_EQUAL (ref[0],
 250                         calc[0],
 251                         fabs (ref[0]) * 1e-4);
 252     CPPUNIT_ASSERT_DOUBLES_EQUAL (ref[1],
 253                         calc[1],
 254                         fabs (ref[1]) * 1e-4);
 255   }
 256 }
 257
 258 void
 259 qa_complex_dotprod_x86::t1_3dnowext ()
 260 {
 261   if (!gr_cpu::has_3dnowext ()){
 262     cerr << "No 3DNow!Ext support; not tested\n";
 263   }
 264   else
 265     t1_base (complex_dotprod_3dnowext);
 266 }
 267
 268 void
 269 qa_complex_dotprod_x86::t2_3dnowext ()
 270 {
 271   if (!gr_cpu::has_3dnowext ()){
 272     cerr << "No 3DNow!Ext support; not tested\n";
 273   }
 274   else
 275     t2_base (complex_dotprod_3dnowext);
 276 }
 277
 278 void
 279 qa_complex_dotprod_x86::t3_3dnowext ()
 280 {
 281   if (!gr_cpu::has_3dnowext ()){
 282     cerr << "No 3DNow!Ext support; not tested\n";
 283   }
 284   else
 285     t3_base (complex_dotprod_3dnowext);
 286 }
 287
 288 void
 289 qa_complex_dotprod_x86::t1_3dnow ()
 290 {
 291   if (!gr_cpu::has_3dnow ()){
 292     cerr << "No 3DNow! support; not tested\n";
 293   }
 294   else
 295     t1_base (complex_dotprod_3dnow);
 296 }
 297
 298 void
 299 qa_complex_dotprod_x86::t2_3dnow ()
 300 {
 301   if (!gr_cpu::has_3dnow ()){
 302     cerr << "No 3DNow! support; not tested\n";
 303   }
 304   else
 305     t2_base (complex_dotprod_3dnow);
 306 }
 307
 308 void
 309 qa_complex_dotprod_x86::t3_3dnow ()
 310 {
 311   if (!gr_cpu::has_3dnow ()){
 312     cerr << "No 3DNow! support; not tested\n";
 313   }
 314   else
 315     t3_base (complex_dotprod_3dnow);
 316 }
 317
 318 void
 319 qa_complex_dotprod_x86::t1_sse ()
 320 {
 321   if (!gr_cpu::has_sse ()){
 322     cerr << "No SSE support; not tested\n";
 323   }
 324   else
 325     t1_base (complex_dotprod_sse);
 326 }
 327
 328 void
 329 qa_complex_dotprod_x86::t2_sse ()
 330 {
 331   if (!gr_cpu::has_sse ()){
 332     cerr << "No SSE support; not tested\n";
 333   }
 334   else
 335     t2_base (complex_dotprod_sse);
 336 }
 337
 338 void
 339 qa_complex_dotprod_x86::t3_sse ()
 340 {
 341   if (!gr_cpu::has_sse ()){
 342     cerr << "No SSE support; not tested\n";
 343   }
 344   else
 345     t3_base (complex_dotprod_sse);
 346 }
 347