CN112036072A - 一种三维示踪粒子匹配方法及速度场测量装置 - Google Patents

一种三维示踪粒子匹配方法及速度场测量装置 Download PDF

Info

Publication number
CN112036072A
CN112036072A CN202010663981.3A CN202010663981A CN112036072A CN 112036072 A CN112036072 A CN 112036072A CN 202010663981 A CN202010663981 A CN 202010663981A CN 112036072 A CN112036072 A CN 112036072A
Authority
CN
China
Prior art keywords
particle
objective function
matching
particles
dimensional space
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010663981.3A
Other languages
English (en)
Other versions
CN112036072B (zh
Inventor
潘翀
聂铭远
王晋军
王隆威
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beihang University
Original Assignee
Beihang University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beihang University filed Critical Beihang University
Priority to CN202010663981.3A priority Critical patent/CN112036072B/zh
Publication of CN112036072A publication Critical patent/CN112036072A/zh
Application granted granted Critical
Publication of CN112036072B publication Critical patent/CN112036072B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/25—Design optimisation, verification or simulation using particle-based methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/004—Artificial life, i.e. computing arrangements simulating life
    • G06N3/006—Artificial life, i.e. computing arrangements simulating life based on simulated virtual individual or collective life forms, e.g. social simulations or particle swarm optimisation [PSO]

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Evolutionary Computation (AREA)
  • General Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Computer Hardware Design (AREA)
  • Computing Systems (AREA)
  • Software Systems (AREA)
  • Mathematical Physics (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Molecular Biology (AREA)
  • General Health & Medical Sciences (AREA)
  • Biophysics (AREA)
  • Computational Linguistics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Medical Informatics (AREA)
  • Data Mining & Analysis (AREA)
  • Biomedical Technology (AREA)
  • Algebra (AREA)
  • Fluid Mechanics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Indicating Or Recording The Presence, Absence, Or Direction Of Movement (AREA)
  • Image Analysis (AREA)

Abstract

本申请提供一种三维示踪粒子匹配方法及速度场测量装置，该方法包括：获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标；根据至少一个第一粒子的三维空间坐标和至少一个第二粒子的三维空间坐标，生成全局目标函数，所述全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成，根据全局目标函数，得到第一粒子和第二粒子之间的匹配结果，该匹配结果用于反映流体速度场，从而实现了长比位移情况下对流动速度场的高精度测量。该方法的有效的提高了示踪粒子匹配结果的精度；同时，提高了全局目标函数的收敛速度，显著提升了计算效率。

Description

一种三维示踪粒子匹配方法及速度场测量装置
技术领域
本申请涉及流体力学的测量技术领域，尤其涉及一种三维示踪粒子匹配方法及速度场测量装置。
背景技术
粒子追踪测速技术(Particle Tracking Velocimetry，简称PTV)是一种流动速度场光学测量技术，该技术追踪由激光照明、由数字相机拍摄得到的示踪粒子图像，将多台相机拍摄得到的平面粒子图像还原为粒子的三维空间坐标，得到单个粒子在相邻两张照片帧之间的跨帧位移，实现对示踪粒子所表征的流场运动速度的测量。由于该技术对单一示踪粒子的运动进行追踪，因此速度场测量的空间分辨率高达像素级。该技术的一个核心问题是如何确定相邻两帧图像上的示踪粒子之间的匹配关系。
现有技术中，可以使用最小位移法或粒子群相似度法来确定相邻两帧图像上的示踪粒子之间的匹配关系。其中，最小位移法是将后一帧中与前一帧的目标粒子间隔最近的备选粒子作为匹配粒子，粒子群相似度法是将前一帧目标粒子与临近粒子的空间位置关系作为匹配准则在后一帧中进行查找匹配。
然而，以上两种方法在长比位移情况下，均会出现较大的匹配误差，导致示踪粒子匹配结果的精度较低。
发明内容
本申请提供一种三维示踪粒子匹配方法及速度场测量装置，用于解决现有技术中示踪粒子匹配精度较低的问题，从而完成长比位移流动的高精度速度场测量。
第一方面，本申请提供一种示踪粒子匹配检测方法，包括：
获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标，所述第一粒子和所述第二粒子分别为相邻两帧图像中的示踪粒子；根据所述至少一个第一粒子的三维空间坐标和所述至少一个第二粒子的三维空间坐标，生成全局目标函数，所述全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成；根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果。
在一种可选的实施方式中，所述最小位移法函数为无量纲化处理后的函数。
在一种可选的实施方式中，所述根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果，包括：确定所述全局目标函数的全局最优解；将所述全局最优解对应的匹配结果作为所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果。
在一种可选的实施方式中，所述确定所述全局目标函数的全局最优解，包括：使用匹配蚁群算法求解所述全局目标函数，得到所述全局目标函数的全局最优解，所述匹配蚁群算法在迭代过程中，匹配目标从后一帧粒子随机跳回前一帧粒子。
在一种可选的实施方式中，所述使用目标蚁群算法求解所述全局目标函数，得到所述全局目标函数的全局最优解，包括：根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，确定所述第一粒子和第二粒子组成的粒子对的信息素和启发素；A、根据所述第一粒子和第二粒子组成的粒子对的信息素和启发素，计算第一粒子和第二粒子组成的粒子对的匹配概率；B、根据所述第一粒子和所述第二粒子的匹配概率，得到所述第一粒子和所述第二粒子的一个粒子对的匹配结果；C、将匹配目标随机跳回前一帧图像中未配对过的粒子，将所述未配对过的粒子作为新的第一粒子；迭代执行步骤B-C，直至得到所有所述第一粒子与所述第二粒子的匹配结果；D、根据所述信息素、所有所述第一粒子与所述第二粒子的匹配结果，更新所述信息素，得到新的信息素；迭代执行步骤A-D，直至达到预设的迭代次数或预设的匹配精度，得到所述全局目标函数的全局最优解。
在一种可选的实施方式中，所述根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，确定所述第一粒子和第二粒子组成的粒子对的信息素和启发素，包括：根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，计算所述第一粒子和第二粒子组成的粒子对的启发素；根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数、信息素的遗留速率和蒸发速率，计算所述信息素。所述信息素的遗留速率为归一化处理后的遗留速率。
在一种可选的实施方式中，根据泰森多面体边界面极径形式的相关函数得到所述粒子群相似度法函数。
第二方面，本申请实施例提供一种使用第一方面所述任一项方法的三维速度场测量装置，包括：
获取模块，用于获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标，所述第一粒子和所述第二粒子分别为相邻两帧图像中的示踪粒子。
处理模块，用于根据所述至少一个第一粒子的三维空间坐标和所述至少一个第二粒子的三维空间坐标，生成全局目标函数，所述全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成；以及，根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果，该匹配结果用于反映流体速度场。
在一种可选的实现方式中，所述最小位移法函数为无量纲化处理后的函数。
在一种可选的实现方式中，所述处理模块具体用于：
确定所述全局目标函数的全局最优解；将所述全局最优解对应的匹配结果作为所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果。
在一种可选的实现方式中，所述处理模块具体用于：
使用目标蚁群算法求解所述全局目标函数，得到所述全局目标函数的全局最优解，所述目标蚁群算法在迭代过程中，匹配目标从后一帧粒子随机跳回前一帧粒子。
在一种可选的实现方式中，所述处理模块具体用于：
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，确定所述第一粒子和第二粒子组成的粒子对的信息素和启发素；
A、根据所述第一粒子和第二粒子组成的粒子对的信息素和启发素，计算第一粒子和第二粒子组成的粒子对的匹配概率；
B、根据所述第一粒子和所述第二粒子的匹配概率，得到所述第一粒子和所述第二粒子的一个粒子对的匹配结果；
C、将匹配目标随机跳回前一帧图像中未配对过的粒子，将所述未配对过的粒子作为新的第一粒子；
迭代执行步骤B-C，直至得到所有所述第一粒子与所述第二粒子的匹配结果；
D、根据所述信息素、所有所述第一粒子与所述第二粒子的匹配结果，更新所述信息素，得到新的信息素；
迭代执行步骤A-D，直至达到预设的迭代次数或预设的匹配精度，得到所述全局目标函数的全局最优解。
在一种可选的实现方式中，所述处理模块具体用于：
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，计算所述启发素；
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标、所述全局目标函数、信息素的遗留速率和蒸发速率，计算所述信息素；
归一化处理所述信息素的遗留速率。
在一种可选的实现方式中，所述处理模块具体用于：
根据泰森多面体边界面极径形式的相关函数得到所述粒子群相似度法函数。
本申请实施例提供的一种三维示踪粒子匹配方法及速度场测量装置，获取流场中相邻两帧示踪粒子的三维空间坐标，可以由最小位移法函数和粒子群相似度法函数混合生成全局目标函数，进而基于全局目标函数得到粒子的匹配结果，该匹配结果进一步反映了流动速度场。由于全局目标函数由最小位移法函数和粒子群相似度法函数以预设权重混合得到，因此，全局目标函数能够同时表征最小位移法函数和粒子群相似度法函数的特征，通过最小位移法的辅助可以限制粒子群相似度法中粒子群的范围，从而保证粒子群相似法在全局下的随机性得到抑制，极大的提高了示踪粒子匹配的精度；同时，通过引入最小位移法，还能够使得跨帧相邻粒子匹配的其他因子得到聚集，从而提高了全局目标函数的收敛速度，显著提升了计算效率。
附图说明
为了更清楚地说明本申请或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
图1为本申请实施例提供的一种基于三维PTV测速系统的硬件实现示意图；
图2为本申请实施例提供的一种示踪粒子匹配方法的流程示意图；
图3为本申请实施例提供的一种泰森多面体的示意图；
图4为本申请实施例提供的一种生成全局目标函数的流程示意图；
图5为本申请实施例提供的一种生成基于最小位移法的目标函数的流程示意图；
图6为本申请实施例提供的一种生成粒子群相似度法函数和基于粒子群相似度法的目标函数的流程示意图；
图7为本申请实施例提供的一种匹配目标随机跳回前一帧粒子的示意图；
图8为本申请实施例提供的一种求解全局目标函数最优解的流程示意图；
图9为本申请实施例提供的一种三种不同比位移情况的示意图；
图10为本申请实施例提供的一种在三种不同比位移情况下本申请全局目标函数、最小位移法函数、粒子群相似度法函数在多次迭代运行中的匹配精度曲线示意图；
图11为本申请实施例提供的一种在两种不同比位移情况下参数与匹配精度关系的示意图；
图12为本申请实施例提供的一种速度场测量装置的模块结构图；
图13为本申请实施例提供的一种使用改进蚁群算法处理全局目标函数的完整流程示意图。
具体实施方式
为使本申请的目的、技术方案和优点更加清楚，下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本申请保护的范围。
为了便于理解，首先，对本申请所涉及的概念进行说明。
1、粒子追踪测速技术(Particle Tracking Velocimetry，简称PTV)
是指通过追踪单个粒子的运动轨迹计算单个粒子的速度，属于拉格朗日方法，即把质点作为研究对象，着眼于流体质点，设法描述出每一个流体质点自始至终的运动过程，因此能够实现较高的空间精度。粒子追踪测速技术的硬件设备一般包括：流体通道、激光光源、示踪粒子、图像采集设备、图像处理设备和通信设备等。
2、比位移
是指粒子在两帧图像之间的平均跨帧位移大小除以单帧图像上粒子的平均间距的比值，用以表征流动速度与粒子浓度的相对大小关系。其中，在流体特性基本平稳的情况下，上述单帧图像的粒子平均间距既可以通过前一帧图像，也可以通过后一帧图像以计算。在本申请中，当比位移值处于[0.5,1.5)的区间时，定义其为长比位移。在本申请中，短比位移也可以称为短位移，中比位移也可以称为中位移，长比位移也可以称为长位移。
3、跨帧位移
跨帧位移属于矢量概念，通过同一坐标系下至少两帧图像获得，其大小为某个粒子在前一帧图像中与该粒子在后一帧图像中的三维坐标间的距离，其方向为该粒子在前一帧图像中的三维坐标点指向其在后一帧图像中的三维坐标点。
4、粒子平均间距
是指在一帧图像中，每一个粒子的质心到其他粒子质心的三维距离的均值，用以表征流体中粒子的浓度。
5、跨帧时间
是指相机在连续拍摄过程中，获取相邻两帧图像的时间间隔。
6、对极几何
是计算机视觉理论的几何学基础，用于描述同一景物两幅图像之间的视觉几何关系。
7、轮盘算法
是一种依据概率大小对不同对象进行选择的计算机算法。例如：存在事件1、事件2和事件3，三者互相独立且三者概率和为1，其中事件1发生概率p1＝0.2，事件2发生概率p2＝0.3，事件3发生概率p3＝0.5。若计算机程序要依据概率对此三个事件进行选择，则计算机先随机生成一个处于(0,1]区间的数字，若该随机数处于(0,0.2]区间，则对应着选择了事件1；若该随机数处于(0.2,0.5]区间，则对应事件2；若该随机数处于(0.5,1]区间，则对应事件3。
在粒子追踪测速技术中，主要包括粒子识别过程和粒子匹配过程。现有技术中，在粒子匹配过程上，可以使用最小位移法和粒子群相似度法对不同图像帧之间的粒子进行匹配，即对不同图像帧之间的同一个粒子进行匹配确认。其中，最小位移法属于纯位移测度的方法，粒子群相似度法属于纯粒子群相似测度的方法。
在实际粒子测速的应用场景中，首先，对于最小位移法，其仅就粒子之间的距离对粒子进行匹配，在长比位移情况下，粒子跨帧位移大小超过平均粒子间距，使得可能出现该粒子在前一帧与其相邻粒子在后一帧的三维距离，较该粒子在前一帧与后一帧的跨帧距离更为接近，因此，造成难以区分跨帧相邻粒子之间匹配关系，导致最小位移法的匹配精度显著降低。其次，对于粒子群相似度法，在长比位移情况下，粒子群的变形程度增加，使得同一粒子群在前一帧与后一帧的相似程度降低，且搜索半径的增大会使得参与候选的粒子群增多，因此，难以区分搜索半径内粒子之间匹配关系，导致粒子群相似度法的匹配精度显著降低。
如图10所示，在短比位移、中比位移和长比位移情况下不同方法多次迭代处理后的精度结果，能够体现最小位移法与粒子群相似度法的特性：在短比位移、中比位移情况下，粒子群相似度法在区分粒子的匹配关系上非常有效，而最小位移法的收敛速度较快。
考虑到最小位移法和粒子群相似度法的上述特性，以及在长比位移情况下二者匹配精度不高的问题，本申请实施例通过混合最小位移法和粒子群相似度法构造全局目标函数并使用特定算法求解该全局目标函数的全局最优解，从而得到高精度的粒子之间的匹配结果。
本申请实施例的方法及装置可以应用在水洞、风洞或其仿真环境等场景中。以下以水洞仿真环境这一场景为例进行说明。
图1为本申请实施例方法在水洞或其仿真环境下的硬件测量系统示意图，如图1所示：在多功能水洞1中均匀混入荧光示踪粒子2，使用Nd-YAG双腔激光器3和透镜组发出双脉冲激光照明待测试区域；使用四台跨帧CCD相机5，从不同角度对被双曝光激光照亮的荧光粒子进行跨帧拍照，相机和激光器的同步工作由专用的同步器6控制，每台相机可得至少一对跨帧图像，并传入计算机7存储。其中，水洞1用以产生均匀且稳定的水流，荧光粒子2在测量区域上游均匀散布于水中；Nd-YAG双腔激光器3用以产生跨帧时间满足需求的、波长为532nm的两次激光脉冲，并经过配套的扩束透镜组产生具有一定厚度的体激光，通过折射镜4照亮待测区域，其方向平行于水洞前观测壁；四台CCD相机5在镜头前可以安装有高通滤光片，以截断532nm的照射激光、只让由荧光粒子受激发出的更大波长的荧光信息被CCD相机记录；四台CCD相机以不同角度同时对观测区域进行拍照，在每一曝光时刻产生四张荧光强度照片，并采样成数字形式；用一个专用的同步控制器连接中控计算机、激光器和四台相机，以确保激光器和相机可以在同一时刻协同工作；相机采集的照片由数据线传输到中控计算机7中，处理过程或结果由显示器8呈现。
在水洞、风洞或其仿真环境等场景下，可以通过本申请实施例的速度场测量装置，由四台跨帧CCD相机采集粒子在不同角度的图像，基于这些图像能够确定前一帧图像和后一帧图像中粒子质心的三维空间坐标。在确定出前一帧图像和后一帧图像中粒子质心的三维空间坐标后，首先基于这些粒子质心的三维空间坐标生成最小位移法函数和粒子群相似度法函数，再对上述两种函数按照预设的权重进行混合，得到本申请的全局目标函数，最后通过特定算法求解该混合生成的全局目标函数的全局最优解，得到前一帧图像与后一帧图像中粒子之间的匹配关系，并基于该匹配关系得到前后两帧图像相对应粒子的跨帧位移。进一步，用该跨帧位移除以跨帧时间得到示踪粒子的速度，从而表征该流动的速度场。
图2为本申请实施例提供的一种流动速度场测量方法的流程示意图，如图2所示，该方法包括以下步骤：
S201、获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标，该第一粒子和该第二粒子分别为相邻两帧图像中的粒子。
可选的，上述第一粒子和第二粒子，可以分别是水洞或其仿真环境中相邻两帧图像的示踪粒子，也可以分别是风洞或其仿真环境中相邻两帧图像的示踪粒子。
本申请中，粒子的三维空间坐标可以指粒子质心的三维空间坐标。示例性的，如图1所示，可以将四台相机在同一时刻、以不同视角测量获得的二维平面粒子图片组，通过使用标准PTV粒子识别算法并配合空间标定重构出三维粒子云。其中，标准PTV粒子识别算法将传输到计算机的数字化粒子图像使用基于高斯分布模版的标准粒子识别算法进行处理，获得每一个粒子的质心。对于四台相机中的任一台而言，均可以识别出同一个粒子的质心二维平面坐标。空间标定应用了对极几何理论，即任意两台相机之间拍摄的粒子图像点与相机空间位置点和被拍摄粒子三维空间位置点存在五点共面约束，从而能确定在其余相机中相同粒子的投影。当在全部四台相机的粒子图像上均能识别出同一个粒子时，则可通过空间三角定位法确定该粒子质心的三维空间坐标。
可选的，上述三维空间坐标，可以是在空间直角坐标系中的坐标，可以是在球坐标系中的坐标，可以是在柱坐标系中的坐标，也可以是其他空间坐标系中的坐标；相应的，上述二维平面坐标，可以是在平面直角坐标系中的坐标，可以是平面极坐标系中的坐标，也可以是其他平面坐标系中的坐标，本申请对此不作具体限定。
S202、根据至少一个第一粒子的三维空间坐标和至少一个第二粒子的三维空间坐标，生成全局目标函数，全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成。
本申请实施例中，首先通过至少一个第一粒子的三维空间坐标和至少一个第二粒子的三维空间坐标，生成最小位移法函数和粒子群相似度法函数。其中，生成最小位移法函数的方法，包括但不限于以下任一项：四帧最小加速度(4Minimal Acceleration，简称4MA)PTV算法、四帧最优估计(4Best Estimation，简称4BE)算法。生成粒子群相似度法函数的方法，包括但不限于以下任一项：松弛算法、基于三角嵌铺的粒子追踪测速算法(DelaunayTessellation Particle Tracking Velocimetry，简称DT-PTV)、泰森图粒子追踪测速算法(Voronoi Diagram Particle Tracking Velocimetry，简称VD-PTV)。
全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成。示例性的，通过至少一个第一粒子的三维空间坐标和至少一个第二粒子的三维空间坐标生成最小位移法函数和粒子群相似度法函数，再对最小位移法函数和粒子群相似度法函数以预设权重混合，可以得到上述全局目标函数。
S203、根据全局目标函数，得到至少一个第一粒子和至少一个第二粒子之间的匹配结果，该匹配结果进一步反映了流动速度场。
本申请实施例中，基于匹配关系得到前后两帧图像相对应粒子的跨帧位移，并用该跨帧位移除以跨帧时间得到示踪粒子的速度，从而表征该流动的速度场。
由于全局目标函数基于最小位移法函数和粒子群相似度法函数按照预设权重混合生成，因此，该全局目标函数能够同时表征最小位移法函数的特征以及粒子群相似度法函数的特征。
本申请实施例中，通过全局目标函数得到至少一个第一粒子和至少一个第二粒子之间的匹配结果，可以是以特定方法求解全局目标函数的全局最优解。可选的，粒子之间的匹配结果可以指粒子之间的匹配关系。求解全局目标函数的全局最优解的方法包括但不限于以下方法中一种或多种方法的结合，同样也包括对其改进后的方法：蚁群算法、模拟退火算法、粒子群算法、遗传算法、禁忌搜索算法、爬山算法。
本实施例中，基于相邻两帧中粒子的三维空间坐标，可以由最小位移法函数和粒子群相似度法函数混合生成全局目标函数，进而基于全局目标函数得到粒子的匹配结果。由于全局目标函数由最小位移法函数和粒子群相似度法函数以预设权重混合得到，因此，全局目标函数能够同时表征最小位移法函数和粒子群相似度法函数的特征，通过最小位移法的辅助可以限制粒子群相似度法中粒子群的范围，从而保证粒子群相似度法在全局下的随机性得到抑制，极大的提高了粒子匹配的精度；同时，通过引入最小位移法，还能够使得跨帧相邻粒子匹配的其他因子得到聚集，从而提高了全局目标函数的收敛速度，显著提升了计算效率。
以下说明上述步骤S202中基于最小位移法函数和粒子群相似度法函数按照预设权重混合生成全局目标函数的一种可选方式。
在该方式中，将最小位移法函数进行无量纲化处理后，再与粒子群相似度法函数按照预设权重进行混合，以生成全局目标函数。
图4为生成全局目标函数的流程示意图，如图4所示，生成全局目标函数的一种可选方式包括如下步骤：
S301、基于前一帧中的任意粒子i和后一帧中的任意粒子j的质心的三维空间坐标，生成基于最小位移法函数。
步骤1、通过下述公式(1)，生成最小位移法函数DFij，两两配对计算前一帧粒子质心的三维空间坐标(xi,yi,zi)到后一帧粒子质心三维空间坐标(xj,yj,zj)的位移的二范数。二范数又名Euclid范数，即向量元素绝对值的平方和再开方，是一种计算向量长度的方法。
Figure BDA0002579667650000101
在执行上述步骤1之后，还可以执行下述步骤2，得到基于最小位移法的目标函数min(DF)。
步骤2、还可以执行本步骤，得到基于最小位移法的目标函数min(DF)，以用于后续的分析过程。当存在一个匹配关系使得在关系内的DFij总和最小，则求解粒子间匹配关系的问题可转化为对应的求解全局最优解的问题，故纯基于最小位移法生成的目标函数可定义为公式(2)。
Figure BDA0002579667650000102
其中，(i,j)∈match表示该目标函数的某一个解确立的相邻两帧中一组粒子的匹配关系，β1是在后述生成全局目标函数DPF过程中，调节最小位移法函数和粒子群相似度法函数之间相对权重的加权系数之一。
S302、基于前一帧中的任意粒子i和后一帧中的任意粒子j的质心的三维空间坐标，计算得到粒子群相似度法函数PFij。
示例性的，测量一个粒子与其邻域之间相对空间位置的方法可以是泰森图法，它将离散粒子所占据的空间分割为一组三维泰森多面体或二维泰森多边形，每个多面体(或多边形)由连接两个相邻粒子的线段的垂直平分面形成，并且在每个泰森多面体的中心仅包含一个粒子。从泰森多面体内的所有点到中心粒子质心的距离总是比到其他任何粒子质心的距离都短，可以通过下述公式(3)所表达。
Figure BDA0002579667650000103
其中，P是图像中待分割的粒子，Voronoi(pi)是指以粒子pi的质心为中心的泰森多面体，dist(x,pi)表示Voronoi(pi)内的点x到粒子质心坐标pi的距离。对于给定的粒子分布，泰森图是确定的；因此，每个泰森多面体的几何结构可用于描述相应的中心粒子的分布模式。
计算得到粒子群相似度法函数PFij的过程可以包括如下7个步骤：
步骤1、基于相邻两帧图像中粒子质心的三维空间坐标，生成泰森多面体。
可选的，生成泰森多面体的方式，可以使用图像处理软件或数据分析软件中的已有代码库生成泰森多面体，或者，也可以使用自定义的代码库生成泰森多面体。
示例性的，本实施例中粒子群相关性函数的构造遵循与在二维VD-PTV算法中相同的思想，即两个粒子的泰森多边形之间的相似性是通过多边形边缘球坐标系函数ρ(θ)的互相关值来测量的。
步骤2、求解三维泰森多面体其边界平面S在球坐标系下的形式。
可以通过如下公式(4)求解。
ρ＝ρS(θ,φ) (4)
其中，极径ρ是中心粒子到边界面S上一点(x,y,z)沿(θ,φ)角方向的距离，如图3所示的一种泰森多面体的示意图，以下述公式(5)的形式表示。
Figure BDA0002579667650000111
步骤3、将边界平面S转换为空间直角坐标系下的形式，可以下述公式(6)的形式表示。
Ax+By+Cz+D＝0 (6)
其中，A,B,C,D为可以由边界面S上的顶点坐标确定的几何系数。
将公式(4)至公式(6)进行组合，可得到衡量泰森多面体边界面的球坐标函数公式(7)：
Figure BDA0002579667650000112
步骤4、通过上述公式(7)，对每一个粒子的每个上述离散化后的角度自变量(θ,φ)进行求解，得到每个粒子对应的至少一个极径值。
步骤5、判断每一个极径所在射线与泰森多面体边界平面的交点是否落在其有限的边界面里。
其中，极径所在射线指的是由某粒子质心出发，沿着各个离散化的角度(θ,φ)发出的射线，且当该射线与该粒子质心所对应的泰森多面体边界平面不平行时，两者之间最多只有一个交点。
其中，有限的边界面指的是，对于某一个泰森多面体的每一个边界平面来说，其数学公式表达的是一个数学意义上的平面，具有无限性。但泰森多面体的每一个边界平面都是有限的，都属于一个平面上的具体的多边形，其有限的边界面能够以该多边形在平面内所有顶点的三维坐标等一些要素进行限定或勾勒。
步骤6、基于步骤1-步骤5所得到的信息对前一帧中的粒子和后一帧中的粒子，两两计算泰森多面体边界面极径函数的互相关性性。
计算前一帧的粒子i和后一帧的粒子j的泰森多面体边界面极径函数的互相关函数Rij，可以下述公式(8)进行计算。
Figure BDA0002579667650000121
上述公式(8)中的参数含义与前述步骤中的参数含义相同，不再赘述。
步骤7、生成粒子群相似度法函数PFij。
可以使用下述公式(9)生成粒子群相似度法函数。
Figure BDA0002579667650000122
在上述公式(9)中，粒子群相似度法函数PFij被定义为Rij的倒数，由于相关系数Rij是一个处于数值区间[0,1]的量，因此粒子群相似度法函数PFij处于[1,+∞)的数值区间。
在上述步骤7之后，还可以执行下述步骤8，得到基于粒子群相似度的目标函数，以用于后续的分析过程。
步骤8、生成基于粒子群相似度的目标函数min(PF)。
与前述最小位移法函数DFij的处理方法类似，如果存在一个匹配关系使得在关系内的PFij总和最小，则求解匹配关系的问题可转化为对应求解最优解的问题，故基于粒子群相似度的目标函数可定义为公式(10)：
Figure BDA0002579667650000123
其中，(i,j)∈match表示该目标函数的某一个解确立的相邻两帧中一组粒子的匹配关系，β2是在后述生成全局目标函数DPF过程中，调节最小位移法函数和粒子群相似度法函数之间相对权重的加权系数之一。
S303，对最小位移法函数进行无量纲化处理，得到无量纲化处理后的函数DF。
由于PF的量纲为1，而DF的量纲为L，其中L用于表征长度，在混合二者之前，需要对DF进行无量纲化处理，可以下述公式(11)对其进行处理。
Figure BDA0002579667650000124
其中,分母是所有可能路径(p,q)的最小距离，p和q分别指代第一帧和第二帧的所有粒子。在混合二者时，用上述公式(11)取代前述公式(3)，这样使得DFij≥1。结合PFij≥1的条件，DPFij的最小值为1，且属于[1,+∞)的数值区间范围。
S304，以预设权重混合对最小位移法函数DFij无量纲化处理后的函数和粒子群相似度法函数，得到全局目标函数min(DPF)。
可以通过如下公式(12)和公式(13)计算全局目标函数min(DPF)。
Figure BDA0002579667650000131
其中，β1和β2是调节最小位移法函数DF和粒子群相似度法函数PF之间相对权重的加权系数。且当β1＝2且β2＝0时，本算法退化为一种最小位移法函数，当β1＝0且β2＝6时退化为一种粒子群相似度法函数。
Figure BDA0002579667650000132
其中，(i,j)∈match表示该全局目标函数的某一个解确立的相邻两帧中一组粒子的匹配关系。
本实施例中，由于PF的量纲为1，而DF的量纲为L，而函数在量纲不同时是不能直接混合的。通过公式(11)对DF函数进行无量纲化处理，使得PF和DF的量纲得到统一，因此，两者能够通过公式(13)进行混合，生成全局目标函数，即上述公式(13)。同时，本申请通过公式(7)和公式(8)，在业界首次提出了使用泰森多面体边界面的极径形式相关性函数，以刻画相邻粒子群的相似度，并得到本申请中的粒子群相似度法函数PF。因此，本申请中的粒子群相似度法函数PF能够在DF函数的辅助下，生成结合两种方法优点的全局目标函数DPF，使粒子匹配精度得到了极大提升。
以下说明上述步骤S202中根据全局目标函数，得到至少一个第一粒子和至少一个第二粒子之间的匹配结果的一种可选方式，即确定全局目标函数的全局最优解，将全局最优解对应的匹配结果作为该至少一个第一粒子和至少一个第二粒子之间的匹配结果。
可选的，可以首先确定全局目标函数的全局最优解，进而，将全局最优解对应的匹配结果作为至少一个第一粒子和至少一个第二粒子之间的匹配结果。
其中，求解某些函数的全局最优解方法一般包括：蚁群算法、遗传算法、粒子群算法、爬山算法、模拟退火算法、禁忌搜索算法。
本申请中，在确定全局目标函数的全局最优解时，可以基于匹配蚁群算法来确定，其中，匹配蚁群算法是对传统的蚁群算法进行改进所得到的一种蚁群算法。在匹配蚁群算法的迭代过程中，考虑到了蚂蚁或匹配目标从位于后一帧的粒子随机跳回前一帧粒子的步骤。
其中，传统的蚁群算法可以适用于求解传统的旅行商问题(Travelling SalesmanProblem，简称TSP)函数的全局最优解，其中，旅行商问题可以理解为代表了一类能够转化为某些排列组合问题的具体问题集合。但在本申请中，粒子匹配问题由于客观存在着蚂蚁或匹配目标从后一帧粒子随机跳回前一帧粒子的步骤，因此不能直接使用传统的蚁群算法求解得到其全局目标函数的全局最优解。因此，本申请对传统的蚁群算法进行了改进，克服了前述的问题，使得匹配蚁群算法能够求解得到本申请中全局目标函数的全局最优解。
以下对使用匹配蚁群算法确定全局目标函数的全局最优解的过程进行说明。
可选的，可以使用匹配蚁群算法求解全局目标函数，得到全局目标函数的全局最优解，匹配蚁群算法在迭代过程中，蚂蚁或称为匹配目标从后一帧粒子随机跳回前一帧粒子。
其中，使用匹配蚁群算法求解全局目标函数，得到全局目标函数的全局最优解的方式可以是通过匹配蚁群算法在通过信息素和启发素正向反馈的基础上多次迭代，逐步逼近并最终生成S304中的全局目标函数min(DPF)的全局最优解。
其中，图7是匹配蚁群算法中蚂蚁或称为匹配目标从后一帧粒子随机跳回前一帧粒子的示意图，匹配蚁群算法在迭代过程中，蚂蚁或称为匹配目标从后一帧粒子随机跳回前一帧粒子，可以是在蚁群算法中，若前一帧中的有N个待匹配粒子i，在第k次迭代过程中，当前一帧中的粒子i按照本次迭代生成的的匹配概率与后一帧中的粒子j完成了一次匹配，此处不同于传统的蚁群算法，蚂蚁或匹配目标会随机跳回前一帧图像中，从以前未配对过的粒子中随机选择一个作为新的粒子i，与后一帧中未匹配过的粒子继续进行匹配处理。
本实施例中，由于传统的蚁群算法不能适用于求解本申请中全局目标函数的全局最优解，因此，通过对蚁群算法进行改进，使得匹配蚁群算法在迭代过程中，位于后一帧粒子的蚂蚁随机跳回前一帧粒子，因而能够实现对全局目标函数全局最优解的求解。
以下对使用匹配蚁群算法求解全局目标函数得到全局目标函数的全局最优解的具体过程进行说明。
图8为匹配蚁群算法求解全局目标函数得到全局目标函数的全局最优解的过程示意图，如图8所示，使用匹配蚁群算法求解全局目标函数得到全局目标函数的全局最优解的过程包括：
首先，根据第一粒子的三维空间坐标、第二粒子的三维空间坐标和上述的全局目标函数，确定第一粒子和第二粒子组成的粒子对的信息素和启发素。
可选的，可以根据第一粒子的三维空间坐标、第二粒子的三维空间坐标和上述全局目标函数，计算第一粒子和第二粒子组成的粒子对的启发素，其中启发素一旦得到确定，在多次迭代过程中并不发生改变。并且，根据第一粒子的三维空间坐标、第二粒子的三维空间坐标、上述的全局目标函数、信息素的初始值、信息素的遗留速率、信息素的蒸发速率和信息素的更新公式，计算信息素。具体执行过程可以参照下述的公式(14)和公式(15)。
进而，执行下述步骤。
S401、根据第一粒子和第二粒子组成的粒子对的信息素和启发素，计算第一粒子和第二粒子组成的粒子对的匹配概率；
S402、根据一个第一粒子和第二粒子的归一化处理后的匹配概率，使用轮盘算法，得到一个粒子对的匹配结果。
S403、匹配目标随机跳回前一帧图像中一个未配对过的粒子，将该粒子作为新的第一粒子，与后一帧图像中若干个未配对过的粒子进行配对。
迭代执行步骤S402-S403，直至得到所有所述第一粒子与所述第二粒子的匹配结果；。
S404、根据信息素、所有第一粒子与第二粒子的匹配结果，更新信息素，得到新的信息素；
迭代执行步骤S401-S404，直至达到预设的迭代次数或达到预设的全局目标函数收敛条件，此时即可得到全局目标函数的全局最优解。
以下基于具体的公式说明上述步骤S401-S404的第k次迭代中的实施方式。
首先，基于前述的公式得到DFij、PFij和DPFij。
其次，基于DFij、PFij和DPFij，计算得到粒子i和粒子j组成的粒子对的启发素和信息素。
启发素ηij是匹配概率中与路径本身信息相关的部分，启发素不会随着迭代次数而改变，反映了粒子之间的几何学性质，用于表征不同的粒子对(i,j)或不同选择路径(i,j)的收益，代表着蚂蚁的个体智慧，启发素ηij的定义通过下述公式(14)表达。
Figure BDA0002579667650000151
其中，由于DFij和PFij均属于[1,+∞)的数值区间，因此，DPFij也属于[1,+∞)的数值区间，启发素ηij则属于[0,1]的数值区间范围。越大的DPFij，则对应的启发素ηij越小；反之亦然。
信息素
Figure BDA00025796676500001613
是匹配概率中由之前迭代过程中蚂蚁经验所贡献的部分，在每一次迭代中，每只蚂蚁都会在其行进的路径上留下一小部分信息素，因此，路径的距离越小，信息素遗留的越多，这正代表了蚁群的集体智慧。其中，k表示当前处于第k次迭代过程中，且信息素的初始值为：
Figure BDA0002579667650000161
与启发素不同，信息素在每次迭代后都有所不同，加之考虑到在迭代中没有蚂蚁走过的道路的信息素会挥发，因此，信息素
Figure BDA0002579667650000162
可通过初始值、信息素增量公式(15)等以信息素更新公式(18)所计算。
Figure BDA0002579667650000163
其中，k表示当前处于第k次迭代过程。M是一个整数常数，即蚁群中蚂蚁的总数量；m表示第m只蚂蚁，且m＝1,2,...,M。
Figure BDA0002579667650000164
表示第m只蚂蚁在第k次迭代中的路径(i,j)。
Figure BDA0002579667650000165
为第k次迭代过程中第m只蚂蚁的总行程，如图7所示，总行程表示在N对粒子全部匹配完成后，对应的N条粗黑实线的长度之和，
Figure BDA0002579667650000166
可以通过下述公式(16)计算。Q为一常数，表示信息素的遗留速率，为了使得
Figure BDA0002579667650000167
在每次迭代后都能保持
Figure BDA0002579667650000168
公式(16)中的常数Q可以为通过公式(17)进行归一化处理后的Q；其中，M表示蚁群中蚂蚁的个数，min(DPF)表示该次迭代中最优化的DPF值。
Figure BDA0002579667650000169
Figure BDA00025796676500001610
在第k次迭代基础上的第k+1次迭代中信息素更新的计算方式，通过下述公式(18)进行。
Figure BDA00025796676500001611
其中，γ为一个常数，表示信息素的蒸发速率，此处γ＝0.5。
综合上述公式(12)、(13)、(15)、(16)、(17)和(18)，得到下述公式(19)。
Figure BDA00025796676500001612
其中，
Figure BDA0002579667650000171
表示本轮迭代结束后各个蚂蚁的总行程中的最小值；m'表示蚁群中的任意蚂蚁。
同时，生成公式(19)处于[0,1]的数值区间范围内，通过下述公式(20)可以计算得到第k次迭代中蚁群中的一只蚂蚁选择(i,j)的路径行进的概率
Figure BDA0002579667650000172
Figure BDA0002579667650000173
其中，
Figure BDA0002579667650000174
被定义为信息素和启发素的幂指数乘积，由加权参数α和β共同决定信息素和启发素对
Figure BDA0002579667650000175
的相对贡献。由于信息素
Figure BDA0002579667650000176
和启发素ηij都处于[0,1]的数值区间范围，因此，
Figure BDA0002579667650000177
也处于[0,1]的数值区间范围。由于启发素ηij公式(14)中已经存在权重系数β1和β2，此处设置β＝1。
由于在匹配蚁群算法迭代过程中，对于第k次迭代中位于前一帧中粒子i处的一个蚂蚁m，可能存在多条匹配概率偏大的路径，这表明前一帧中的粒子i，可能与多个后一帧中的粒子j存在较大的匹配概率
Figure BDA0002579667650000178
所以有可能使得总的概率之和大于1，因此，需要对匹配概率
Figure BDA0002579667650000179
进行归一化处理，通过下述公式(21)生成归一化处理后的匹配概率
Figure BDA00025796676500001710
Figure BDA00025796676500001711
其中，
Figure BDA00025796676500001712
表示在第k次迭代过程中，对于第m只蚂蚁，在该蚂蚁已经匹配了一部分粒子的情况下，第二帧中所有未配对过的粒子的集合；l表示对该集合中的任意粒子的遍历。
图14是本申请实施例提供的一种使用改进蚁群算法处理全局目标函数的完整流程示意图，通过上述公式(21)，对于第k次迭代中位于前一帧中粒子i处的一个蚂蚁m，对位于后一帧中多个可被接受的j，得到归一化处理后的匹配概率
Figure BDA00025796676500001713
即前一帧中粒子i与若干个后一帧的粒子j匹配概率之和
Figure BDA00025796676500001714
为1，然后使用轮盘算法基于概率大小选择路径(i,j)，得到一对粒子的匹配结果。
然后，蚂蚁或匹配目标随机跳回前一帧图像中一个未配对过的粒子，将该粒子作为新的第一粒子，继续与后一帧图像中若干个未配对过的粒子进行配对。直到得到第一粒子与第二粒子所有的匹配结果。
最后，判断是否达到预设的迭代次数，或，匹配结果达到预设的匹配精度，若是，则流程结束，以本轮迭代后所有蚂蚁的总路程中最小值作为全局最优解，该全局最优解即对应着一组最优的所有的第一粒子与第二粒子的匹配结果。若否，则进入步骤S404，更新信息素，并进入下一轮迭代过程。
本实施例中，通过公式(17)，对信息素的遗留速率Q进行归一化处理，能够防止信息素更新公式(18)
Figure BDA0002579667650000181
中的右侧，在当信息素
Figure BDA0002579667650000182
与信息素增量
Figure BDA0002579667650000183
出现较大的数值差异时(例如数量级式的差异时)，部分信息被忽略因而造成的迭代过程中信息丢失问题。
以下说明本申请的方案在具体应用场景中的实施过程。该应用场景为水洞仿真环境。在该应用场景中，以标准三维粒子图像测速技术(Particle Image Velocimetry，简称PIV)流场的VSJ-PIV(the Visualization Society of Japan，简称VSJ)数据库中的#351作为测试流场，该测试流场采用大涡模拟(Large Eddy Simulation，简称LES)方法求解了射流雷诺数为3000的二维平面射流冲击壁面流场中的三维粒子位移。因此，图像帧中的示踪粒子质心的三维坐标以及其正确的匹配关系已经包含在VSJ-PIV数据库中。在下文中，#351的第0帧和第1帧被归为短比位移测试的图像对，而同一图像集的第0帧和第2帧被用于中比位移测试，第0帧和第4帧被用于长比位移测试。测试图像平均每帧包含2000个粒子。为了评估算法的不确定性，随机选取1000个无粒子丢失的跨帧粒子对子集作为输入。图9表示了三种情况下(短比位移、中比位移和长比位移)粒子对的一个子集。
其中，精度衡量可以用如下公式表示：
Figure BDA0002579667650000184
其中，分子Ncorrect是正确匹配的粒子对数量，分母N是输入粒子对的总数。对于这三种情况，在独立且随机选择的子集上运行30次，将匹配精度的平均值以及均方根(RootMean Square，简称RMS)值与其他算法获得的结果进行比较。
下述表一中示例性的给出本申请中可能使用的默认参数。
Figure BDA0002579667650000185
表一
图10显示了在三种不同比位移情况下本申请全局目标函数、最小位移法函数、粒子群相似度法函数在多次迭代运行中的匹配精度曲线。
图11显示了两种位移情况下采用了本申请的全局目标函数DPF(其中β1＝2且β2＝6)和将β2或β1设置为零，全局目标函数DPF退化为纯粹的最小位移法函数DF或纯粹的粒子群相似度函数PF在30次独立运行中的平均短比位移情况下的精度收敛曲线。所有算法都采用了匹配蚁群算法实现作为全局优化解决方案。结果表明，虽然随着比位移的增加，所有测试算法的匹配精度都会降低，但本申请中本申请全局目标函数不仅具有最高的粒子匹配精度，而且具有最快的函数收敛速度。
图11为短比位移、长比位移情况下的测试结果示意图，如图11所示，在本申请实施例中，短比位移测试中，三种测度方法均能达到95％以上精度，但本申请全局目标函数在第1轮迭代就能达到98％精度，最小位移法和粒子群相似度法在该测试中达到相同精度分别需要迭代10次和18次；但在长比位移的情况下，经过10次迭代后，本申请全局目标函数仍然具有约94％的匹配精度，相比之下，最小位移法的粒子匹配精度降低到约50％，粒子群相似度法的粒子匹配精度也降低到80％。因此，本申请全局目标函数在匹配精度、函数收敛速度和运算效率上都具备显著的提升或改善。
以下结合具体应用场景描述一种本申请推荐参数的选取方式。在本实施例中，以标准三维PIV流场VSJ-PIV数据库中的#351作为测试流场，与图9所示例的测试场景基本相同，采用大涡模拟方法求解了射流雷诺数为3000的二维平面射流冲击壁面流场中的三维粒子位移。因此，图像帧中的示踪粒子质心的三维坐标以及其正确的匹配关系已经包含在VSJ-PIV数据库中。在下文中，#351的第0帧和第1帧被归为短比位移测试的图像对，第0帧和第4帧被用于长比位移测试。
参照前述公式(14)和公式(20)，本申请中可以有三个自由参数，即α、β1和β2。权重α的影响已在先前蚁群算法的研究中得到详细说明，该参数是全局最小化的权重，大的α将减慢收敛速度，但可以避免陷入局部最优解。因此，本实施例只关注β1和β2。对于短比位移和长比位移情况，这里，算法中的所有参数，除了β1和β2，都采用实施例二中所示的默认值，匹配精度取30次独立运行的平均值。
继续参照图11，对于短比位移流场，本申请中匹配蚁群算法对β1和β2不敏感。只要β1和β2足够大，匹配精度就可以达到99％以上。然而，在长比位移流场中，匹配精度对β1比β2更敏感。对于固定的β2，增加β1将大大降低性能；同时，一旦β2＞2，匹配精度相对于β2的变化相对较弱。在本发明中，β1和β2分别调节距离测度和粒子群相似性测度相对于混合测度函数DPF的相对权重。本实施例表明，距离测度应在DPF中占据较少的权重，这与实施例二中的结果一致，即纯粒子群相似测度的性能优于纯位移测度。另外，如图11所示，β1比β2均不能为零；否则，匹配精度会严重下降。可以说明，需要将这两项措施作为一个整体结合起来，特别是在比位移较大的情况下。根据图11，β1的建议范围是[1,2]。对于β2，由于较大的β2抑制了蚁群迭代过程的随机性，从而增加了陷入局部最优解的概率，降低了鲁棒性。因此，β2建议在[3,8]的范围内。本实施例中，选择的默认值为β1＝2.0且β2＝6.0。
在通过上述实施例中的第一粒子与第二粒子的匹配结果之后，可以基于该匹配结果以及第一粒子和第二粒子在至少一对相邻帧中的高精度三维空间坐标，等信息，得到前后两帧图像相对应粒子的跨帧位移。进一步，用该跨帧位移除以跨帧时间得到示踪粒子的速度，表征该流动的速度场，从而完成了长比位移流动的高精度速度场测量。
图12为本申请实施例提供的速度场测量装置的模块结构图，如图12所示，该装置1200包括：
获取模块1201，用于获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标，所述第一粒子和所述第二粒子分别为相邻两帧图像中的示踪粒子。
处理模块1202，用于根据所述至少一个第一粒子的三维空间坐标和所述至少一个第二粒子的三维空间坐标，生成全局目标函数，所述全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成；以及，根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果，该匹配结果用于反映流体速度场。
在一种可选的实施方式中，所述最小位移法函数为无量纲化处理后的函数。
在一种可选的实施方式中，处理模块1202具体用于：
确定所述全局目标函数的全局最优解；将所述全局最优解对应的匹配结果作为所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果。
在一种可选的实施方式中，处理模块1202具体用于：
使用匹配蚁群算法求解所述全局目标函数，得到所述全局目标函数的全局最优解，所述匹配蚁群算法在迭代过程中，位于后一帧粒子的蚂蚁随机跳回前一帧粒子。
在一种可选的实施方式中，处理模块1202具体用于：
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述目标函数，确定所述第一粒子和第二粒子组成的粒子对的信息素和启发素；
A、根据所述第一粒子和第二粒子组成的粒子对的信息素和启发素，计算第一粒子和第二粒子组成的粒子对的匹配概率；
B、根据所述第一粒子和所述第二粒子的匹配概率，得到所述第一粒子和所述第二粒子的一个粒子对的匹配结果；
C、将匹配目标随机跳回前一帧图像中未配对过的粒子，将所述未配对过的粒子作为新的第一粒子；
迭代执行步骤B-C，直至得到所有所述第一粒子与所述第二粒子的匹配结果；
D、根据所述信息素、所有所述第一粒子与所述第二粒子的匹配结果，更新所述信息素，得到新的信息素；
迭代执行步骤A-D，直至达到预设的迭代次数或预设的匹配精度，得到所述全局目标函数的全局最优解。
在一种可选的实施方式中，处理模块1202具体用于：
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，计算所述第一粒子和第二粒子组成的粒子对的启发素。
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标、所述全局目标函数、信息素的遗留速率和蒸发速率，计算所述信息素。
本申请实施例提供的速度场测量装置，可以执行上述方法实施例中的方法步骤，其实现原理和技术效果类似，在此不再赘述。
需要说明的是，应理解以上装置的各个模块的划分仅仅是一种逻辑功能的划分，实际实现时可以全部或部分集成到一个物理实体上，也可以物理上分开。且这些模块可以全部以软件通过处理元件调用的形式实现；也可以全部以硬件的形式实现；还可以部分模块通过处理元件调用软件的形式实现，部分模块通过硬件的形式实现。例如，确定模块可以为单独设立的处理元件，也可以集成在上述装置的某一个芯片中实现，此外，也可以以程序代码的形式存储于上述装置的存储器中，由上述装置的某一个处理元件调用并执行以上确定模块的功能。其它模块的实现与之类似。此外这些模块全部或部分可以集成在一起，也可以独立实现。这里所述的处理元件可以是一种集成电路，具有信号的处理能力。在实现过程中，上述方法的各步骤或以上各个模块可以通过处理器元件中的硬件的集成逻辑电路或者软件形式的指令完成。
在本申请实施例中，“至少一个”是指一个或者多个，“多个”是指两个或两个以上。“和/或”，描述关联对象的关联关系，表示可以存在三种关系，例如，A和/或B，可以表示：单独存在A，同时存在A和B，单独存在B的情况，其中A，B可以是单数或者复数。字符“/”一般表示前后关联对象是一种“或”的关系；在公式中，字符“/”，表示前后关联对象是一种“相除”的关系。“以下至少一项(个)”或其类似表达，是指的这些项中的任意组合，包括单项(个)或复数项(个)的任意组合。例如，a，b，或c中的至少一项(个)，可以表示：a，b，c，a-b，a-c，b-c，或a-b-c，其中，a，b，c可以是单个，也可以是多个。
可以理解的是，在本申请实施例中涉及的各种数字编号仅为描述方便进行的区分，并不用来限制本申请实施例的范围。
可以理解的是，在本申请的实施例中，上述各过程的序号的大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本申请实施例的实施过程构成任何限定。
最后应说明的是：以上各实施例仅用以说明本申请的技术方案，而非对其限制；尽管参照前述各实施例对本申请进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本申请各实施例技术方案的范围。

Claims (8)

1.一种示踪粒子匹配方法，其特征在于，包括：
获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标，所述第一粒子和所述第二粒子分别为相邻两帧图像中的示踪粒子；
根据所述至少一个第一粒子的三维空间坐标和所述至少一个第二粒子的三维空间坐标，生成全局目标函数，所述全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成；
根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果。
2.根据权利要求1所述的方法，其特征在于，所述最小位移法函数为无量纲化处理后的函数。
3.根据权利要求1所述的方法，其特征在于，所述根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果，包括：
确定所述全局目标函数的全局最优解；
将所述全局最优解对应的匹配结果作为所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果。
4.根据权利要求3所述的方法，其特征在于，所述确定所述全局目标函数的全局最优解，包括：
使用匹配蚁群算法求解所述全局目标函数，得到所述全局目标函数的全局最优解，所述匹配蚁群算法在迭代过程中，匹配目标从后一帧粒子随机跳回前一帧粒子。
5.根据权利要求4所述的方法，其特征在于，所述使用匹配蚁群算法求解所述全局目标函数，得到所述全局目标函数的全局最优解，包括：
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，确定所述第一粒子和第二粒子组成的粒子对的信息素和启发素；
A、根据所述第一粒子和第二粒子组成的粒子对的信息素和启发素，计算第一粒子和第二粒子组成的粒子对的匹配概率；
B、根据所述第一粒子和所述第二粒子的匹配概率，得到所述第一粒子和所述第二粒子的一个粒子对的匹配结果；
C、将匹配目标随机跳回前一帧图像中未配对过的粒子，将所述未配对过的粒子作为新的第一粒子；
迭代执行步骤B-C，直至得到所有所述第一粒子与所述第二粒子的匹配结果；
D、根据所述信息素、所有所述第一粒子与所述第二粒子的匹配结果，更新所述信息素，得到新的信息素；
迭代执行步骤A-D，直至达到预设的迭代次数或预设的匹配精度，得到所述全局目标函数的全局最优解。
6.根据权利要求5所述的方法，其特征在于，所述根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，确定所述第一粒子和第二粒子组成的粒子对的信息素和启发素，包括：
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标和所述全局目标函数，计算所述启发素；
根据所述第一粒子的三维空间坐标、所述第二粒子的三维空间坐标、所述全局目标函数、信息素的遗留速率和蒸发速率，计算所述信息素；
所述信息素的遗留速率为归一化处理后的遗留速率。
7.根据权利要求1-6所述任一项方法，其特征在于，根据泰森多面体边界面极径形式的相关函数得到所述粒子群相似度法函数。
8.一种速度场测量装置，其特征在于，包括：
获取模块，用于获取至少一个第一粒子的三维空间坐标以及至少一个第二粒子的三维空间坐标，所述第一粒子和所述第二粒子分别为相邻两帧图像中的示踪粒子；
处理模块，用于根据所述至少一个第一粒子的三维空间坐标和所述至少一个第二粒子的三维空间坐标，生成全局目标函数，所述全局目标函数由最小位移法函数和粒子群相似度法函数按照预设权重混合生成；以及，根据所述全局目标函数，得到所述至少一个第一粒子和所述至少一个第二粒子之间的匹配结果，所述匹配结果用于反映流体速度场。
CN202010663981.3A 2020-07-10 2020-07-10 一种三维示踪粒子匹配方法及速度场测量装置 Active CN112036072B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010663981.3A CN112036072B (zh) 2020-07-10 2020-07-10 一种三维示踪粒子匹配方法及速度场测量装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010663981.3A CN112036072B (zh) 2020-07-10 2020-07-10 一种三维示踪粒子匹配方法及速度场测量装置

Publications (2)

Publication Number Publication Date
CN112036072A true CN112036072A (zh) 2020-12-04
CN112036072B CN112036072B (zh) 2023-05-12

Family

ID=73579022

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010663981.3A Active CN112036072B (zh) 2020-07-10 2020-07-10 一种三维示踪粒子匹配方法及速度场测量装置

Country Status (1)

Country Link
CN (1) CN112036072B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113446957A (zh) * 2021-05-08 2021-09-28 北京航空航天大学 基于神经网络标定和散斑追踪的三维轮廓测量方法及装置
CN114859072A (zh) * 2022-05-11 2022-08-05 北京航空航天大学 一种体视粒子追踪测速方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8793200B1 (en) * 2009-09-22 2014-07-29 Hrl Laboratories, Llc Method for particle swarm optimization with random walk
CN108593958A (zh) * 2017-12-15 2018-09-28 北京航空航天大学 同步获取气固两相流速度场的方法和装置
CN109633199A (zh) * 2018-12-11 2019-04-16 天津大学 一种基于相关匹配的粒子聚焦两点像测速方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8793200B1 (en) * 2009-09-22 2014-07-29 Hrl Laboratories, Llc Method for particle swarm optimization with random walk
CN108593958A (zh) * 2017-12-15 2018-09-28 北京航空航天大学 同步获取气固两相流速度场的方法和装置
CN109633199A (zh) * 2018-12-11 2019-04-16 天津大学 一种基于相关匹配的粒子聚焦两点像测速方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
周云龙等: "基于粒子群优化Hopfield网络匹配的稀相颗粒速度测量", 《化工学报》 *
张伟等: "粒子图像测速技术互相关算法研究进展", 《力学进展》 *
聂铭远，潘翀: "改进的蚁群3D-PTV 匹配算法", 《北京力学会第26届学术年会论文集》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113446957A (zh) * 2021-05-08 2021-09-28 北京航空航天大学 基于神经网络标定和散斑追踪的三维轮廓测量方法及装置
CN114859072A (zh) * 2022-05-11 2022-08-05 北京航空航天大学 一种体视粒子追踪测速方法

Also Published As

Publication number Publication date
CN112036072B (zh) 2023-05-12

Similar Documents

Publication Publication Date Title
CN111968229B (zh) 高精度地图制图方法及装置
Chen et al. Transforming a 3-d lidar point cloud into a 2-d dense depth map through a parameter self-adaptive framework
CN113160068B (zh) 基于图像的点云补全方法及系统
Kitzhofer et al. Generation and visualization of volumetric PIV data fields
CN109900274B (zh) 一种图像匹配方法及系统
CN111311663B (zh) 一种实时大场景三维语义建模的方法
CN112036072A (zh) 一种三维示踪粒子匹配方法及速度场测量装置
Bourgoin et al. Using ray-traversal for 3D particle matching in the context of particle tracking velocimetry in fluid mechanics
CA2998676A1 (en) An apparatus and a method for calibrating an optical acquisition system
CN115485718A (zh) 从单个图像确定可穿越空间
Sušanj et al. Effective area coverage of 2D and 3D environments with directional and isotropic sensors
Chen et al. Multi-level scene modeling and matching for smartphone-based indoor localization
He et al. Semi-synthesis: A fast way to produce effective datasets for stereo matching
Panek et al. Visual localization using imperfect 3d models from the internet
Zhang et al. Improved feature point extraction method of ORB-SLAM2 dense map
Han et al. A PTV-based feature-point matching algorithm for binocular stereo photogrammetry
Gautam et al. An experimental comparison of visual SLAM systems
Li et al. 3d fluid flow reconstruction using compact light field piv
Biström Comparative analysis of properties of LiDAR-based point clouds versus camera-based point clouds for 3D reconstruction using SLAM algorithms
Kang et al. 3D urban reconstruction from wide area aerial surveillance video
Yuan et al. Mfusenet: Robust depth estimation with learned multiscopic fusion
Cotra et al. Lidar-based methods for tracking and identification
Frenzel et al. Quantitative comparison between neural network-and SGM-based stereo matching
Irmisch Visual localization in challenging environments
Höhner et al. Particle Filter Based Tracking and Mapping

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant