Home | History | Annotate | Download | only in test 
      1 /*
      2  *  Copyright (c) 2013 The WebRTC project authors. All Rights Reserved.
      3  *
      4  *  Use of this source code is governed by a BSD-style license
      5  *  that can be found in the LICENSE file in the root of the source
      6  *  tree. An additional intellectual property rights grant can be found
      7  *  in the file PATENTS.  All contributing project authors may
      8  *  be found in the AUTHORS file in the root of the source tree.
      9  */
     10 
     11 #include "dl/sp/src/test/gensig.h"
     12 
     13 #include <stdio.h>
     14 #include <math.h>
     15 #include <stdlib.h>
     16 
     17 #define MAX_REAL_SIGNAL_TYPE 3
     18 #define MAX_SIGNAL_TYPE 4
     19 
     20 int MaxSignalType(int real_only) {
     21   return real_only ? MAX_REAL_SIGNAL_TYPE : MAX_SIGNAL_TYPE;
     22 }
     23 
     24 /*
     25  * Generate a test signal and compute the theoretical FFT.
     26  *
     27  * The test signal is specified by |signal_type|, and the test signal
     28  * is saved in |x| with the corresponding FFT in |fft|.  The size of
     29  * the test signal is |size|.  |signalValue| is desired the amplitude
     30  * of the test signal.
     31  *
     32  * If |real_only| is true, then the test signal is assumed to be real
     33  * instead of complex, which is the default.  This is only applicable
     34  * for a |signal_type| of 0 or 3; the other signals are already real-valued.
     35  */
     36 void GenerateTestSignalAndFFT(struct ComplexFloat* x,
     37                               struct ComplexFloat* fft,
     38                               int size,
     39                               int signal_type,
     40                               float signal_value,
     41                               int real_only) {
     42   int k;
     43 
     44   switch (signal_type) {
     45     case 0:
     46       /*
     47        * Signal is a constant signal_value + i*signal_value (or just
     48        * signal_value if real_only is true.)
     49        */
     50       for (k = 0; k < size; ++k) {
     51         x[k].Re = signal_value;
     52         x[k].Im = real_only ? 0 : signal_value;
     53       }
     54 
     55       fft[0].Re = signal_value * size;
     56       fft[0].Im = real_only ? 0 : signal_value * size;
     57 
     58       for (k = 1; k < size; ++k) {
     59         fft[k].Re = fft[k].Im = 0;
     60       }
     61       break;
     62     case 1:
     63       /*
     64        * A real-valued ramp
     65        */
     66       {
     67         double factor = signal_value / (float) size;
     68         double omega = 2 * M_PI / size;
     69 
     70         for (k = 0; k < size; ++k) {
     71           x[k].Re = ((k + 1)*factor);
     72           x[k].Im = 0;
     73         }
     74 
     75         fft[0].Re = factor * size * (size + 1) / 2;
     76         fft[0].Im = 0;
     77         for (k = 1; k < size; ++k) {
     78           double phase;
     79           phase = omega * k;
     80           fft[k].Re = factor * -size / 2;
     81           fft[k].Im = factor * size / 2 * (sin(phase) / (1 - cos(phase)));
     82         }
     83 
     84         /*
     85          * Remove any roundoff for k = N/2 since sin(2*pi/N*N/2) = 0.
     86          */
     87         fft[size / 2].Im = 0;
     88       }
     89       break;
     90     case 2:
     91       /*
     92        * Pure real-valued sine wave, one cycle.
     93        */
     94       {
     95         double omega = 2 * M_PI / size;
     96 
     97         for (k = 0; k < size; ++k) {
     98           x[k].Re = signal_value * sin(omega * k);
     99           x[k].Im = 0;
    100         }
    101 
    102         /*
    103          * Remove any roundoff for k = N/2 since sin(2*pi/N*N/2) = 0.
    104          */
    105         x[size / 2 ].Re = 0;
    106 
    107         for (k = 0; k < size; ++k) {
    108           fft[k].Re = 0;
    109           fft[k].Im = 0;
    110         }
    111 
    112         /*
    113          * When size == 2, x[k] is identically zero, so the FFT is also zero.
    114          */
    115         if (size != 2) {
    116           fft[1].Im = -signal_value * (size / 2);
    117           fft[size - 1].Im = signal_value * (size / 2);
    118         }
    119       }
    120       break;
    121     case 3:
    122       /*
    123        * The source signal is x[k] = 0 except x[1] = x[size-1] =
    124        * -i*signal_value.  The transform is one period of a pure real
    125        * (negative) sine wave.  Only defined when real_only is false.
    126        */
    127       if (!real_only) {
    128         double omega = 2 * M_PI / size;
    129         for (k = 0; k < size; ++k) {
    130           x[k].Re = 0;
    131           x[k].Im = 0;
    132         }
    133         x[1].Im = -signal_value;
    134         x[size-1].Im = signal_value;
    135 
    136         if (size == 2) {
    137           fft[0].Re = 0;
    138           fft[0].Im = signal_value;
    139           fft[1].Re = 0;
    140           fft[1].Im = -signal_value;
    141         } else {
    142           for (k = 0; k < size; ++k) {
    143             fft[k].Re = -2 * signal_value * sin(omega * k);
    144             fft[k].Im = 0;
    145           }
    146 
    147           /*
    148            * Remove any roundoff for k = N/2 since sin(2*pi/N*N/2) = 0.
    149            */
    150           fft[size / 2].Re = 0;
    151         }
    152         break;
    153       }
    154       /* Fall through if real_only */
    155     case MAX_SIGNAL_TYPE:
    156     default:
    157       fprintf(stderr, "invalid signal type: %d\n", signal_type);
    158       exit(1);
    159   }
    160 }
    161