我试图计算使用下面的代码时间序列样本窗口自相关。 我在申请FFT到窗口,然后计算的实部和虚部的幅度和虚部设定为零,最后采取逆它变换,以获得自相关:
DoubleFFT_1D fft = new DoubleFFT_1D(magCnt);
fft.realForward(magFFT);
magFFT[0] = (magFFT[0] * magFFT[0]);
for (int i = 1; i < (magCnt - (magCnt%2)) / 2; i++) {
magFFT[2*i] = magFFT[2*i] * magFFT[2*i] + magFFT[2*i + 1] * magFFT[2*i + 1];
magFFT[2*i + 1] = 0.0;
}
if (magCnt % 2 == 0) {
magFFT[1] = (magFFT[1] * magFFT[1]);
} else {
magFFT[magCnt/2] = (magFFT[magCnt-1] * magFFT[magCnt-1] + magFFT[1] * magFFT[1]);
}
autocorr = new double[magCnt];
System.arraycopy(magFFT, 0, autocorr, 0, magCnt);
DoubleFFT_1D ifft = new DoubleFFT_1D(magCnt);
ifft.realInverse(autocorr, false);
for (int i = 1; i < autocorr.length; i++)
autocorr[i] /= autocorr[0];
autocorr[0] = 1.0;
第一个问题是:可以看出,这个代码的自相关结果映射到[0,1]
范围内,但是相关性被认为是在-1和1。当然,很容易将结果映射到之间[-1,1]
范围内,但我不知道,如果这个映射是正确的。 我们如何解释导致价值观autocorr
阵列?
其次,与此代码我歌厅良好的效果对于一些周期性的系列,这是我根据信号的周期获取特定的自相关指数值越高。 然而,结果去奇怪,当我把它应用到非周期信号:在所有的值autocorr
阵列看上去非常接近1是什么原因呢?
对于基于FFT算法的工作,你必须特别注意的定义,包括你正在使用的库的约定。 你似乎混淆交流的“信号处理”公约和“统计”之一。 再有就是FFT包装和零填充。
下面是的工作为偶数n情况下,代码,信号处理惯例。 它是针对蛮力包裹自相关测试。 该言论表明如何将其转化为信号处理惯例。 对于统计交流,数据的平均值减去。 这可以仅仅通过归零FFT的“为0Hz”组件来完成。 然后将交流的第零元素是变化的,你可以通过除以这个数量正常化。 所得的值将落在-1..1就像你说的。
您的代码似乎是在做分割通过,但不能忽视的数据0 Hz分量。 因此,它的计算某种约定的混搭。
import edu.emory.mathcs.jtransforms.fft.DoubleFFT_1D;
import java.util.Arrays;
public class TestFFT {
void print(String msg, double [] x) {
System.out.println(msg);
for (double d : x) System.out.println(d);
}
/**
* This is a "wrapped" signal processing-style autocorrelation.
* For "true" autocorrelation, the data must be zero padded.
*/
public void bruteForceAutoCorrelation(double [] x, double [] ac) {
Arrays.fill(ac, 0);
int n = x.length;
for (int j = 0; j < n; j++) {
for (int i = 0; i < n; i++) {
ac[j] += x[i] * x[(n + i - j) % n];
}
}
}
private double sqr(double x) {
return x * x;
}
public void fftAutoCorrelation(double [] x, double [] ac) {
int n = x.length;
// Assumes n is even.
DoubleFFT_1D fft = new DoubleFFT_1D(n);
fft.realForward(x);
ac[0] = sqr(x[0]);
// ac[0] = 0; // For statistical convention, zero out the mean
ac[1] = sqr(x[1]);
for (int i = 2; i < n; i += 2) {
ac[i] = sqr(x[i]) + sqr(x[i+1]);
ac[i+1] = 0;
}
DoubleFFT_1D ifft = new DoubleFFT_1D(n);
ifft.realInverse(ac, true);
// For statistical convention, normalize by dividing through with variance
//for (int i = 1; i < n; i++)
// ac[i] /= ac[0];
//ac[0] = 1;
}
void test() {
double [] data = { 1, -81, 2, -15, 8, 2, -9, 0};
double [] ac1 = new double [data.length];
double [] ac2 = new double [data.length];
bruteForceAutoCorrelation(data, ac1);
fftAutoCorrelation(data, ac2);
print("bf", ac1);
print("fft", ac2);
double err = 0;
for (int i = 0; i < ac1.length; i++)
err += sqr(ac1[i] - ac2[i]);
System.out.println("err = " + err);
}
public static void main(String[] args) {
new TestFFT().test();
}
}