带有或不带有窗口的KISS FFT输出 [英] KISS FFT output with or without windowing
问题描述
我目前正在尝试使用Kiss FFT将FFT实现到AVR32微控制器中，以进行信号处理. 而且我的输出有一个奇怪的问题. 基本上，我将ADC样本(使用函数发生器进行测试)传递到fft(真实输入，256 n大小)中，并且检索到的输出对我来说很有意义. 但是，如果我将汉明窗应用于ADC样本，然后将其传递给FFT，则峰值幅度的频率仓是错误的(并且与之前没有开窗的结果不同). ADC样本具有DC偏移，因此我消除了偏移，但仍不适用于加窗样本.
以下是通过rs485输出的前几个输出值. 第一列是不带窗口的fft输出，而第二列是带窗口的输出.从第1列开始，峰值位于第6行(6 x fs(10.5kHz)/0.5N)为我提供了正确的输入频率结果，其中第2列在第2行具有峰值幅度(直流仓除外)，这对我来说没有意义. 任何建议都会有所帮助. 预先感谢.
488 260 //dc bin
5 97
5 41
5 29
4 26
10 35
133 76
33 28
21 6
17 3
kiss_fft_scalar zero;
memset(&zero,0,sizeof(zero));
kiss_fft_cpx fft_input[n];
kiss_fft_cpx fft_output[n];
for(ctr=0; ctr<n; ctr++)
{
fft_input[ctr].r = zero;
fft_input[ctr].i = zero;
fft_output[ctr].r =zero;
fft_output[ctr].i = zero;
}
// IIR filter calculation
for (ctr=0; ctr<n; ctr++)
{
// filter calculation
y[ctr] = num_coef[0]*x[ctr];
y[ctr] += (num_coef[1]*x[ctr-1]) - (den_coef[1]*y[ctr-1]);
y[ctr] += (num_coef[2]*x[ctr-2]) - (den_coef[2]*y[ctr-2]);
//y1[ctr] += y[ctr] - 500;
// hamming window
hamming[ctr] = (0.54-((0.46) * cos(2*PI*ctr/256)));
window[ctr] = hamming[ctr]*y[ctr];
fft_input[ctr].r = window[ctr];
fft_input[ctr].i = 0;
fft_output[ctr].r = 0;
fft_output[ctr].i = 0;
}
kiss_fftr_cfg fftConfig = kiss_fftr_alloc(n,0,NULL,NULL);
kiss_fftr(fftConfig, (kiss_fft_scalar * )fft_input, fft_output);
for (ctr=0; ctr<n; ctr++)
{
fft_mag[ctr] = (sqrt((fft_output[ctr].r * fft_output[ctr].r) + (fft_output[ctr].i * fft_output[ctr].i)))/(0.5*n);
//Usart write
char filtResult[10];
sprintf(filtResult, "%04d %04d\n", (int)x[ctr], (int)fft_mag[ctr]);
//sprintf(filtResult, "%04d %04d\n", (int)x[ctr], (int)window[ctr]);
char c;
char *ptr = &filtResult[0];
do
{
c = *ptr;
ptr++;
usart_bw_write_char(&AVR32_USART2, (int)c);
// sendByte(c);
} while (c != '\n');
}
kiss_fft_cleanup();
free(fftConfig);
频率域输出说明
在频域中，矩形和汉明窗看起来像:
您可能已经知道，时域乘以一个窗口对应于频域中的卷积，这实际上将信号的能量散布在通常称为
最后，请注意，IIR过滤器的实现方式还有一个问题，即当 另一方面，如果您想过滤 然后，您可以在这些缓冲区中使用偏移量.索引0和1代表过去的样本，而索引2中的其余缓冲区则用当前输入数据块填充.这将导致以下经过稍微修改的过滤代码: 最后，您必须在每个块的末尾更新索引为0& 0的过去样本". 1，准备好处理当前块的最后一个样本，以处理下一个输入块: I am currently trying to implement fft into avr32 micro controllers for signal processing purposes using kiss fft.
And having a strange problem with my output.
Basically, I am passing ADC samples (testing with function generator) into fft (real input, 256 n size) and retrieved output makes sense to me.
However, if I apply Hamming window to ADC samples and then pass these to FFT, the frequency bin of the peak magnitude is wrong (and different from previous result without windowing).
ADC samples have DC offset, thus I eliminated the offset but still it doesnt work with windowed samples. The below is the first few output values through rs485.
First column is the fft output without window whereas second column is the output with window. From Column 1 the peak is at row 6 (6 x fs (10.5kHz) / 0.5N) gave me the correct input freq result where column 2 has a peak magnitude at row 2 (except dc bin) which does not make sense to me.
Any suggestion would be helpful.
Thanks in advance.
In the frequency domain, the rectangular and Hamming windows look like: As you may be aware, multiplication in the time domain by a window corresponds to a convolution in the frequency domain, which essentially spreads the energy of the signal over multiple frequency bins in what is often called spectral leakage. For the specific windows that you have chosen (depicted above in the frequency domain), you may notice that the Hamming window spreads much less energy outside the main lobe, but that main lobe is a little bit wider than that of the rectangular window. As a result, the spike of DC energy ends up spreading past bin 0, and into bin 1 when using the Hamming window. It isn't so much that you have a strong peak in bin 1. In fact, if you plot the data you provided, you should see that the spike you were seeing at index 6 is in fact still there at the same frequency with and without the use of the Hamming window: If you want to get rid of the DC spike and leakage in surrounding bins, either remove the bias in your data (essentially applying a notch filter), or you'll have to ignore a few more low-frequency bins when looking for your "first strong spike". Finally, note that there is also an issue with the way the IIR filter is implemented whereby the arrays If on the other hand you want to filter multiple blocks of Then you can use an offset in those buffer. Index 0 & 1 represent past samples, whereas the rest of the buffer from index 2 are filled with the current input data block. This leads to the following slightly modified filtering code: Finally, at the end of each block you then have to update the "past samples" at index 0 & 1, with the last samples of the current block to be ready to process the next input block:
这篇关于带有或不带有窗口的KISS FFT输出的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！过滤问题
ctr==0和ctr==1时，数组x和y的索引将超出范围(除非您已经实际上，某些特殊规定和x& y是与分配的缓冲区的起始位置偏移的指针.无论有无窗口，这都可能影响结果.如果仅过滤单个数据块，则通常的假设是较早的样本为零.在这种情况下，您可以使用以下方法避免超出范围的索引编制:// filter calculation
y[ctr] = num_coef[0]*x[ctr];
if (ctr>=1)
{
y[ctr] += (num_coef[1]*x[ctr-1]) - (den_coef[1]*y[ctr-1]);
}
if (ctr>=2)
{
y[ctr] += (num_coef[2]*x[ctr-2]) - (den_coef[2]*y[ctr-2]);
}
n个样本的多个块，则必须记住前一个块的最后几个样本.这可以通过分配比块大一点的缓冲区来完成:x = malloc((n+2)*sizeof(kiss_fft_scalar));
y = malloc((n+2)*sizeof(kiss_fft_scalar));
// initialize "past samples" for the first block, assuming data was zero
x[0] = x[1] = 0;
y[0] = y[1] = 0;
// filter calculation
y[ctr+2] = num_coef[0]*x[ctr+2];
y[ctr+2] += (num_coef[1]*x[ctr+1]) - (den_coef[1]*y[ctr+1]);
y[ctr+2] += (num_coef[2]*x[ctr]) - (den_coef[2]*y[ctr]);
// hamming window
hamming[ctr] = (0.54-((0.46) * cos(2*PI*ctr/256)));
window[ctr] = hamming[ctr]*y[ctr+2];
// remember last 2 samples of block
x[0] = x[n-2];
x[1] = x[n-1];
y[0] = y[n-2];
y[1] = y[n-1];
488 260 //dc bin
5 97
5 41
5 29
4 26
10 35
133 76
33 28
21 6
17 3
kiss_fft_scalar zero;
memset(&zero,0,sizeof(zero));
kiss_fft_cpx fft_input[n];
kiss_fft_cpx fft_output[n];
for(ctr=0; ctr<n; ctr++)
{
fft_input[ctr].r = zero;
fft_input[ctr].i = zero;
fft_output[ctr].r =zero;
fft_output[ctr].i = zero;
}
// IIR filter calculation
for (ctr=0; ctr<n; ctr++)
{
// filter calculation
y[ctr] = num_coef[0]*x[ctr];
y[ctr] += (num_coef[1]*x[ctr-1]) - (den_coef[1]*y[ctr-1]);
y[ctr] += (num_coef[2]*x[ctr-2]) - (den_coef[2]*y[ctr-2]);
//y1[ctr] += y[ctr] - 500;
// hamming window
hamming[ctr] = (0.54-((0.46) * cos(2*PI*ctr/256)));
window[ctr] = hamming[ctr]*y[ctr];
fft_input[ctr].r = window[ctr];
fft_input[ctr].i = 0;
fft_output[ctr].r = 0;
fft_output[ctr].i = 0;
}
kiss_fftr_cfg fftConfig = kiss_fftr_alloc(n,0,NULL,NULL);
kiss_fftr(fftConfig, (kiss_fft_scalar * )fft_input, fft_output);
for (ctr=0; ctr<n; ctr++)
{
fft_mag[ctr] = (sqrt((fft_output[ctr].r * fft_output[ctr].r) + (fft_output[ctr].i * fft_output[ctr].i)))/(0.5*n);
//Usart write
char filtResult[10];
sprintf(filtResult, "%04d %04d\n", (int)x[ctr], (int)fft_mag[ctr]);
//sprintf(filtResult, "%04d %04d\n", (int)x[ctr], (int)window[ctr]);
char c;
char *ptr = &filtResult[0];
do
{
c = *ptr;
ptr++;
usart_bw_write_char(&AVR32_USART2, (int)c);
// sendByte(c);
} while (c != '\n');
}
kiss_fft_cleanup();
free(fftConfig);
Frequency domain output clarifications
Filtering issue
x and y will be indexed out of bounds when ctr==0 and ctr==1 (unless you've made some special provision and x & y are in fact pointers with an offset from the start of the allocated buffers). This could affect the results for both with and without window. If you are only filtering a single block of data, a common assumption is that earlier samples were zeros. In that case, you can avoid the out-of-bound indexing with:// filter calculation
y[ctr] = num_coef[0]*x[ctr];
if (ctr>=1)
{
y[ctr] += (num_coef[1]*x[ctr-1]) - (den_coef[1]*y[ctr-1]);
}
if (ctr>=2)
{
y[ctr] += (num_coef[2]*x[ctr-2]) - (den_coef[2]*y[ctr-2]);
}
n samples, you'll have to remember the last few samples from the previous block. This could be done by allocating buffers that are slightly larger than the block size:x = malloc((n+2)*sizeof(kiss_fft_scalar));
y = malloc((n+2)*sizeof(kiss_fft_scalar));
// initialize "past samples" for the first block, assuming data was zero
x[0] = x[1] = 0;
y[0] = y[1] = 0;
// filter calculation
y[ctr+2] = num_coef[0]*x[ctr+2];
y[ctr+2] += (num_coef[1]*x[ctr+1]) - (den_coef[1]*y[ctr+1]);
y[ctr+2] += (num_coef[2]*x[ctr]) - (den_coef[2]*y[ctr]);
// hamming window
hamming[ctr] = (0.54-((0.46) * cos(2*PI*ctr/256)));
window[ctr] = hamming[ctr]*y[ctr+2];
// remember last 2 samples of block
x[0] = x[n-2];
x[1] = x[n-1];
y[0] = y[n-2];
y[1] = y[n-1];
