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
MaxSignalType(int real_only)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 */
GenerateTestSignalAndFFT(struct ComplexFloat * x,struct ComplexFloat * fft,int size,int signal_type,float signal_value,int real_only)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