CN103645341B - 全流场3d可视化测速方法 - Google Patents

Abstract

本发明公开了一种微流场的全流场可视化测速方法，包括以下步骤：基于频域光学相干层析技术建立微粒子跟踪全流场可视化测速系统；对流体进行二维扫描，连续采集流体的干涉谱数据；基于傅立叶变换方法获得流体三维图像；利用局部灰度阈值和体积滤波方法搜索流场中各个微粒子，采用平方加权质心法求各个微粒子的三维坐标，实现微粒子的可视化；通过定义代价函数对微粒子进行匹配；利用微粒子的三维坐标变化求运动速度矢量。通过对微粒子成像和跟踪实现流场中流速的可视化；通过跟踪全流场微粒子的运动轨迹实现全流场三维速度矢量测量；具有微米级的空间分辨率，特别适合于复杂微流场的三维速度矢量检测。

Description

全流场3D可视化测速方法

技术领域

本发明涉及一种全流场可视化测速方法。尤其是涉及微流场的全流场可视化测速方法。

背景技术

流体流速是流场最基本的物理量，人们对流体流动特性的认识较大程度上取决于速度场的获取。到目前为止，还无法从理论上彻底解决湍流等问题，只能借助科学实验手段建立或改进流体力学模型，并对实际问题进行更深入的研究。微流动器件在生命科学、化学、微机械加工等领域应用也更加广泛，微流动的特性更加复杂，科学实验的作用越来越突出。测量流体中流速的方法有机械方法、散热率法、动力测压法、激光多普勒测速和粒子成像测速等。其中机械方法测量流速的原理是根据流体中叶轮的旋转角速度与流体流速成正比；散热率法测量流速是基于发热的测速传感器的散热率与流体流速成比例的原理；毕托管(PitotTube)为典型的动力测压法，基于流体力学Bernoulli能量方程进行平均速度测量和流量测量；以上流速测量为单点、接触式测量方法，对流场的影响也较大。激光多普勒测速原理是利用流体中运动微粒散射光的多普勒频移获得流体的速度信息，采用多光束测量系统可获得三维流速。激光多普勒测速方法有较高的空间分辨率，但为单点测速技术，无法测量湍流等复杂流场。粒子成像测速(PIV)原理是在流体中散布跟随性好的固体粒子(或流体中已经存在)作为示踪粒子，通过观察流场中微粒子图像，对连续多个瞬时速度场的相关性分析，直观反映全场瞬时流动信息，可以深入研究复杂流场瞬时过程，具有单点测速技术无可比拟的优点。粒子图像处理算法主要有相关算法、最小二次差方算法和粒子跟踪算法等。Brody等[″BiotechnologyatlowReynoldsnumbers"：BiophysicalJournal，3430-3441，1996]使用超荧光显微镜对直径900nm的荧光示踪粒子进行长时间照明和观察，获得粒子运动轨迹图像，这种方法获得的速度场分辨率、精度都很低。SANTIAGO等[″Aparticleimagevelocimetrysystemformicrofluidics"，ExperimentsinFluids，316-319，1998]使用汞弧灯照明对直径300nm的荧光示踪粒子进行照明，用CCD采集示踪粒子图像，利用相关算法得到流场速度，这种方法获得的空间分辨率可以得到10um。二维测速的空间分辨率已经达到100nm[″Micro-ParticleImageVelocimetry(μPIV)：Recentdevelopments，applications，andguidelines"，LabonaChip，2551-2567，2009]，二维平面流场测速技术已经相对成熟。全流场测量可以通过多二维平面扫描、数字全息技术、立体粒子图像测速、散焦数字图像测速和共焦荧光显微技术等。但一般结构复杂，而且全流场信息需要通过物镜焦平面的移动获得，一般只适合于定常流动和周期性流动[″显微粒子图像测速技术——微流场可视化测速技术及应用综述″，机械工程学报，155-168，2012]。

发明内容

为了克服背景技术的不足，本发明的目的是提供一种全流场可视化测速方法。本发明包括以下步骤：

(1)基于光学低相干原理，建立频域光学相干层析微粒子跟踪全流场可视化测速系统；

(2)利用激光束对渗入微粒子或本身存在微粒子的流体进行二维扫描，连续实时采集流体的干涉谱数据；

(3)利用每个扫描位置的干涉谱数据到空间域转换算法，得到每个扫描位置在深度方向的结构图像，一组二维扫描可以获得流体的三维结构图像；

(4)利用各个微粒子图像具有较深灰度原理，并根据图像块的体积大小搜索流场中各个微粒子三维图像，实现全流场微粒子的可视化；

(5)利用平方加权质心法求取各个微粒子的中心三维坐标；

(6)通过定义代价函数并使其最小化方法对连续扫描得到的两个流体体积中的微粒子进行匹配；

(7)根据每一对已匹配的微粒子的三维坐标值求所对应微粒子的运动方向和位移，即流体中此处的速度矢量，对流场中所有的已匹配的微粒子进行同样的处理即可得到全流场的速度矢量；

(8)对连续采样的相邻三维流体图像中的微粒子进行匹配，跟踪流场中微粒子的运动轨迹，检测全流场中速度的实时变化情况。

本发明由于采取以上技术方案，其具有以下优点：

通过对流体中的微粒子成像和跟踪过程实现了流场中流速的可视化；利用流体中全流场分布的微粒子实现全流场速度检测；通过跟踪微粒子的运动轨迹实现三维速度矢量测量；具有微米级的空间分辨率，特别适合于复杂微流场的三维速度矢量检测。

附图说明

图1是根据本发明实施例的基于频域光学相干层析技术的微粒子跟踪全流场可视化测速系统示意图；

图2是根据本发明实施例的全流场可视化测速系统数据采集与处理流程图；

图3是根据本发明实施例的典型的流场中微粒子可视化二维图像；

图4是根据本发明实施例的典型的流场中微粒子可视化三维图像；

图5是根据本发明实施例的三维速度矢量检测原理示意图；

图6是根据本发明实施例的对流流体二维切面速度矢量示意图；

具体实施方式

下面结合附图详细说明本发明的实施例。

图1是根据本发明实施例的基于频域光学相干层析技术的微粒子跟踪全流场可视化测速系统示意图。

参照图1，超辐射发光二极管(SLD)(101)宽带光源发出的弱相干光进入2×2光纤耦合器(102)后分为两束光，一束光进入参考臂，并经过准直透镜(103)和平面镜(104)返回，另外一束光进入样品臂，经由准直透镜(105)、二维振镜系统(Y扫描仪106和X扫描仪107)、物镜(108)入射到被测流体上(114)。两束光经原路返回后在光纤耦合器内发生干涉，其干涉谱由准直透镜(109)、衍射光栅(110)、消色差透镜(111)和CMOS线阵摄像机(112)组成的光谱仪(115)接收，并通过与线阵摄像机数据接口相对应的接口板将数据读入计算机(113)。物镜的数值孔径决定了测量系统空间的横向分辨率：其中R_l为横向分辨率，λ₀为光源中心波长，NA_obj为物镜数值孔径。测量系统空间的纵向分辨率取决于光源的带宽，其中R_a为纵向分辨率，λ₀为光源中心波长，Δλ为光源带宽。

图2是根据本发明实施例的全流场可视化测速系统数据采集与处理流程图；

在步骤201中，通过二维振镜系统(106和107)将聚焦光束在渗入微粒子的流体的被测范围内(amm*bmm)进行横向快速扫描，获取被测流场的干涉谱数据。可包括以下步骤：

(1)X扫描仪(107)在锯齿波的驱动下将聚焦光束进行X方向的扫描，一个锯齿波周期完成一个B扫描(快速扫描)，在一个B扫描内(amm)进行M次采样，每次采样即可得到此扫描位置对应的干涉谱数据；

(2)Y扫描仪(106)在斜波信号的驱动下使聚焦光束在Y方向扫描(慢速扫描)，将Y方向的扫描长度(bmm)分为N等份，即设置N个B扫描位置，完成一次被测流体的二维扫描共得到M*N个扫描位置的干涉谱数据；

(3)重复以上步骤(1)和(2)K次，连续采集被测流体的K组干涉谱数据。

在步骤202中，通过干涉谱域到空间域转换关系得到每一个扫描位置的沿光轴方向结构信息，将M*N个扫描位置的结构信息组合起来即可得到流体三维结构图像。典型的干涉谱域到空间域转换关系为：

W(z)＝IFFT(I(λ))

其中I(λ)为干涉谱强度随波长变化函数，W(z)反映流体沿光轴方向结构信息，IFFT为逆傅立叶变换。此处可以采用均匀或非均匀傅立叶逆变换。由于加入流体中微粒子的直径在1微米以下，肉眼根本无法看到，其尺寸也小于光学系统的分辨率。但在光学系统传递函数的作用下，微粒子形成的图像大小要大于其实际尺寸，便于观察，实现了流体中各种复杂情况的可视化。典型的流场中微粒子可视化二维图像和三维图像分别如图3和图4所示。

在步骤203中，根据灰度值搜索流体中微粒子，并求取每一个微粒子中心的三维坐标。可包括以下步骤：

(1)为减小各种随机噪声对可靠提取微粒子图像的影响，采用中值滤波方法对原始三维图像进行处理；

(2)利用最大类间方差原理设定灰度阈值，若灰度值大于阈值，认为是微粒子图像部分，对流体三维图像进行分割；

(3)流体三维图像中除了大小比较均匀的微粒子图像外，还有更小的杂质形成的图像，以及微粒子粘连图像，为去除它们的影响设置体积最大阈值V_max和最小阈值V_min,若分割后小块体积大于V_max，则认为出现了微粒子粘连情况，若分割后小块体积小于V_min,则认为是更微小粒子形成的图像，予以去除，微粒子形成的图像大小V_particle由光学系统的点扩展函数决定，一般可以选择V_max＝1.5V_particle，V_min＝0.5V_particle；

(4)利用上面的图像二值化过程可以初步搜索出各个微粒子图像所在位置。为求取微粒子运动速度，需要从二值化前的三维灰度图像中精确计算各个微粒子图像中心的三维坐标，因此，采用以下的平方加权质心法求取微粒子中心三维坐标。

$\left\{\begin{array}{c}{x}_0=\frac{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}{F}^2(x,y,z)x}{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}{F}^2(x,y,z)},\\ y_0=\frac{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}{F}^2(x,y,z)y}{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}{F}^2(x,y,z)},\\ z_0=\frac{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}{F}^2(x,y,z)z}{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}{F}^2(x,y,z)}\end{array}\right.$

其中，F(x，y，z)为灰度值，m，n，p为考虑的微粒子图像窗口范围。

在步骤204中，通过定义代价函数并使其最小化方法对连续扫描得到的两个流体体积中的微粒子进行匹配，可包括以下步骤：

(1)微粒子匹配就是从连续扫描得到的两个体积图像中搜寻来自同一个微粒子的图像。假定第一个三维图像中搜寻到K₁个微粒子p_i，i=1，2，…，K₁，第二个三维图像中搜寻到K₂个微粒子q_j，j＝1，2，…，K₂，p_i和q_j的关系用矩阵r_ij表示，如果p_i和q_j为同一个微粒子，则r_ij＝1，否则r_ij＝0。为了得到最优关系r_ij，定义了代价函数并使其最小化。定义的代价函数为：

$\mathrm{\ψ}=\underset{i=1}{\overset{K_1}{\mathrm{\Σ}}}\underset{j=1}{\overset{K_2}{\mathrm{\Σ}}}{\mathrm{\φ}}_{\mathrm{ij}}{r}_{\mathrm{ij}}$

其中，φ_ij表示将p_i和q_j关联的代价，定义为p_i和q_j之间距离的平方与二阶矩之差的平方之和，即

${\mathrm{\φ}}_{\mathrm{ij}}={(x_{p_i}-x_{q_j})}^2+{(y_{p_i}-y_{q_j})}^2+{(z_{p_i}-z_{q_j})}^2+{({\mathrm{moment}}_{p_i}-{\mathrm{moment}}_{q_j})}^2$

其中，moment为二阶矩，表示为

$\mathrm{moment}=\frac{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}[{(x-x_0)}^2+{(y-y_0)}^2+{(z-z_0)}^2]F(x,y,z)}{\underset{x=1}{\overset{m}{\mathrm{\Σ}}}\underset{y=1}{\overset{n}{\mathrm{\Σ}}}\underset{z=1}{\overset{p}{\mathrm{\Σ}}}F(x,y,z)}$

其中，F(x，y，z)为灰度值，m，n，p为考虑的微粒子图像窗口范围，x₀，y₀，z₀为微粒子中心三维坐标。

(2)求连续采集到的两个三维流体中微粒子中心三维坐标(和)并对微粒子进行匹配后，用以下公式计算微粒子三维运动速度矢量。

$v_x=\frac{x_{q_j}-x_{p_i}}{T},v_y=\frac{y_{q_j}-y_{p_i}}{T},v_z=\frac{z_{q_j}-z_{p_i}}{T}$

T为采样时间间隔。其原理如图5所示，图6所示为纵向二维切面上的速度场分布，其中箭头方向表示在此二维平面上流速方向，长度表示流速大小。

在步骤205中，对连续采样的相邻三维流体图像中的微粒子进行匹配，跟踪流场中微粒子的运动轨迹，检测全流场中速度的实时变化情况。

Claims (5)

1.一种全流场可视化测速方法，其特征在于包括以下步骤：

(1)基于频域光学相干层析原理，建立微粒子跟踪全流场可视化测速系统；

(2)利用激光束对渗入微粒子或本身存在微粒子的流体进行二维扫描，连续实时采集流体的干涉谱数据；

(3)基于频域光学相干层析成像原理将每个扫描位置的干涉谱数据转换到空间域，得到每个扫描位置在深度方向的结构图像，一组二维扫描可以获得流体的三维结构图像；

(4)根据各个微粒子图像具有较深灰度特点，并利用体积滤波方法搜索流场中各个微粒子三维图像，利用平方加权质心法求取各个微粒子的中心三维坐标，实现全流场微粒子的可视化；

(5)通过定义代价函数并使其最小化方法对连续扫描得到的两个流体体积中的微粒子进行匹配；根据每一对已匹配的微粒子的三维坐标值求所对应微粒子的运动方向和位移，即流体中此处的速度矢量，对流场中所有的已匹配的微粒子进行同样的处理即可得到全流场的速度矢量；

(6)对连续采样的相邻三维流体图像中的微粒子进行匹配，跟踪流场中微粒子的运动轨迹，检测全流场中速度的实时变化情况。

2.根据权利要求1所述的全流场可视化测速方法，其中，将每个扫描位置的干涉谱数据转换到空间域步骤包括均匀傅立叶逆变换或非均匀傅立叶逆变换方法。

3.根据权利要求1所述的全流场可视化测速方法，其中，全流场微粒子的可视化步骤包括：

(1)采用中值滤波方法减小各种随机噪声的影响；

(2)利用最大类间方差方法设定灰度阈值，对微粒子图像进行分割；

(3)为去除更小的杂质以及微粒子粘连的影响，设置体积最大阈值V_max和最小阈值V_min，若分割后小块体积大于V_max，则认为出现了微粒子粘连情况，若分割后小块体积小于V_min，则认为是更微小粒子形成的图像，予以去除，微粒子形成的图像大小V_particle由光学系统的点扩展函数决定，选择V_max＝1.5V_particle，V_min＝0.5V_particle；

(4)利用上面的图像二值化过程可以初步搜索出各个微粒子图像所在位置，为求取微粒子运动速度，需要从二值化前的三维灰度图像中精确计算各个微粒子图像中心的三维坐标，因此，采用平方加权质心法求取微粒子中心三维坐标。

4.根据权利要求1所述的全流场可视化测速方法，其特征在于，通过将微粒子之间距离的平方与二阶矩之差的平方之和定义为代价函数，通过对其优化，实现微粒子之间的可靠匹配。

5.根据权利要求1所述的全流场可视化测速方法，其特征在于，对连续采样的相邻三维流体图像中的微粒子进行匹配，跟踪流场中微粒子的运动轨迹，检测全流场中速度矢量的实时变化情况。