CN109308712B - 一种从液滴流视频中计算液滴运动频率的方法 - Google Patents

一种从液滴流视频中计算液滴运动频率的方法 Download PDF

Info

Publication number
CN109308712B
CN109308712B CN201710620334.2A CN201710620334A CN109308712B CN 109308712 B CN109308712 B CN 109308712B CN 201710620334 A CN201710620334 A CN 201710620334A CN 109308712 B CN109308712 B CN 109308712B
Authority
CN
China
Prior art keywords
frequency
droplet
drop
motion
similarity vector
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.)
Active
Application number
CN201710620334.2A
Other languages
English (en)
Other versions
CN109308712A (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.)
Beijing Targeting One Biotechnology Co ltd
Tsinghua University
Original Assignee
Beijing Targeting One Biotechnology Co ltd
Tsinghua 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 Beijing Targeting One Biotechnology Co ltd, Tsinghua University filed Critical Beijing Targeting One Biotechnology Co ltd
Priority to CN201710620334.2A priority Critical patent/CN109308712B/zh
Publication of CN109308712A publication Critical patent/CN109308712A/zh
Application granted granted Critical
Publication of CN109308712B publication Critical patent/CN109308712B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/20Analysis of motion
    • G06T7/246Analysis of motion using feature-based methods, e.g. the tracking of corners or segments
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/14Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/22Matching criteria, e.g. proximity measures
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10016Video; Image sequence

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • General Engineering & Computer Science (AREA)
  • Mathematical Optimization (AREA)
  • Artificial Intelligence (AREA)
  • Multimedia (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Evolutionary Computation (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Image Analysis (AREA)

Abstract

本发明涉及一种从液滴流视频中计算液滴运动频率的方法,所述方法包括以下步骤:载入液滴流视频及其采集频率;选取所述液滴流视频中的一帧作为参考帧;依次计算所述参考帧与所述液滴流视频中的每一帧之间的相似度,构造相似度向量;和根据所述相似度向量,计算所述液滴流视频的液滴运动频率。本发明的液滴运动频率计算方法能广泛地应用于液滴微流控技术中,特别是并行化、阵列化的液滴微流控技术中,液滴的通量和均一性的表征。

Description

一种从液滴流视频中计算液滴运动频率的方法
技术领域
本发明涉及液滴微流控技术领域,尤其涉及从一种从液滴流视频中计算液滴运动频率的方法。
背景技术
液滴微流控技术能够高通量地制备尺寸高度均一的液滴序列,在高通量均一材料的制备和高通量生化分析等领域有广泛的应用。为了进一步提高制备或分析通量,液滴微流控单元往往还会并行化、阵列化使用。这些液滴序列流过这些并行液滴微流控单元中每个通道的频率,即“液滴运动频率”,代表该器件在单位时间内的最高分析通量,通道间液滴运动频率的变异系数也是表征液滴均一性的关键指标。
目前,能用于液滴运动频率检测的方法主要有:光学、电学和图像处理等方法。光学方法主要包括光散射等方法。虽然在使用CCD等成像装置的条件下,光散射方法能够进行多通道并行的液滴运动频率检测,但是这种方法需要在液滴微流控技术平台上定制额外的设备,且相对检测误差普遍在15%到25%的范围内,不够准确。电学方法主要包括电阻检测和电容检测等方法。电阻检测方法需要在微流控芯片内增加与液滴接触的电极,电容检测方法可以在微流控芯片内增加与液滴不接触的电极,虽然这两种方法分别能将相对检测误差分别降低到约3.9%和约5.1%,但是这两种方法都需要在液滴微流控芯片内增加电极,且对于每一个液滴流过的通道,都需要一套单独的电路系统来进行检测,不便于多通道并行的液滴运动频率检测。图像处理法包括微粒子测速等方法。微粒子测速方法在计算频率的过程中需要使用相邻液滴的间距,当某帧图像中只有1个液滴时,微粒子测速方法将无法计算液滴的运动频率;当某帧图像中有3个或更多个液滴时,微粒子测速方法可能会造成相邻液滴间隔的误判,计算会出现错误(计算出液滴运动频率的一个分频);此外,微粒子测速法只能在液滴运动方向相同、液滴运动频率一致的条件下进行多通道并行液滴运动频率的计算。
总之,现有的各个方法难以同时做到以下三点:第一,使用常规液滴生成平台,不在液滴微流控平台或液滴微流控芯片上增加额外的定制装备或定制结构;第二,准确检测或计算液滴的运动频率,目前液滴微流控技术制备的液滴,直径的变异系数一般约为1%,在分散相流量不变的情况下(一般液滴微流控的应用场合都满足该条件),可以折合成液滴运动频率的变异系数约为3%,液滴运动频率的检测或计算方法的相对误差应小于该值,才能可靠地表征液滴序列在均一性上的差异;第三,随着液滴微流控的使用趋于并行化、阵列化,液滴运动频率的检测方法应当适应这种趋势,具备对多通道和液滴尺寸、多液滴运动频率和多液滴运动方向的并行液滴运动进行频率检测或计算的能力。
发明内容
本发明的目的在于提供从液滴流视频中计算液滴运动频率的方法,旨在解决现有技术在液滴运动频率计算中存在的应用不方便、结果不准确和难以进行多通道、多尺寸、多频率的并行液滴运动频率检测等问题。
在一种实施方式中,本发明一种从液滴流视频中计算液滴运动频率的方法,所述方法包括以下步骤:步骤1,载入液滴流视频及其采集频率;步骤2,选取所述液滴流视频中的一帧作为参考帧;步骤3,依次计算所述参考帧与所述液滴流视频中的每一帧之间的相似度,构造相似度向量;和步骤4,根据所述相似度向量,计算所述液滴流视频的液滴运动频率。
在一种实施方式中,当所述相似度向量的周期性明显,则从所述相似度向量中提取峰值,通过相邻峰值间的平均时间差计算相似度向量的周期,计算所述液滴流视频的液滴运动频率。
在一种实施方式中,当所述相似度向量的周期性不明显,则对所述相似度向量进行时域分析,消除所述相似度向量中的非周期分量,保留和放大周期分量,然后通过相邻峰值间的平均时间差计算相似度向量的相关函数的周期,计算所述液滴流视频的液滴运动频率。
在一种实施方式中,所述时域分析包括基线校正、幅值变换和时域噪声消除的一种或多种。
在一种实施方式中,当相似度向量中存在多个频率分量,则计算所述相似度向量或其相关函数的频谱,进行频域分析,计算相似度向量或其相关函数的频谱,计算所述液滴流视频的液滴运动频率。
在一种实施方式中,所述频域分析包括:对所述相似度向量进行时域-频域的变换;和对所述变换的结果进行频谱分析。
在一种实施方式中,所述频谱分析包括频谱噪声消除、峰值分析、基频提取、频谱截断和频谱密度加权平均中的一种或多种。
在一种实施方式中,当所述相似度向量的频谱中存在分布集中的频率分量,则对所述频率分量的峰值进行提取,所提取的频率分量直接作为所述液滴流视频的液滴运动频率。
在一种实施方式中,当所述相似度向量的频谱分量连续分布于一定的频率范围中,则在设定的频率范围内对强度不小于阈值的频率分量进行加权平均,将加权平均得到的频率作为所述液滴流视频的液滴运动频率。
在一种实施方式中,步骤1所述液滴流视频中的所有液滴的中心沿着至少1条确定的轨迹运动。
在一种实施方式中,当所述液滴流视频中的所有液滴中心沿着至少2条确定的轨迹运动,则任意2条不同的所述运动轨迹均不相交。
在一种实施方式中,步骤1中所述液滴流视频的采集频率不小于液滴运动频率的2倍。
本发明提供的从液滴流视频中计算液滴运动频率的方法,该方法使用液滴流运动的周期性进行液滴运动频率的计算,液滴流运动的周期性表现为每隔一段时间液滴就会出现在同样的位置。使用商业化的高速图像采集设备对液滴流的运动进行视频采集(采集频率不小于液滴运动频率的2倍),则视频中各帧之间的相似度会随着液滴的运动过程呈现周期性变化,计算出视频中各帧之间相似度,即可用其变化的周期表征液滴流运动的周期,进而通过时域和/或频域的分析方法,准确地从中计算出液滴的运动频率。这种方法不需要进行图像识别,对液滴并行流动的通道个数,微流道和液滴的尺寸,以及各通道内液滴的运动方向等没有要求,只要求液滴的中心在各自的通道内沿固定的轨迹运动,且任意两条不同的轨迹互不相交即可。在绝大多数液滴微流控的应用场合,一般都存在一段微流道,其中液滴的中心沿着该微流道的中心线运动,在绝大多数并行液滴微流控的应用场合,一般也都存在一系列的微流道片段,其中的任意两条微流道互不相交,因此在绝大多数液滴微流控的应用场合,上述两个条件都容易得到满足。有鉴于此,本发明提出的液滴运动频率计算方法能广泛地应用于液滴微流控技术中,特别是并行化、阵列化的液滴微流控技术中,液滴的通量和均一性的表征。
本发明的方法只使用商业化的高速图像采集设备采集的液滴流视频,无需在液滴微流控技术平台上定制额外的设备或部件,即可计算液滴流中液滴的运动频率。当不同尺寸的液滴以不同频率流过单一通道,则在横跨2个数量级的液滴运动频率下,所有计算误差均小于0.3%;当不同尺寸的液滴以不同频率和不同方向流过多个宽度不同的通道,则所有计算误差均小于1.6%。因此,本发明的方法是一种便于应用、准确可靠的液滴运动频率计算方法。
附图说明
为了更清楚地说明本申请实施例中的技术方案,下面将对实施例中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本申请中记载的一些实施例,对于本领域普通技术人员来说,在不付出创造性劳动的前提下,还可以根据这些附图获得其它的附图。
图1是本发明的一种从液滴流视频中的液滴运动频率计算方法的装置示意图;
图2是本发明的一种从液滴流视频中的液滴运动频率计算方法的实施流程图;
图3是本发明的一种从液滴流视频中的液滴运动频率计算方法所载入的单一通道液滴流视频中的图像序列;
图4是依据图3所示的液滴流视频计算的相似度向量;
图5是对图4所示的相似度向量的一种时域分析方法:计算相似度向量的循环自相关函数的结果;
图6是对图4所示的相似度向量的一种频域分析方法:计算相似度向量的循环自功率谱密度函数的结果;
图7是对图6所示的相似度向量的频域分析方法中的一种频谱分析方法:基于降噪阈值和主频带内功率谱密度加权的平均频率计算结果,其中7a和7b两条线之间的频段为主频带,7c线为降噪阈值,7d线为所计算的平均频率的位置;
图8是本发明根据表2中的数据进行绘制的在不同液滴运动频率下的计算误差;
图9是本发明的方法所载入的多流量条件、多通道和液滴尺寸、多运动方向的并行液滴流视频(由多个不同流量组合的液滴流视频经平移、旋转和缩放拼合而成)中的一帧图像;
图10是依据图9所示的液滴流视频计算的相似度向量;
图11是对图10所示的相似度向量的一种时域分析方法:计算相似度向量的循环自相关函数的结果;和
图12是对图10所示的相似度向量的一种频域分析方法:计算相似度向量的循环自功率谱密度函数的结果,以及其中的一种频谱分析方法:带降噪的局部峰值提取方法的计算结果,其中12a线为降噪阈值。
具体实施方式
为了使本领域技术领域人员更好地理解本申请中的技术方案,下面将结合下面结合实施例对本发明作进一步说明,显然,所描述的实施例仅仅是本申请一部分实施例,而不是全部的实施例。基于本申请中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其它实施例,都应当属于本申请保护的范围。
实施例一:给定流量条件下单一通道内液滴运动频率的计算
请参阅图1所示的装置示意图,此例中,注射泵141和142(也可以用气压源)驱动连续相102和分散相103进入流动聚焦微流控芯片10,产生大小均一的液滴流101。经过显微镜13放大后,该运动过程被高速图像采集设备12采集,采集频率不低于液滴运动频率的2倍。采集后的液滴流视频16和高速图像采集设备12的采集频率被输入计算机15,使用本发明中所述的液滴运动频率计算方法进行计算,最终输出液滴的运动频率。
在本实施例中,连续相102和分散相103分别以12000μL/h和2000μL/h流动,流动聚焦微流控芯片10的喷口特征尺寸为75μm×75μm、微流道深度为75μm。液滴运动频率的估算方法如下:在流动聚焦微流控芯片中,液滴的体积一般不小于喷口区域的体积,由此可以估算液滴的运动频率。根据香农采样定理,采集频率不应小于该频率的2倍,因此最低采集频率为2634Hz,为保险起见,将采集频率设定为12000Hz。本段所述的实验条件和估算方法,其目的是确定采集频率,采集频率是本实施例中所述液滴运动频率计算方法的一个输入,因此本段所述估算过程不属于液滴运动频率计算方法中的步骤。
从液滴流视频中计算液滴运动频率有多种实施途径,具体可参阅图2显示的流程图。在图2中,液滴流运动视频S201.1和相应的采集频率S201.2作为液滴运动频率计算方法的输入。在步骤S202中,从液滴流运动视频S201.1中取一帧作为参考帧。在步骤S203中,依次计算步骤S202中的参考帧与液滴流运动视频S201.1中各帧的相似度,并结合采集频率S201.2计算液滴流运动视频S201.1中各帧的时间位置,构造相似度向量。在后续的计算方法中,将有不同的实施途径可供选择:
实施途径1:当相似度向量的周期性明显,不需要进行校正基线、增强幅值或消除噪声等处理,则可以进行步骤S206,从相似度向量中提取峰值,通过相邻峰值间的平均时间差计算相似度向量的周期,求倒数后得到液滴的运动频率,在步骤S209中进行输出;
实施途径2:当相似度向量的周期性不明显,需要进行校正基线、增强幅值和消除噪声等处理,则可以进行步骤S204进行时域分析,计算相似度向量的相关函数,消除相似度向量中的非周期分量(含直流分量),保留和放大周期分量,从而达到校正基线、增强幅值和消除噪声的效果。然后可以进行步骤S206,从相似度向量的相关函数中提取峰值,通过相邻峰值间的平均时间差计算相似度向量的相关函数的周期,求倒数后得到液滴的运动频率,在步骤S209中进行输出;
当相似度向量中存在多个频率分量,难以直接求解相似度向量或其相关函数的周期,则可以进行步骤S205进行频域分析,计算相似度向量或其相关函数的频谱,将相似度向量或其相关函数的频率分量直观地展现出来。为求取液滴的运动频率,需要进行频谱分析,频谱分析有两种实施途径可供选择:
实施途径3:当相似度向量的频谱中存在分布集中的频率分量,则可以进行步骤S207,对这些频率分量的峰值进行提取,所提取的频率分量可以直接作为液滴的运动频率,在步骤S209中进行输出;
实施途径4:当相似度向量的频谱分量连续分布于一定的频率范围中,则可以进行步骤S208,在设定的频率范围内对强度不小于阈值的频率分量进行加权平均,综合考虑该频率范围内各频率分量对频谱的贡献,将该平均频率作为液滴的运动频率,在步骤S209中进行输出。
下面对图2中的各个步骤进行详细的解释说明。
步骤S201:图3为本实施例中高速图像采集设备12采集到的液滴流视频中的一段帧序列。从图3中可以看出,所有液滴的中心均沿着微流道的中心线运动,且沿该中心线呈等距排列(见图中虚线),大约每10帧产生一个液滴并沿着微流道从左向右运动到前一个液滴的位置(见图中实线),因此液滴的运动频率
Figure BDA0001361529790000051
步骤S202:选择视频的第一帧作为参考帧。当然,由于本实施例中的液滴流视频每间隔10帧就会出现一帧相似的图像,比如第1帧和第11帧、第21帧……第111帧……等相似,第2帧和第12帧、第22帧……第112帧……等相似,以此类推,因此选择视频中的任意一帧作为参考帧,都不会明显地影响计算结果。
步骤S203:依次计算视频参考帧与各帧之间相似度,这里所述的相似度可以采用两帧之间对应位置的颜色向量的余弦相似度,该相似度描述两帧中颜色分布的整体线性相关程度,当所有颜色向量的分量全部不小于0(比如灰度值、RGB或CMYK等颜色空间)时,该相似度的值在0~1之间(含0和1),当两帧各像素的颜色向量相同或成比例,该相似度为1,表示完全相似,该相似度的计算方法如下:
Figure BDA0001361529790000061
其中R(x,y)表示第x帧和第y帧之间对应位置的颜色向量的余弦相似度,D(j,i,k,x)表示液滴流视频中第x帧在坐标(i,j)处的第k个颜色分量。
参考帧与视频中每一帧之间的相似度构成相似度向量,其每个分量的定义如下:
S(iΔt)=R(i,iref),i=1,2,…,N
其中S表示相似度向量,iref表示参考帧的编号,Δt表示采样间隔(即采样频率的倒数),N表示液滴流视频的总帧数。
计算所得相似度向量如图4所示,图中横坐标是各帧在液滴流视频采集过程中的时间位置(即iΔt,i=1,2,…,N),纵坐标是通过计算得到的相似度向量的各个分量的值:S(iΔt),i=1,2,…,N。在液滴流视频中,第1帧与第11帧、第21帧……第111帧……等相似,因此取第1帧作为参考帧时,第11帧、第21帧……第111帧……等应与参考帧的相似度较高,这在图中相应的时间位置能够得到印证,比如第11帧(对应时间位置
Figure BDA0001361529790000062
)和第111帧(对应时间位置
Figure BDA0001361529790000063
)等位置出现了较高的相似度,因此当液滴流视频中各帧变化的周期与相似度变化的周期完全对应。图4中的两个X和Y代表第1个峰和第11个峰对应的时间位置X和该时间位置对应的相似度向量的分量的值Y;其值通过峰值提取得到。
步骤S204:当相似度向量的周期性表现得不明显时,则可以对相似度向量进行时域分析,去除其中的非周期分量。相似度向量中的非周期分量包括:基线的无规则漂移,直流分量和毛刺噪声等。这种时域分析可以使用相关函数进行实现。考虑到液滴流的运动是无穷周期运动,高速图像采集设备采集的只是其中的一段,因此所述相关函数可以采用循环自相关函数,定义如下:
Figure BDA0001361529790000064
其中rx(mΔt)表示循环自相关函数,S*(iΔt)表示扣除直流分量的相似度向量,其定义如下:
Figure BDA0001361529790000065
本实施例中,相似度向量(图4)的循环自相关函数如图5所示。虽然本实施例中相似度向量的周期性表现得非常明显,但是从图5中还是可以看出,经过相关分析,相似度向量(图4)的基线漂移和直流分量全部被消除,周期分量的幅值被放大了约100倍,且相似度向量中的毛刺噪声(表现为图4波谷中的微小凹凸)也被成功消除;与此同时,相似度向量中的周期分量被保留下来,且周期没有发生明显的变化,这在步骤S206(使用相关函数)中可以加以证明。图5中的两个X和Y代表第1个峰和第16个峰对应的时间位置X和该时间位置对应的相似度向量的分量的值Y;其值通过峰值提取得到。
步骤S206:在计算出相似度向量或相似度向量的相关函数后,可以采用周期分析的方法计算液滴的运动频率,这里所述的周期可以采用峰值的时间差与峰值数之差的比例进行计算:
Figure BDA0001361529790000071
其中T表示液滴运动的周期,ti和tj分别表示第i个峰值和第j个峰值对应的时间位置。
无论使用相似度向量还是使用相似度向量的相关函数,计算出的周期是完全一致的。
使用相似度向量:根据图4中的结果,可以根据第1个峰值和第11个峰值计算出相似度向量的周期为:
Figure BDA0001361529790000072
相应的液滴运动频率为:
Figure BDA0001361529790000073
(作为步骤S209中的输出)。
使用相似度向量的相关函数:根据图5中的结果,可以根据第1个峰值和第16个峰值计算出相似度向量的周期为:
Figure BDA0001361529790000074
相应的液滴运动频率为:
Figure BDA0001361529790000075
(作为步骤S209中的输出),与使用相似度向量时的计算结果完全一致。
步骤S205:有时相似度向量中存在多种不同的周期分量,难以直接使用基于峰值提取的周期分析方法计算液滴的运动频率。这时可以使用频域分析方法,将相似度向量或相似度向量的相关函数的各个周期分量以频谱的方式进行直观的展现,方便液滴运动频率的计算。
由于相似度向量与相似度向量的相关函数具有相同的周期分量,因此两者的频域分析过程没有差别。以下使用相似度向量的相关函数为例进行频域分析。
频域分析方法分为时域-频域的变换和频谱分析两个步骤。时域-频域的变换在步骤S205中完成,使用的变换可以是离散傅里叶变换:
Figure BDA0001361529790000076
其中
Figure BDA0001361529790000077
表示循环自相关函数rx(mΔt)的离散傅里叶变换结果,rx(mΔt)表示循环自相关函数,Δt表示采样间隔(即采样频率的倒数),N表示液滴流视频的总帧数。
对相似度向量的循环自相关函数进行离散傅里叶变换,结果为相似度向量的循环自功率谱密度函数,如图6所示。图6中的两个X和Y代表图中最大的循环自功率谱密度Y和对应的频率X,该点可以通过峰值提取获得。频谱分析将在步骤S207或S208中完成。
步骤S207:本步骤是实现频谱分析的一种手段。当频谱中的频率分量集中且不为0时(频率为0的分量为直流分量而非周期分量),可以直接采用主频率表示液滴的运动频率。在单一通道的情况下,主频率可以通过峰值提取的方法进行计算:
Figure BDA0001361529790000081
其中f主要表示频谱的主频率,p(f)表示频谱密度函数,f表示频率。
在图6的频谱中,频率分量单一且不为0,因此可以通过峰值提取的方法,计算出主频率为:f主要=1997.62Hz,即液滴的运动频率为:fS207=1197.62Hz(作为步骤S209中的输出)。
步骤S208:本步骤是实现频谱分析的另一种手段。当频谱中的频率分量不单一时,可以采用主频带内的平均频率表示液滴的运动频率。主频带是以主频率为中心、具有一定宽度的频带:
F主要=[f主要-Δf,f主要+Δf],Δf>0
其中F主要表示主频带,Δf表示主频带的半带宽。
在本实施例中,取
Figure BDA0001361529790000082
因此主频带为[998.81,2996.43]Hz,对应于图7中7a和7b两条线之间的频率范围。
主频带内的平均频率可以采用基于频谱强度加权平均的方法进行计算。考虑到频谱噪声对均值可能造成的影响,计算时只考虑频谱强度不小于特定阈值的频率分量:
Figure BDA0001361529790000083
其中
Figure BDA0001361529790000084
表示主频带内的平均频率,p0表示降噪阈值。
在本实施例中,取
Figure BDA0001361529790000085
对应于图7中的7c线。
在本实施例中,主频带内平均频率的计算结果为:
Figure BDA0001361529790000086
对应于图7中的7d线,因此液滴的运动频率为:fS208=1197.62Hz(作为步骤S209中的输出),与步骤S207的结果一致。
步骤S209:本步骤完成对步骤S206的计算结果fS206,S或fS206,r、步骤S207的计算结果fS207或步骤S208的计算结果fS208的输出。
以上介绍了基于本发明的单一通道内液滴运动频率计算的4种实施途径。为了评价它们的准确性,可以对单一通道内的液滴运动进行手动计时和计数,获得准确的液滴运动频率:
Figure BDA0001361529790000091
其中n液滴表示从第x帧到与该帧相似度最高的另一帧(第y帧)之间,通过通道内某一截面的液滴个数。
在本实施例的液滴流视频中,从第1帧到第4990帧之间共有499个液滴通过图3所示通道内左侧实线所示的截面,因此可以计算出液滴的运动频率为:
Figure BDA0001361529790000092
各实施途径的相对计算误差如表1所示,所有的相对计算误差均<0.3%。
表1:实施例一中,液滴运动频率的测量值,通过不同实施途径求解的液滴运动频率的计算值,以及相对计算误差(采集频率为12000Hz)。
Figure BDA0001361529790000093
综上所述,本实施例证明,根据本发明揭露一种从液滴流视频中计算液滴运动频率的方法,可以在给定流量条件下进行单通道内液滴运动频率的计算,相对计算误差<0.3%。
实施例二:不同流量条件下单一通道内液滴运动频率的计算
请参阅图1所示的装置示意图,此例中,使用实施例一的流动聚焦微流控芯片10,改变连续相102和分散相103的流量,可以得到不同运动频率的液滴序列。使用高速图像采集设备12分别在合适的采集频率(请参见实施例一中的液滴运动频率估算方法)下对这些液滴的运动过程进行记录,就可以得到运动频率不同的液滴流视频16。使用本发明所述的液滴运动频率计算方法分别对这些液滴流视频和它们的采集频率进行液滴运动频率的计算,就可以在不同液滴运动频率下评价该液滴运动频率计算方法的准确性。
考虑到需要对不同的液滴流视频实施液滴运动频率的计算,使用较通用的实施途径进行液滴运动频率的计算,根据实施例一中对图2的各个实施途径的分析,采用“频域分析-使用主频带内的平均频率”(S201→S202→S203→S204→S205→S208→S209)这一实施途径进行计算。相应步骤的计算方法请参考实施例一中的说明。
表2显示不同连续相和分散相流速下,液滴运动频率的测量值、计算值和相对误差。从表中可以看出,调节连续相和分散相的流速,可以实现40~5000Hz范围内的液滴运动频率,所有液滴运动频率下的相对计算误差均<0.3%。
表2:实施例二中,不同连续相流量、分散相流量和相应的采集频率条件下,液滴运动频率的测量值,液滴运动频率的计算值,以及相对计算误差。
Figure BDA0001361529790000101
图8显示不同的液滴运动频率测量值与相对计算误差之间的关系,从图中可以看出,所有相对计算误差<0.3%,,且分布与液滴运动频率的测量值无关,因此液滴运动频率的计算误差与液滴的运动频率无关,为随机误差。
综上所述,本实施例证明,本发明揭露一种从液滴流视频中计算液滴运动频率的方法,可以在不同流量条件下进行单通道内液滴运动频率的计算,所有的相对计算误差均<0.3%。
实施例三:多通道和液滴尺寸、多液滴运动频率和多液滴运动方向的并行液滴流的运动频率计算
此例中,首先在三种不同流量条件(见表3)下,按图1所示的方法,使用流动聚焦微流控芯片10得到不同运动频率的液滴序列。根据表2中的数据,这三种不同流量条件下,液滴的最高运动频率<2000Hz,因此设定高速图像采集设备12的采集频率为20000Hz,依次对这些液滴的运动过程进行记录,得到三个运动频率不同的液滴流视频16。使用实施例二中所述的方法,就可以分别计算出这三个液滴流视频中的液滴运动频率,计算结果和相对计算误差如表3所示。
表3:实施例三中,参与拼合多通道和液滴尺寸、多液滴运动频率、多液滴运动方向的并行液滴流视频的单一通道液滴流视频的液滴运动频率的测量值,液滴运动频率的计算值,以及相对计算误差(采集频率均为20000Hz)。
Figure BDA0001361529790000111
将这三个液滴流视频按不同的缩放比例和方向进行画面拼合(参数见表4),就可以制作成多通道和液滴尺寸、多液滴运动频率和多液滴运动方向的并行液滴流视频。图9为拼合后的并行液滴流视频中的一帧,图中的箭头系额外添加,用于指示液滴的运动方向,不属于视频帧的内容。在图中,所有液滴的中心沿着3条确定的轨迹运动,这3条确定的轨迹分别为图中3条微流道的中心线,其中任意2条不同的轨迹均不相交,沿着同一条轨迹运动的液滴的中心沿着该轨迹呈等距排列。
表4:实施例三中,使用单一通道液滴流视频制作多通道和液滴尺寸、多液滴运动频率、多液滴运动方向的并行液滴流视频时,各单一通道液滴流视频使用的平移、旋转和缩放参数。
液滴流视频编号 中心位置 旋转角度(deg) 缩放比例
1 50%,30.47% 3 200%
2 50%,58.59% -13 110%
3 73.44%,59.77% 147 150%
注:中心位置指各液滴流视频的中心在拼合后画面中距离左上角的位置,逗号前后的数值分别以拼合后画面的宽和高作为100%。
可以采用“频域分析-使用频谱峰值”(S201→S202→S203→S204→S205→S207→S209)这一实施途径对该多流量条件、多尺寸、多方向的并行液滴流视频进行运动频率计算:取第1帧为参考帧,计算参考帧与视频中各帧之间的相似度,构成相似度向量,结果如图10所示;计算相似度向量的循环自相关函数,结果如图11所示。从图10和图11来看,相似度向量及其循环自相关函数似乎具有一定的周期性,但由于两者都表现出明显的拍频现象,因此其中存在多种频率分量。根据实施例一的分析,这种情况适合采用频域分析方法。这里可以使用离散傅里叶变换计算相似度向量的循环自功率谱密度函数,如图12所示,图中的三个X和Y代表图中循环自功率谱密度的三个局部峰值Y和对应的三个频率X,这些点可以通过局部峰值提取获得。从图12中可以看出,相似度向量及其循环自相关函数中存在3个主要的频率分量,因此适合采用带降噪的局部峰值提取的方法进行多种液滴运动频率的计算:
Figure BDA0001361529790000121
其中符号的定义与实施例一中的对应符号相同。
可以取降噪阈值:
Figure BDA0001361529790000122
对应于图12中的12a线。经过计算,可以得出图12的频谱中共有3个主要分量:399.2Hz,1577Hz和1916Hz,分别与表3中3个液滴流视频的运动频率接近。在计算得出的所有3个频率中,最大相对计算误差<1.6%。
综上所述,本实施例证明,根据本发明揭露一种从液滴流视频中计算液滴运动频率的方法,可以进行多通道和液滴尺寸、多液滴运动频率和多液滴运动方向的并行液滴流的运动频率计算,最大相对计算误差<1.6%。
综上所述,根据本发明揭露一种从液滴流视频中计算液滴运动频率的方法,主要解决使用常规液滴微流控平台和器件的情况下,准确地提取液滴运动频率的问题。该方法使用液滴流视频及其采集频率,计算选定的参考帧与各帧的相似度,然后可以在不同的实际情况和应用需求下,选用4种不同的实施途径计算液滴的运动频率。在单通道液滴运动频率计算的应用中,在超过两个数量级的液滴运动频率范围内,相对计算误差均<0.3%;在多通道和液滴尺寸、多液滴运动频率和多液滴运动方向的并行液滴流视频的运动频率计算的应用中,最大相对计算误差<1.6%。
应该理解到披露的本发明不仅仅限于描述的特定的方法、方案和物质,因为这些均可变化。还应理解这里所用的术语仅仅是为了描述特定的实施方式方案的目的,而不是意欲限制本发明的范围,本发明的范围仅受限于所附的权利要求。
本领域的技术人员还将认识到,或者能够确认使用不超过常规实验,在本文中所述的本发明的具体的实施方案的许多等价物。这些等价物也包含在所附的权利要求中。

Claims (4)

1.一种从液滴流视频中计算液滴运动频率的方法,其特征在于,所述方法包括以下步骤:
步骤1,载入液滴流视频及其采集频率;
步骤2,选取所述液滴流视频中的一帧作为参考帧;
步骤3,依次计算所述参考帧与所述液滴流视频中的每一帧之间的相似度,构造相似度向量;和
步骤4,根据所述相似度向量,计算所述液滴流视频的液滴运动频率;
当所述相似度向量的周期性明显,则从所述相似度向量中提取峰值,通过相邻峰值间的平均时间差计算相似度向量的周期,计算所述液滴流视频的液滴运动频率;
当所述相似度向量的周期性不明显,则对所述相似度向量进行时域分析,消除所述相似度向量中的非周期分量,保留和放大周期分量,然后通过相邻峰值间的平均时间差计算相似度向量的相关函数的周期,计算所述液滴流视频的液滴运动频率;
所述时域分析包括基线校正、幅值变换和时域噪声消除的一种或多种;
当相似度向量中存在多个频率分量,则计算所述相似度向量或其相关函数的频谱,进行频域分析,计算相似度向量或其相关函数的频谱,计算所述液滴流视频的液滴运动频率所述频域分析包括:对所述相似度向量进行时域-频域的变换;和对所述变换的结果进行频谱分析;
所述频谱分析包括频谱噪声消除、峰值分析、基频提取、频谱截断和频谱密度加权平均中的一种或多种;
当所述相似度向量的频谱中存在分布集中的频率分量,则对所述频率分量的峰值进行提取,所提取的频率分量直接作为所述液滴流视频的液滴运动频率;当所述相似度向量的频谱分量连续分布于一定的频率范围中,则在设定的频率范围内对强度不小于阈值的频率分量进行加权平均,将加权平均得到的频率作为所述液滴流视频的液滴运动频率。
2.根据权利要求1所述的方法,其特征在于,步骤1所述液滴流视频中的所有液滴的中心沿着至少1条确定的轨迹运动。
3.根据权利要求1所述的方法,其特征在于,当所述液滴流视频中的所有液滴中心沿着至少2条确定的轨迹运动,则任意2条不同的运动轨迹均不相交。
4.根据权利要求1所述的方法,其特征在于,步骤1中所述液滴流视频的采集频率不小于液滴运动频率的2倍。
CN201710620334.2A 2017-07-26 2017-07-26 一种从液滴流视频中计算液滴运动频率的方法 Active CN109308712B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710620334.2A CN109308712B (zh) 2017-07-26 2017-07-26 一种从液滴流视频中计算液滴运动频率的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710620334.2A CN109308712B (zh) 2017-07-26 2017-07-26 一种从液滴流视频中计算液滴运动频率的方法

Publications (2)

Publication Number Publication Date
CN109308712A CN109308712A (zh) 2019-02-05
CN109308712B true CN109308712B (zh) 2021-10-26

Family

ID=65202781

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710620334.2A Active CN109308712B (zh) 2017-07-26 2017-07-26 一种从液滴流视频中计算液滴运动频率的方法

Country Status (1)

Country Link
CN (1) CN109308712B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117455922B (zh) * 2023-12-26 2024-04-05 墨卓生物科技(浙江)有限公司 一种基于液滴运动图像的自动计数分析方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6128433A (en) * 1997-12-31 2000-10-03 Indigita Corporation Controlling frequency drift of a digital video tape drive
CN102792145A (zh) * 2010-03-09 2012-11-21 贝克曼考尔特公司 计算流式细胞仪的液滴延迟时间的系统和方法
CN104826673A (zh) * 2015-03-16 2015-08-12 中国科学院微生物研究所 写入式二维微流控液滴阵列化装置、用途及其使用方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6128433A (en) * 1997-12-31 2000-10-03 Indigita Corporation Controlling frequency drift of a digital video tape drive
CN102792145A (zh) * 2010-03-09 2012-11-21 贝克曼考尔特公司 计算流式细胞仪的液滴延迟时间的系统和方法
CN104826673A (zh) * 2015-03-16 2015-08-12 中国科学院微生物研究所 写入式二维微流控液滴阵列化装置、用途及其使用方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
《用于精子筛选及体外受精的微流控芯片的构建及应用》;谢兰;《中国博士学位论文全文数据库 医药卫生科技辑》;20120615;全文 *

Also Published As

Publication number Publication date
CN109308712A (zh) 2019-02-05

Similar Documents

Publication Publication Date Title
Biferale et al. Particle trapping in three-dimensional fully developed turbulence
CN104680483B (zh) 图像的噪声估计方法、视频图像去噪方法及装置
CN105550694B (zh) 一种度量人脸图像模糊程度的方法
CN103054569B (zh) 基于可见光图像测量人体心率的方法、装置及手持设备
Morháč et al. High-resolution boosted deconvolution of spectroscopic data
CN104272347B (zh) 去除包含在静止图像中的雾的图像处理装置及其方法
CN109308712B (zh) 一种从液滴流视频中计算液滴运动频率的方法
Chenouard et al. Curvelet analysis of kymograph for tracking bi-directional particles in fluorescence microscopy images
Vinista et al. A novel modified sobel algorithm for better edge detection of various images
Gu et al. Region sampling for robust and rapid autofocus in microscope
Wang et al. Spatiotemporal saliency model for small moving object detection in infrared videos
US9721179B2 (en) Line segment and arc detection apparatus
CN110288564B (zh) 基于功率谱分析的二值化散斑质量评价方法
US20170295314A1 (en) Focusing method
CN108335310A (zh) 一种便携式粒形粒度检测方法及系统
Priya et al. Retrospective non-uniform illumination correction techniques in images of tuberculosis
CN109211138A (zh) 一种外形检测系统和方法
WO2015114750A1 (ja) フローサイトメータ
US10152774B2 (en) Method and system for estimating point spread function
Zantow et al. Image-based analysis of droplets in microfluidics
CN103914695A (zh) 一种用于微电泳图像识别装置及其方法
Choudhary et al. A novel approach for edge detection for blurry images by using digital image processing
US8000499B2 (en) Automated determination of cross-stream wind velocity from an image series
Janoriya et al. Critical review on edge detection techniques in spatial domain on low illumination images
Taguchi et al. Spectrum imaging measurements with semi-parallel detection using an AES apparatus

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
CB02 Change of applicant information
CB02 Change of applicant information

Address after: 100084 Tsinghua Yuan, Beijing, Haidian District

Applicant after: Tsinghua University

Applicant after: Beijing Xinyi Biotechnology Co., Ltd.

Address before: 100084 Tsinghua Yuan, Beijing, Haidian District

Applicant before: Tsinghua University

Applicant before: Beijing Tianjian Wellcome Biotechnology Co. Ltd.

SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant