49 #pragma warning (disable: 4244)
54 #include "sphinxbase/byteorder.h"
55 #include "sphinxbase/fixpoint.h"
56 #include "sphinxbase/fe.h"
61 #include "fe_internal.h"
67 #define FLOAT2COS(x) FLOAT2FIX_ANY(x,30)
68 #define COSMUL(x,y) FIXMUL_ANY(x,y,30)
70 #define FLOAT2COS(x) (x)
71 #define COSMUL(x,y) ((x)*(y))
86 static const unsigned char fe_logadd_table[] = {
87 177, 177, 176, 176, 175, 175, 174, 174, 173, 173,
88 172, 172, 172, 171, 171, 170, 170, 169, 169, 168,
89 168, 167, 167, 166, 166, 165, 165, 164, 164, 163,
90 163, 162, 162, 161, 161, 161, 160, 160, 159, 159,
91 158, 158, 157, 157, 156, 156, 155, 155, 155, 154,
92 154, 153, 153, 152, 152, 151, 151, 151, 150, 150,
93 149, 149, 148, 148, 147, 147, 147, 146, 146, 145,
94 145, 144, 144, 144, 143, 143, 142, 142, 141, 141,
95 141, 140, 140, 139, 139, 138, 138, 138, 137, 137,
96 136, 136, 136, 135, 135, 134, 134, 134, 133, 133,
97 132, 132, 131, 131, 131, 130, 130, 129, 129, 129,
98 128, 128, 128, 127, 127, 126, 126, 126, 125, 125,
99 124, 124, 124, 123, 123, 123, 122, 122, 121, 121,
100 121, 120, 120, 119, 119, 119, 118, 118, 118, 117,
101 117, 117, 116, 116, 115, 115, 115, 114, 114, 114,
102 113, 113, 113, 112, 112, 112, 111, 111, 110, 110,
103 110, 109, 109, 109, 108, 108, 108, 107, 107, 107,
104 106, 106, 106, 105, 105, 105, 104, 104, 104, 103,
105 103, 103, 102, 102, 102, 101, 101, 101, 100, 100,
106 100, 99, 99, 99, 98, 98, 98, 97, 97, 97,
107 96, 96, 96, 96, 95, 95, 95, 94, 94, 94,
108 93, 93, 93, 92, 92, 92, 92, 91, 91, 91,
109 90, 90, 90, 89, 89, 89, 89, 88, 88, 88,
110 87, 87, 87, 87, 86, 86, 86, 85, 85, 85,
111 85, 84, 84, 84, 83, 83, 83, 83, 82, 82,
112 82, 82, 81, 81, 81, 80, 80, 80, 80, 79,
113 79, 79, 79, 78, 78, 78, 78, 77, 77, 77,
114 77, 76, 76, 76, 75, 75, 75, 75, 74, 74,
115 74, 74, 73, 73, 73, 73, 72, 72, 72, 72,
116 71, 71, 71, 71, 71, 70, 70, 70, 70, 69,
117 69, 69, 69, 68, 68, 68, 68, 67, 67, 67,
118 67, 67, 66, 66, 66, 66, 65, 65, 65, 65,
119 64, 64, 64, 64, 64, 63, 63, 63, 63, 63,
120 62, 62, 62, 62, 61, 61, 61, 61, 61, 60,
121 60, 60, 60, 60, 59, 59, 59, 59, 59, 58,
122 58, 58, 58, 58, 57, 57, 57, 57, 57, 56,
123 56, 56, 56, 56, 55, 55, 55, 55, 55, 54,
124 54, 54, 54, 54, 53, 53, 53, 53, 53, 52,
125 52, 52, 52, 52, 52, 51, 51, 51, 51, 51,
126 50, 50, 50, 50, 50, 50, 49, 49, 49, 49,
127 49, 49, 48, 48, 48, 48, 48, 48, 47, 47,
128 47, 47, 47, 47, 46, 46, 46, 46, 46, 46,
129 45, 45, 45, 45, 45, 45, 44, 44, 44, 44,
130 44, 44, 43, 43, 43, 43, 43, 43, 43, 42,
131 42, 42, 42, 42, 42, 41, 41, 41, 41, 41,
132 41, 41, 40, 40, 40, 40, 40, 40, 40, 39,
133 39, 39, 39, 39, 39, 39, 38, 38, 38, 38,
134 38, 38, 38, 37, 37, 37, 37, 37, 37, 37,
135 37, 36, 36, 36, 36, 36, 36, 36, 35, 35,
136 35, 35, 35, 35, 35, 35, 34, 34, 34, 34,
137 34, 34, 34, 34, 33, 33, 33, 33, 33, 33,
138 33, 33, 32, 32, 32, 32, 32, 32, 32, 32,
139 32, 31, 31, 31, 31, 31, 31, 31, 31, 31,
140 30, 30, 30, 30, 30, 30, 30, 30, 30, 29,
141 29, 29, 29, 29, 29, 29, 29, 29, 28, 28,
142 28, 28, 28, 28, 28, 28, 28, 28, 27, 27,
143 27, 27, 27, 27, 27, 27, 27, 27, 26, 26,
144 26, 26, 26, 26, 26, 26, 26, 26, 25, 25,
145 25, 25, 25, 25, 25, 25, 25, 25, 25, 24,
146 24, 24, 24, 24, 24, 24, 24, 24, 24, 24,
147 23, 23, 23, 23, 23, 23, 23, 23, 23, 23,
148 23, 23, 22, 22, 22, 22, 22, 22, 22, 22,
149 22, 22, 22, 22, 21, 21, 21, 21, 21, 21,
150 21, 21, 21, 21, 21, 21, 21, 20, 20, 20,
151 20, 20, 20, 20, 20, 20, 20, 20, 20, 20,
152 19, 19, 19, 19, 19, 19, 19, 19, 19, 19,
153 19, 19, 19, 19, 18, 18, 18, 18, 18, 18,
154 18, 18, 18, 18, 18, 18, 18, 18, 18, 17,
155 17, 17, 17, 17, 17, 17, 17, 17, 17, 17,
156 17, 17, 17, 17, 16, 16, 16, 16, 16, 16,
157 16, 16, 16, 16, 16, 16, 16, 16, 16, 16,
158 16, 15, 15, 15, 15, 15, 15, 15, 15, 15,
159 15, 15, 15, 15, 15, 15, 15, 15, 14, 14,
160 14, 14, 14, 14, 14, 14, 14, 14, 14, 14,
161 14, 14, 14, 14, 14, 14, 14, 13, 13, 13,
162 13, 13, 13, 13, 13, 13, 13, 13, 13, 13,
163 13, 13, 13, 13, 13, 13, 13, 12, 12, 12,
164 12, 12, 12, 12, 12, 12, 12, 12, 12, 12,
165 12, 12, 12, 12, 12, 12, 12, 12, 12, 11,
166 11, 11, 11, 11, 11, 11, 11, 11, 11, 11,
167 11, 11, 11, 11, 11, 11, 11, 11, 11, 11,
168 11, 11, 11, 10, 10, 10, 10, 10, 10, 10,
169 10, 10, 10, 10, 10, 10, 10, 10, 10, 10,
170 10, 10, 10, 10, 10, 10, 10, 10, 10, 9,
171 9, 9, 9, 9, 9, 9, 9, 9, 9, 9,
172 9, 9, 9, 9, 9, 9, 9, 9, 9, 9,
173 9, 9, 9, 9, 9, 9, 9, 9, 8, 8,
174 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
175 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
176 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
177 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
178 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
179 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
180 7, 7, 7, 7, 7, 7, 7, 7, 6, 6,
181 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
182 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
183 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
184 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
185 6, 5, 5, 5, 5, 5, 5, 5, 5, 5,
186 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
187 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
188 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
189 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
190 5, 5, 5, 4, 4, 4, 4, 4, 4, 4,
191 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
192 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
193 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
194 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
195 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
196 4, 4, 4, 4, 4, 4, 4, 4, 3, 3,
197 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
198 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
199 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
200 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
201 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
202 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
203 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
204 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
205 3, 3, 3, 3, 2, 2, 2, 2, 2, 2,
206 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
207 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
208 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
209 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
210 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
211 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
212 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
213 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
214 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
215 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
216 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
217 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
218 2, 2, 2, 2, 2, 2, 1, 1, 1, 1,
219 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
220 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
221 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
222 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
223 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
224 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
225 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
226 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
227 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
228 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
229 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
230 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
231 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
232 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
233 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
234 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
235 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
236 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
237 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
238 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
239 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
240 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
241 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
242 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
243 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
244 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
245 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
246 1, 1, 1, 1, 1, 1, 1, 0
248 static const int fe_logadd_table_size =
sizeof(fe_logadd_table) /
sizeof(fe_logadd_table[0]);
251 fe_log_add(fixed32 x, fixed32 y)
256 d = (x - y) >> (DEFAULT_RADIX - 8);
260 d = (y - x) >> (DEFAULT_RADIX - 8);
263 if (d > fe_logadd_table_size)
266 r += ((fixed32)fe_logadd_table[d] << (DEFAULT_RADIX - 8));
283 return FLOAT2FIX(log(x));
289 fe_mel(
melfb_t *mel, float32 x)
291 float32 warped = fe_warp_unwarped_to_warped(mel, x);
293 return (float32) (2595.0 * log10(1.0 + warped / 700.0));
297 fe_melinv(
melfb_t *mel, float32 x)
299 float32 warped = (float32) (700.0 * (pow(10.0, x / 2595.0) - 1.0));
300 return fe_warp_warped_to_unwarped(mel, warped);
304 fe_build_melfilters(
melfb_t *mel_fb)
306 float32 melmin, melmax, melbw, fftfreq;
310 mel_fb->spec_start =
ckd_malloc(mel_fb->num_filters *
sizeof(*mel_fb->spec_start));
311 mel_fb->filt_start =
ckd_malloc(mel_fb->num_filters *
sizeof(*mel_fb->filt_start));
312 mel_fb->filt_width =
ckd_malloc(mel_fb->num_filters *
sizeof(*mel_fb->filt_width));
316 melmin = fe_mel(mel_fb, mel_fb->lower_filt_freq);
317 melmax = fe_mel(mel_fb, mel_fb->upper_filt_freq);
320 melbw = (melmax - melmin) / (mel_fb->num_filters + 1);
321 if (mel_fb->doublewide) {
324 if ((fe_melinv(mel_fb, melmin) < 0) ||
325 (fe_melinv(mel_fb, melmax) > mel_fb->sampling_rate / 2)) {
327 (
"Out of Range: low filter edge = %f (%f)\n",
328 fe_melinv(mel_fb, melmin), 0.0);
330 (
" high filter edge = %f (%f)\n",
331 fe_melinv(mel_fb, melmax), mel_fb->sampling_rate / 2);
332 return FE_INVALID_PARAM_ERROR;
337 fftfreq = mel_fb->sampling_rate / (float32) mel_fb->fft_size;
341 for (i = 0; i < mel_fb->num_filters; ++i) {
345 for (j = 0; j < 3; ++j) {
346 if (mel_fb->doublewide)
347 freqs[j] = fe_melinv(mel_fb, (i + j * 2) * melbw + melmin);
349 freqs[j] = fe_melinv(mel_fb, (i + j) * melbw + melmin);
351 if (mel_fb->round_filters)
352 freqs[j] = ((int)(freqs[j] / fftfreq + 0.5)) * fftfreq;
356 mel_fb->spec_start[i] = -1;
358 for (j = 0; j < mel_fb->fft_size/2+1; ++j) {
359 float32 hz = j * fftfreq;
362 else if (hz > freqs[2] || j == mel_fb->fft_size/2) {
364 mel_fb->filt_width[i] = j - mel_fb->spec_start[i];
366 mel_fb->filt_start[i] = n_coeffs;
367 n_coeffs += mel_fb->filt_width[i];
370 if (mel_fb->spec_start[i] == -1)
371 mel_fb->spec_start[i] = j;
376 mel_fb->filt_coeffs =
ckd_malloc(n_coeffs *
sizeof(*mel_fb->filt_coeffs));
380 for (i = 0; i < mel_fb->num_filters; ++i) {
384 for (j = 0; j < 3; ++j) {
385 if (mel_fb->doublewide)
386 freqs[j] = fe_melinv(mel_fb, (i + j * 2) * melbw + melmin);
388 freqs[j] = fe_melinv(mel_fb, (i + j) * melbw + melmin);
390 if (mel_fb->round_filters)
391 freqs[j] = ((int)(freqs[j] / fftfreq + 0.5)) * fftfreq;
394 for (j = 0; j < mel_fb->filt_width[i]; ++j) {
395 float32 hz, loslope, hislope;
397 hz = (mel_fb->spec_start[i] + j) * fftfreq;
398 if (hz < freqs[0] || hz > freqs[2]) {
399 E_FATAL(
"WTF, %f < %f > %f\n", freqs[0], hz, freqs[2]);
401 loslope = (hz - freqs[0]) / (freqs[1] - freqs[0]);
402 hislope = (freqs[2] - hz) / (freqs[2] - freqs[1]);
403 if (mel_fb->unit_area) {
404 loslope *= 2 / (freqs[2] - freqs[0]);
405 hislope *= 2 / (freqs[2] - freqs[0]);
407 if (loslope < hislope) {
409 mel_fb->filt_coeffs[n_coeffs] = fe_log(loslope);
411 mel_fb->filt_coeffs[n_coeffs] = loslope;
416 mel_fb->filt_coeffs[n_coeffs] = fe_log(hislope);
418 mel_fb->filt_coeffs[n_coeffs] = hislope;
430 fe_compute_melcosine(
melfb_t * mel_fb)
441 freqstep = M_PI / mel_fb->num_filters;
444 for (i = 0; i < mel_fb->num_cepstra; i++) {
445 for (j = 0; j < mel_fb->num_filters; j++) {
448 cosine = cos(freqstep * i * (j + 0.5));
449 mel_fb->mel_cosine[i][j] = FLOAT2COS(cosine);
454 mel_fb->sqrt_inv_n = FLOAT2COS(sqrt(1.0 / mel_fb->num_filters));
455 mel_fb->sqrt_inv_2n = FLOAT2COS(sqrt(2.0 / mel_fb->num_filters));
458 if (mel_fb->lifter_val) {
459 mel_fb->lifter = calloc(mel_fb->num_cepstra,
sizeof(*mel_fb->lifter));
460 for (i = 0; i < mel_fb->num_cepstra; ++i) {
461 mel_fb->lifter[i] = FLOAT2MFCC(1 + mel_fb->lifter_val / 2
462 * sin(i * M_PI / mel_fb->lifter_val));
470 fe_pre_emphasis(int16
const *in, frame_t * out, int32 len,
471 float32 factor, int16 prior)
476 int16 fxd_alpha = (int16)(factor * 0x8000);
479 tmp1 = (int32)in[0] << 15;
480 tmp2 = (int32)prior * fxd_alpha;
481 out[0] = (int16)((tmp1 - tmp2) >> 15);
482 for (i = 1; i < len; ++i) {
483 tmp1 = (int32)in[i] << 15;
484 tmp2 = (int32)in[i-1] * fxd_alpha;
485 out[i] = (int16)((tmp1 - tmp2) >> 15);
487 #elif defined(FIXED_POINT)
488 fixed32 fxd_alpha = FLOAT2FIX(factor);
489 out[0] = ((fixed32)in[0] << DEFAULT_RADIX) - (prior * fxd_alpha);
490 for (i = 1; i < len; ++i)
491 out[i] = ((fixed32)in[i] << DEFAULT_RADIX)
492 - (fixed32)in[i-1] * fxd_alpha;
494 out[0] = (frame_t) in[0] - (frame_t) prior * factor;
495 for (i = 1; i < len; i++)
496 out[i] = (frame_t) in[i] - (frame_t) in[i-1] * factor;
501 fe_short_to_frame(int16
const *in, frame_t * out, int32 len)
506 memcpy(out, in, len *
sizeof(*out));
507 #elif defined(FIXED_POINT)
508 for (i = 0; i < len; i++)
509 out[i] = (int32) in[i] << DEFAULT_RADIX;
511 for (i = 0; i < len; i++)
512 out[i] = (frame_t) in[i];
517 fe_create_hamming(window_t * in, int32 in_len)
522 for (i = 0; i < in_len / 2; i++) {
524 hamm = (0.54 - 0.46 * cos(2 * M_PI * i /
525 ((float64) in_len - 1.0)));
527 in[i] = (int16)(hamm * 0x8000);
529 in[i] = FLOAT2COS(hamm);
535 fe_hamming_window(frame_t * in, window_t * window, int32 in_len, int32 remove_dc)
546 for (i = 0; i < in_len; i++)
549 for (i = 0; i < in_len; i++)
550 in[i] -= (frame_t)mean;
554 for (i = 0; i < in_len/2; i++) {
557 tmp1 = (int32)in[i] * window[i];
558 tmp2 = (int32)in[in_len-1-i] * window[i];
559 in[i] = (int16)(tmp1 >> 15);
560 in[in_len-1-i] = (int16)(tmp2 >> 15);
563 for (i = 0; i < in_len/2; i++) {
564 in[i] = COSMUL(in[i], window[i]);
565 in[in_len-1-i] = COSMUL(in[in_len-1-i], window[i]);
571 fe_spch_to_frame(
fe_t *fe,
int len)
574 if (fe->pre_emphasis_alpha != 0.0) {
575 fe_pre_emphasis(fe->spch, fe->frame, len,
576 fe->pre_emphasis_alpha, fe->prior);
577 if (len >= fe->frame_shift)
578 fe->prior = fe->spch[fe->frame_shift - 1];
580 fe->prior = fe->spch[len - 1];
583 fe_short_to_frame(fe->spch, fe->frame, len);
586 memset(fe->frame + len, 0,
587 (fe->fft_size - len) *
sizeof(*fe->frame));
590 fe_hamming_window(fe->frame, fe->hamming_window, fe->frame_size,
597 fe_read_frame(
fe_t *fe, int16
const *in, int32 len)
601 if (len > fe->frame_size)
602 len = fe->frame_size;
605 memcpy(fe->spch, in, len *
sizeof(*in));
608 for (i = 0; i < len; ++i)
609 SWAP_INT16(&fe->spch[i]);
611 for (i = 0; i < len; ++i)
612 fe->spch[i] += (int16) ((!(s3_rand_int31() % 4)) ? 1 : 0);
614 return fe_spch_to_frame(fe, len);
618 fe_shift_frame(
fe_t *fe, int16
const *in, int32 len)
622 if (len > fe->frame_shift)
623 len = fe->frame_shift;
624 offset = fe->frame_size - fe->frame_shift;
627 memmove(fe->spch, fe->spch + fe->frame_shift,
628 offset *
sizeof(*fe->spch));
629 memcpy(fe->spch + offset, in, len *
sizeof(*fe->spch));
632 for (i = 0; i < len; ++i)
633 SWAP_INT16(&fe->spch[offset + i]);
635 for (i = 0; i < len; ++i)
637 += (int16) ((!(s3_rand_int31() % 4)) ? 1 : 0);
639 return fe_spch_to_frame(fe, offset + len);
646 fe_create_twiddle(
fe_t *fe)
650 for (i = 0; i < fe->fft_size / 4; ++i) {
651 float64 a = 2 * M_PI * i / fe->fft_size;
653 fe->ccc[i] = (int16)(cos(a) * 0x8000);
654 fe->sss[i] = (int16)(sin(a) * 0x8000);
655 #elif defined(FIXED_POINT)
656 fe->ccc[i] = FLOAT2COS(cos(a));
657 fe->sss[i] = FLOAT2COS(sin(a));
673 fe_fft_real(
fe_t *fe)
675 int i, j, k, m, n, lz;
684 for (i = 0; i < n - 1; ++i) {
699 for (i = 0; i < n; ++i)
705 for (lz = 0; lz < m; ++lz)
706 if (max & (1 << (15-lz)))
716 for (i = 0; i < n; i += 2) {
717 int atten = (lz == 0);
719 x[i] = xt + (x[i + 1] >> atten);
720 x[i + 1] = xt - (x[i + 1] >> atten);
724 for (k = 1; k < m; ++k) {
727 int atten = (k >= lz);
733 for (i = 0; i < n; i += (1 << n1)) {
739 x[i] = xt + (x[i + (1 << n2)] >> atten);
740 x[i + (1 << n2)] = xt - (x[i + (1 << n2)] >> atten);
748 x[i + (1 << n2) + (1 << n4)] = -x[i + (1 << n2) + (1 << n4)] >> atten;
749 x[i + (1 << n4)] = x[i + (1 << n4)] >> atten;
754 for (j = 1; j < (1 << n4); ++j) {
755 frame_t cc, ss, t1, t2;
759 i2 = i + (1 << n2) - j;
760 i3 = i + (1 << n2) + j;
761 i4 = i + (1 << n2) + (1 << n2) - j;
767 cc = fe->ccc[j << (m - n1)];
768 ss = fe->sss[j << (m - n1)];
774 tmp1 = (int32)x[i3] * cc + (int32)x[i4] * ss;
775 tmp2 = (int32)x[i3] * ss - (int32)x[i4] * cc;
776 t1 = (int16)(tmp1 >> 15) >> atten;
777 t2 = (int16)(tmp2 >> 15) >> atten;
780 x[i4] = (x[i2] >> atten) - t2;
781 x[i3] = (-x[i2] >> atten) - t2;
782 x[i2] = (x[i1] >> atten) - t1;
783 x[i1] = (x[i1] >> atten) + t1;
793 fe_fft_real(
fe_t *fe)
804 for (i = 0; i < n - 1; ++i) {
822 for (i = 0; i < n; i += 2) {
824 x[i] = (xt + x[i + 1]);
825 x[i + 1] = (xt - x[i + 1]);
829 for (k = 1; k < m; ++k) {
836 for (i = 0; i < n; i += (1 << n1)) {
842 x[i] = (xt + x[i + (1 << n2)]);
843 x[i + (1 << n2)] = (xt - x[i + (1 << n2)]);
851 x[i + (1 << n2) + (1 << n4)] = -x[i + (1 << n2) + (1 << n4)];
852 x[i + (1 << n4)] = x[i + (1 << n4)];
857 for (j = 1; j < (1 << n4); ++j) {
858 frame_t cc, ss, t1, t2;
862 i2 = i + (1 << n2) - j;
863 i3 = i + (1 << n2) + j;
864 i4 = i + (1 << n2) + (1 << n2) - j;
870 cc = fe->ccc[j << (m - n1)];
871 ss = fe->sss[j << (m - n1)];
875 t1 = COSMUL(x[i3], cc) + COSMUL(x[i4], ss);
876 t2 = COSMUL(x[i3], ss) - COSMUL(x[i4], cc);
878 x[i4] = (x[i2] - t2);
879 x[i3] = (-x[i2] - t2);
880 x[i2] = (x[i1] - t1);
881 x[i1] = (x[i1] + t1);
892 fe_spec_magnitude(
fe_t *fe)
896 int32 j, scale, fftsize;
900 scale = fe_fft_real(fe);
905 fftsize = fe->fft_size;
908 scale = fe->fft_order - scale;
913 spec[0] = fixlog(abs(fft[0]) << scale) * 2;
914 #elif defined(FIXED_POINT)
915 spec[0] = FIXLN(abs(fft[0]) << scale) * 2;
917 spec[0] = fft[0] * fft[0];
921 for (j = 1; j <= fftsize / 2; j++) {
923 int32 rr = fixlog(abs(fft[j]) << scale) * 2;
924 int32 ii = fixlog(abs(fft[fftsize - j]) << scale) * 2;
925 spec[j] = fe_log_add(rr, ii);
926 #elif defined(FIXED_POINT)
927 int32 rr = FIXLN(abs(fft[j]) << scale) * 2;
928 int32 ii = FIXLN(abs(fft[fftsize - j]) << scale) * 2;
929 spec[j] = fe_log_add(rr, ii);
931 spec[j] = fft[j] * fft[j] + fft[fftsize - j] * fft[fftsize - j];
937 fe_mel_spec(
fe_t * fe)
940 powspec_t *spec, *mfspec;
946 for (whichfilt = 0; whichfilt < fe->mel_fb->num_filters; whichfilt++) {
947 int spec_start, filt_start, i;
949 spec_start = fe->mel_fb->spec_start[whichfilt];
950 filt_start = fe->mel_fb->filt_start[whichfilt];
953 mfspec[whichfilt] = spec[spec_start] + fe->mel_fb->filt_coeffs[filt_start];
954 for (i = 1; i < fe->mel_fb->filt_width[whichfilt]; i++) {
955 mfspec[whichfilt] = fe_log_add(mfspec[whichfilt],
956 spec[spec_start + i] +
957 fe->mel_fb->filt_coeffs[filt_start + i]);
960 mfspec[whichfilt] = 0;
961 for (i = 0; i < fe->mel_fb->filt_width[whichfilt]; i++)
963 spec[spec_start + i] * fe->mel_fb->filt_coeffs[filt_start + i];
969 fe_mel_cep(
fe_t * fe, mfcc_t *mfcep)
977 for (i = 0; i < fe->mel_fb->num_filters; ++i) {
980 mfspec[i] = log(mfspec[i]);
991 if (fe->log_spec == RAW_LOG_SPEC) {
992 for (i = 0; i < fe->feature_dimension; i++) {
993 mfcep[i] = (mfcc_t) mfspec[i];
997 else if (fe->log_spec == SMOOTH_LOG_SPEC) {
999 fe_dct2(fe, mfspec, mfcep, 0);
1000 fe_dct3(fe, mfcep, mfspec);
1001 for (i = 0; i < fe->feature_dimension; i++) {
1002 mfcep[i] = (mfcc_t) mfspec[i];
1005 else if (fe->transform == DCT_II)
1006 fe_dct2(fe, mfspec, mfcep, FALSE);
1007 else if (fe->transform == DCT_HTK)
1008 fe_dct2(fe, mfspec, mfcep, TRUE);
1010 fe_spec2cep(fe, mfspec, mfcep);
1016 fe_spec2cep(
fe_t * fe,
const powspec_t * mflogspec, mfcc_t * mfcep)
1022 mfcep[0] = mflogspec[0] / 2;
1023 for (j = 1; j < fe->mel_fb->num_filters; j++)
1024 mfcep[0] += mflogspec[j];
1025 mfcep[0] /= (frame_t) fe->mel_fb->num_filters;
1027 for (i = 1; i < fe->num_cepstra; ++i) {
1029 for (j = 0; j < fe->mel_fb->num_filters; j++) {
1034 mfcep[i] += COSMUL(mflogspec[j],
1035 fe->mel_fb->mel_cosine[i][j]) * beta;
1040 mfcep[i] /= (frame_t) fe->mel_fb->num_filters * 2;
1045 fe_dct2(
fe_t * fe,
const powspec_t * mflogspec, mfcc_t * mfcep,
int htk)
1051 mfcep[0] = mflogspec[0];
1052 for (j = 1; j < fe->mel_fb->num_filters; j++)
1053 mfcep[0] += mflogspec[j];
1055 mfcep[0] = COSMUL(mfcep[0], fe->mel_fb->sqrt_inv_2n);
1057 mfcep[0] = COSMUL(mfcep[0], fe->mel_fb->sqrt_inv_n);
1059 for (i = 1; i < fe->num_cepstra; ++i) {
1061 for (j = 0; j < fe->mel_fb->num_filters; j++) {
1062 mfcep[i] += COSMUL(mflogspec[j],
1063 fe->mel_fb->mel_cosine[i][j]);
1065 mfcep[i] = COSMUL(mfcep[i], fe->mel_fb->sqrt_inv_2n);
1070 fe_lifter(
fe_t *fe, mfcc_t *mfcep)
1074 if (fe->mel_fb->lifter_val == 0)
1077 for (i = 0; i < fe->num_cepstra; ++i) {
1078 mfcep[i] = MFCCMUL(mfcep[i], fe->mel_fb->lifter[i]);
1083 fe_dct3(
fe_t * fe,
const mfcc_t * mfcep, powspec_t * mflogspec)
1087 for (i = 0; i < fe->mel_fb->num_filters; ++i) {
1088 mflogspec[i] = COSMUL(mfcep[0], SQRT_HALF);
1089 for (j = 1; j < fe->num_cepstra; j++) {
1090 mflogspec[i] += COSMUL(mfcep[j],
1091 fe->mel_fb->mel_cosine[j][i]);
1093 mflogspec[i] = COSMUL(mflogspec[i], fe->mel_fb->sqrt_inv_2n);
1098 fe_write_frame(
fe_t * fe, mfcc_t * fea)
1100 fe_spec_magnitude(fe);
1102 fe_mel_cep(fe, fea);
1109 fe_create_2d(int32 d1, int32 d2, int32 elem_size)
1115 fe_free_2d(
void *arr)