基于数字粒子图像处理的细水雾雾场速度测量(5)
单纯把“中心”相关计算改为“偏心”相关计算只是减少了一定相关窗口的操作，并没有减少计算量。而对于同一位置的Fw和Gw的所有点的相关计算可以根据卷积分定理和傅里叶变换（FFT）的相关定理得到：
(2.8)
其中，*为取共轭，F(U,V)是f(m,n)的傅里叶变换，G(U,V)是g(m,n)的傅里叶变换，只需对 作傅里叶逆变换，得到相关函数 (m,n)，然后检测 (m,n)的峰值，即可得到对应粒子的位移。
基于FFT的互相关算法的计算量可以分为两个方面计算：
1） 一幅速度矢量场图像需要算出的位移s(i,j)的数目为： (2.9)
2） 因为快速傅里叶变换的运算量为O(Nlog2N)次。所以用式(2.8)获得计算一次位移所需的相关函数的计算量为：
(2.10)
这样，总的计算次数为C=Cs×CF，上节所述的图像所需的计算次数为O(108) ，又降低了两个数量级。所以采用基于FFT的互相关算法比普通算法快很多，所需的限制条件只是流体的最大速度在横向和纵向都不能超过wx/3和wy/3，这一点在实际操作中是可以满足的。因此，基于以上原因，作者选择此算法完成DPIV系统处理软件，以完成细水雾雾场的测量。
2．5 亚像素拟合
如果计算仅到此为止，则所获得的数据均为整数值，也即此时算法对小于一个像素的位移不响应。因此系统的位移的动态范围限制在1~wx(或wy)像素这一狭窄的范围内，而且得到的速度值皆为整数，量化效应明显，这种精度往往在实际应用中是不够的。为了获得“亚像素”精度，有必要对结果进一步拟合，得到最大相关峰值在两像素间的精确位置。考虑到示踪粒子成像时亮度的分布为高斯分布：
(2.11)
所以本文将采用高斯曲线拟合峰值。设峰值所在点的x轴上的坐标为xp，相关值为1，则其横向相邻两像素在x轴上的坐标分别为xp-1和xp+1，对应点的相关值为0、2，则在横向拟合，可得极值点的位置为：
(2.12)
对最大相关峰所在像素点进行横向、纵向、斜向（此时分母中的2改为 ）四个方向进行拟合，最终确定最大峰值的精确位置。
2．6 数据后处理-误差矢量剔除
由于此算法中流体的最大速度在横向和纵向都不能超过wx/3和wy/3，所以在测量流场速度时，不仅要注意窗口大小的选取，还需在测量之后将超出限制范围的矢量剔除。
2．7 算法的计算机实现
本文采用矩阵预算语言MATLAB7.7.0作为软件开发和运行平台实现了上述算法，结合2.3和2.4两节的内容，FFT的互相关算法的步骤如图2.3(其中t为时间间隔):
图2.3 基于FFT的相关算法的步骤
下面介绍一下具体的运算过程（图2.4）：首先输入需要进行互相关计算的连续两帧图像，然后选择两幅图的查询计算窗口大小和移动步长，接着对两幅图进行FFT快速傅里叶变换并得到结果，然后对结果在频率进行复数共轭相乘计算，对结果进行反傅里叶变化，求出相关函数值，并将零点平移到谐线中心，就是将FFT变换得到的相关谱结果的图像调整一下，因为FFT变换的结果是将低频分量放在了4个角上，而高频分量放在了正中间，在matlab中采用fftshift函数即可。我们对FFT的结果进行寻找，查找出模最大的点，该点与中心点的距离就是位移量。其中傅里叶变换及反变换可以通过matlab中的fft2函数和ifft2获得，共轭可通过conj函数完成。接下来用高斯曲线拟合精确确定相关峰的位置及位移，然后重复进行相关运算，直至全幅图像计算完毕，根据时间间隔获得流场速度，画出速度矢量图。
上一篇：普通油料挥发性研究+文献综述