使用vDSP函数进行自相关

Question

给定一个一维的浮点或双精度向量，如何使用加速框架中的vDSP库函数计算该向量的自相关性？虽然vDSP_acor()和vDSP_acorD()函数被认为可以执行此计算，但是vDSP_Library.pdf文档（可在这里找到）并没有很好地解释如何使用函数参数。同样，vDSP_conv()和vDSP_convD()函数提到了可以在两个向量之间执行相关性和卷积的能力，但是它们没有提供足够的解释或示例代码，以便我能够成功使用它们。例如，如果使用滤波器核对二维矩阵进行卷积，则我想可能需要调用两次vDSP_convD()，并且signalStride的值不同，但是文档中省略了这一点。另一个省略的内容是滤波器中的数据必须如何打包。如果用零填充，那么零是否要先出现，后出现，或者它们需要均匀分布在非零条目的两侧？是否有滤波器长度、结果长度和输入长度的要求？

我希望您能提供一些有用的示例：使用vDSP_acor()和vDSP_conv()实现向量自身的自相关，二进制乘法是两个在频率域中作为实数据打包并进行前向FFT处理的数组，在IFT返回未归一化答案之前将用于计算自相关函数。在1D和2D数组上实现高斯核卷积。总的来说，这是一个很棒的库（可以说非常快！），但我发现这些特定的函数有点难以理解，上述示例可能会被广泛使用，因为它们在信号处理和图像分析中非常常见。

对于vDSP_Library参考文档的维护者，我的建议是：我假设“空间域”和“时间域”在整个文档中是等效的。如果不是，请确保进行区分。此外，请检查任何公式是否具有与所讨论的函数中声明的参数名称匹配的明确定义的参数。

注：这里所指的自相关性定义如下：A[T] = <(X[t]-m)(X[t-T]-m)>/v，其中A[T]是T时滞的自相关性，t是信号X的索引，m是X在所有t上的平均值，v是X在所有t上的方差，尖括号<>表示对所有间隔为T的X对进行平均。

- oldmanjank

Answer 1

< p>对我来说，vDSP_acorD的文档似乎还算清晰，但略显简洁，并且参数描述中似乎有一个错别字。

void vDSP_acorD (double * A,
    double * C, int N, int M);

A是输入信号，C是自相关输出，N是A中的样本数，M是您在C中需要的输出值的数量（即C将保存从0到M-1的滞后值）。

如果不可用vDSP_acorD，则可以尝试使用vDSP_conv，因为卷积是与一个输入信号反转的相关操作。

或者，您可以使用自己的快速自相关，因为自相关等效于功率谱的反向FFT，因此：

auto_correlation = IFFT(MAG(FFT(x)))

在哪里

FFT = forward FFT
IFFT = inverse FFT
MAG = magnitude of complex spectrum (sqrt(re * re + im * im))

然后使用来自vDSP的FFT例程。

- Paul R

---

谢谢，保罗。你对vDSP_acorD的解释似乎是合理的。对于卷积，似乎默认情况下它实际上可能会执行相关操作，并且您必须指定步幅为-1并将滤波器的最后一个元素传递给函数才能实际执行卷积。但是，卷积文档中还有其他一些使其在一般情况下更难使用（例如，用于执行高斯核的卷积）的事情。 - oldmanjank

哎呀！这可真让人头疼：vDSP_acor和vDSP_acorD在去年发布的vecLib-268.0版本中实际上没有在任何vecLib头文件中声明。 - oldmanjank

哦，但在这里有用的是，而且远远超出职责范围的是，一个自相关函数，它可以操作从调用vDSP_fft2d_zrip()得到的棘手打包数据，这是2D实时DFT。因为1D应用两次（水平然后垂直），但奈奎斯特频率存储在第二个元素中，所以进行逐元素乘法并不容易，特别是当数据结构是DSPSplitComplex时！ - oldmanjank

@oldmanjank：抱歉，我的意思是使用vDSP中的FFT例程来滚动快速自相关例程-请参见上文。 - Paul R

@oldmanjank：啊-你想要2D自相关-我从问题中没有猜到，但可能应该有-这时候生活会变得有点棘手。由于vDSP_acor例程似乎只支持1D，所以无论如何都不能帮助你。 - Paul R