本发明涉及数字图像处理技术领域，具体涉及一种基于扑翼飞行系统的实时电子稳像方法。
背景技术：
近年来，国内外学者和科研人员对仿生扑翼式飞行器进行大量的研究，取得了一定成就。扑翼飞行器的气动性能理论比固定翼和旋翼飞行器要复杂。尤其在低雷诺数下，扑翼飞行器能充分利用非定常空气动力学特性获取更高的升力，在相同能量的消耗下，其续航时间、里程都能更长、能够广泛运用于表演、航拍、监控、测绘、救灾等场景。
视频图像是扑翼飞行器自身导航以及获取任务执行目标信息最主要的手段之一，然而相比于固定翼或者四旋翼飞行器，扑翼飞行器的理论研究不够深入、完善，设计制作出的扑翼飞行器远不如固定翼或者旋翼飞行器稳定。扑翼飞行器自重轻，飞行时翼展扑动带来的高频震动，加上外界气流的影响，所以扑翼飞行器机载设备拍摄的视频抖动剧烈。以平行于地面的直线线路为目标线路的飞行过程为例，由于扑翼飞行器上扑行程和下扑行程不对称，所以扑翼飞行器总是围绕着目标线路上下波动，即扑翼飞行器在y轴方向上的抖动更加剧烈。但是目前几乎没有专门针对扑翼飞行器航拍视频而设计电子稳像算法，相关实验研究也很少。
目前主流的稳像方式有三种：机械式稳像、光学稳像、电子稳像。机械稳像体积大、重量大，根本不适用于低负载的扑翼飞行器，光学稳像效果好但是成本昂贵，只能满足特定方向上的防抖，电子稳像技术以其体积小、重量轻、能耗低、价格低、可移植性高、维护简单和易于升级等众多优点成为小型飞行器的首选，对于扑翼飞行器更是如此。
发明专利202010328381.1公开了一种随机抖动干扰下视频图像的实时电子稳像方法，利用surf算法检测图像的特征点，但是该算法的实时性不高，并且对于边缘光滑目标的特征点提取能力比较弱。在完成放射变换后，并没有对出现的未定义空白区域进行修复，影响稳像效果。发明专利201710563620.x公开了一种基于改进的klt和卡尔曼滤波的电子稳像方法，但并没有给出图像修复的具体方法，当扑翼飞行器飞行时产生的抖动比较剧烈时，仅仅靠单纯的边缘裁剪以完成对图像未定义区的修复达到的效果是不尽如人意的。发明专利201810921737.5公开了一种水下机器人实时电子稳像方法结合imu单元检测出的姿态信息，来计算图像帧间的运动矢量，对视频进行稳像处理，但是机器人的姿态信息并不能完全代表图像帧间的信息，存在一定误差，此外该方法对imu的精度要求极高，稳像的质量不能得到保证。
技术实现要素：
针对上述状况，本发明提出一种结合光流法和savitzky-golay滤波与kalman滤波的实时电子稳像算法，光流法可以跟踪匹配视频帧序列参考帧和目标帧之间的特征角点，并计算出这些特征点的运动矢量；savitzky-golay滤波可以过滤掉由于扑翼飞行器随机抖动产生的运动矢量，也可看成噪声，保留主动运动矢量；利用加速度传感器，确定扑翼飞行器的飞行状态，自适应选择固定帧补偿或者相邻帧补偿，以完成对视频序列的精确补偿；最后，利用插值算法对图像进行放大后将边缘未定义部分区域进行裁剪，便能得到相对稳定的视频。
一种基于扑翼飞行系统的实时电子稳像方法，其特征在于，
所述扑翼飞行系统包括：扑翼飞行器样机、微处理器、无线图传模块、加速传感器、摄像头、锂电池模块、远程服务终端；
所述微处理器分别与所述的摄像头、加速传感器、无线图传模块通过有线方式依次连接；所述无线图传模块与所述远程服务终端通过无线方式依次连接；
所述的微处理器、无线图传模块、加速传感器、摄像头、锂电池模块均放置于所述扑翼飞行器样机上；
所述锂电池模块分别与所述的微处理器、无线图传模块、加速传感器、摄像头通过有线方式依次连接；
所述摄像头用于采集视频帧序列，并视频帧序列传输至所述微处理器；所述微处理器将视频帧序列通过无线图传模块无线传输至所述远程服务终端；
所述加速传感器用于采集x轴加速度、y轴加速度、z轴加速度，并将x轴加速度、y轴加速度、z轴加速度传输至所述微处理器；所述微处理器将x轴加速度、y轴加速度、z轴加速度通过无线图传模块无线传输至所述远程服务终端；
所述实时电子稳像方法，其特征在于，包含以下步骤：
步骤1：所述远程服务终端利用harris角点检测算法检测出视频帧序列的参考帧图像中多个特征角点，利用klt算法在视频帧序列的目标帧图像中跟踪视频帧序列的匹配参考帧图像中的特征角点得到跟踪匹配的视频帧序列的参考帧图像中的特征角点，结合跟踪匹配的视频帧序列的参考帧图像中的特征角点构建视频帧序列的目标帧图像与视频帧序列的参考帧图像的刚性变换模型，通过最小二乘法拟合得到最优解，即为视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量；
步骤1所述视频帧序列为：
{f1,f2，…，fm}
其中，m表示视频帧序列中帧图像的数量，fm为视频帧序列的第m帧图像，定义fp为视频帧序列的第p帧参考帧图像，m∈[1,m]，p∈[1,m)，定义fq为视频帧序列的第q帧目标帧图像，q∈(p,m]；
步骤1所述，若视频帧序列的每帧参考帧图像中多个特征角点，表示为：
{(xp,1，yp,1)，(xp,2，yp,2)，…，(xp,n，yp,n)}
n∈[1,n]
其中，n表示视频帧序列的每帧参考帧图像中特征角点的数量，xp,n表示视频帧序列的第p帧参考帧图像中第n个特征角点的x轴的坐标，yp,n表示视频帧序列的第p帧参考帧图像中第n个特征角点的y轴的坐标；
步骤1所述跟踪匹配的参考帧图像中的特征角点，具体为：
若利用klt算法成功在视频帧序列的第q帧目标帧图像中跟踪匹配到视频帧序列的第p帧参考帧图像中的第k个特征角点(xp,k，yp,k)，则返回布尔值1，以及第k个特征角点在视频帧序列的第q帧目标帧图像中的位置坐标(x'q,k,y'q,k)；
若不成功，则只返回布尔值0；
依次对视频帧序列的第p帧参考帧图像中多个特征角点通过利用klt算法进行跟踪匹配，得到视频帧序列的第p帧参考帧图像中多个特征角点在视频帧序列的第q帧目标帧图像中的对应位置坐标为：
{(x'q，1，y'q，1)，(x'q，2，y'q，2)，…，(x'q,r，y'q,r)}，且r≤n；
步骤1所述目标帧图像与参考帧图像的刚性变换模型，表示为：
其中，(xp,k，yp,k)为视频帧序列的第p帧参考帧图像中的第k个特征角点位置坐标，(x'q,k,y'q,k)为(xp,k，yp,k)在视频帧序列的第q帧目标帧图像中的新位置坐标，δx表示(x'q,k,y'q,k)相对于(xp,k，yp,k)在x轴方向上的位移，δy表示(x'q,k,y'q,k)相对于(xp,k，yp,k)在y轴方向上的位移，α表示(x'q,k,y'q,k)相对于(xp,k，yp,k)沿z轴旋转角度α；
进一步，第k个特征角点的运动矢量可表示为：
pq,k＝[δx,δy,α]t,其中pq,k又称为视频帧序列的第q帧目标帧图像中的第k个特征角点的局部运动矢量；
步骤1所述通过最小二乘法拟合得到最优解为：
r对特征点的位置坐标{(xp,1，yp,1)，(x'q，1，y'q，1)}…{(xq,r，yq,r)，(x'q,r，y'q,r)}，结合上述刚性模型就会组成一个超定方程组am＝b，具体如下：
其最优解为：
m0＝(ata)-1atb
其中，最优解m0中的元素δx0、δy0和α0一起组成步骤1所述的视频帧序列的第q帧目标帧图像相对视频帧序列的第p帧参考帧图像的全局运动矢量，表示为：
vq＝[xt,q,yt,q,zt,q]t；
其中，xt,q＝δx0，xt,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向的位移；
yt,q＝δy0，yt,q表示视频帧序列的第q帧目标帧图像相对视频帧序列的第p帧参考帧图像的在y轴方向上的位移；
zt,q＝α0，zt,q表示视频帧序列的第q帧目标帧图像相对视频帧序列的第p帧参考帧图像的沿z轴方向的旋转角度；
重复vq的计算过程，依次可以得到第(p 1)帧、第(p 2)帧、…、第q帧相对于视频帧序列的第p帧参考帧图像的全局运动矢量的集合{v(p 1),v(p 2),……,vq}；
步骤2：所述远程服务终端将视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量进行累加处理，得到帧序列间的运动轨迹，分别利用savitzky-golay滤波算法和kalman滤波算法对帧序列间的运动轨迹进行平滑处理，得到只包含主动运动的帧序列间的运动轨迹，通过计算求得目标帧相对于参考帧只包含主动运动的运动矢量；
步骤2所述将视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量进行累加处理，具体为：
其中：lq表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像的全局运动矢量和，lx,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向上的矢量和，ly,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴方向上的矢量和，lz,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴方向上的矢量和，q∈(p,m]；
步骤2所述帧序列间的运动轨迹，表示为：
{lx，ly，lz}
lx＝{lx,1，lx,2，…，lx,m}
ly＝{ly,1，ly,2，…，ly,m}
lz＝{lz,1，lz,2，…，lz,m}
其中，lx表示m帧序列在x轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，ly表示m帧序列在y轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，lz表示m帧序列在z轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，lx,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向上的矢量和，ly,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴方向上的矢量和，lz,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴方向上的矢量和，q∈(p,m]；
步骤2所述的对帧序列间的运动轨迹进行平滑处理，具体过程为：
首先利用savitzky-golay滤波算法分别对lx和lz进行滤波平滑后，得到只包含主动运动的x轴运动轨迹即l'x、只包含主动运动的z轴运动轨迹即l'z；
然后利用kalman滤波算法对ly进行滤波平滑处理，得到只包含主动运动的y轴运动轨迹即l'y；
步骤2所述只包含主动运动的帧序列间的运动轨迹，表示为：
{l'x，l'y，l'z}
l'x＝{l'x,1、l'x,2…l'x,m}
l'y＝{l'y,1、l'y,2…l'y,m}
l'z＝{l'z,1、l'z,2…l'z,m}
其中，l'x表示只包含主动运动的x轴运动轨迹，l'z表示只包含主动运动的z轴运动轨迹，l'y只包含主动运动的y轴运动轨迹，l'x,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向上的矢量和，l'y,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴方向上的矢量和，l'z,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴方向上的矢量和，q∈(p,m]；
步骤2所述目标帧相对于参考帧只包含主动运动的运动矢量，表示为：
v'q＝[x't,q,y't,q,z't,q]t
x't,q＝l'x,q-l'x,p
y't,q＝l'y,q-l'y,p
z't,q＝l'z,q-l'z,p
其中：v'q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像只包含主动运动的运动矢量，x't,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴上只包含主动运动的运动矢量，y't,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴只包含主动运动的运动矢量；z't,q示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴上只包含主动运动的运动矢量；
步骤3：通过加速度传感器采集x轴加速度、y轴加速度、z轴加速度，分别计算x轴速度、y轴速度、z轴速度，进一步选择对应补偿算法，对视频帧序列中每幅图像依次进行补偿得到稳定的视频帧序列；稳定的视频帧序列中每幅图像利用双线性插值算法进行图像放大得到放大后的视频帧序列，将放大后的视频帧序列中每幅图像边缘的未定义区域裁剪掉，得到去除了未定义区的稳定帧序列；
步骤3所述加速度传感器采集的x轴加速度、y轴加速度、z轴加速度为：
ax(tn)、ay(tn)、az(tn)；
ax(tn)为tn时刻x轴上加速度，ay(tn)为tn时刻y轴上加速度，az(tn)为tn时刻z轴上加速度；
步骤3所述分别计算x轴速度、y轴速度、z轴速度为：
其中：vx(tn)表示tn时刻飞行器在x轴上的瞬时速度，t0为初始时刻，v(t0)表示初始时刻飞行器在x轴上的瞬时速度；同理有，
步骤3所述进一步选择对应补偿算法为：
分别在x轴、y轴和z轴上设定设定运动阈值βx、βy和βz，当条件满足vx(tn)＞βx，vy(tn)＞βy且vz(tn)＞βz时，并且该条件一定时间td内持续满足，认定当前扑翼飞行器的状态为正常飞行运动状态，对应补偿算法选择相邻帧补偿方法；
否则为悬停状态，对应补偿算法为选择固定帧补偿方法；
步骤3所述原视频帧序列{f1,f2，…，fm}经过对应补偿算法得到稳定的帧序列{f′1，f′2，…，f'm}，相邻帧补偿算法中的补偿向量为(x't,i，y't,i，z't,i)，固定帧补偿中的补偿向量为(l'x,i，l'y,i，l'z,i)；
步骤3所述用利用双线性插值算法将稳定的视频帧序列{f′1，f′2，…，f′m}放大s倍后，得到放大后的视频帧序列{sf′1，sf′2，…，sf′m}，其中，一般s的取值为(1～1.2)。
步骤3所述，利用裁剪的方法依次将放大后的视频帧序列{sf′1，sf′2，…，sf′m}中每帧图像的边缘未定义区裁剪掉，得到去除了未定义区的稳定帧序列{f″1，f″2，…，f″m}，其中有帧序列{f″1，f″2，…，f″m}的尺寸与原视频帧序列{f1,f2，…，fm}一致。
本发明的优点在于，有效解决扑翼飞行器飞行过程中由于机身抖动导致所拍摄的视频出现抖动、不连续和畸变等问题。其优点具体表现在：在运动滤波过程中，针对y轴方向抖动大的特点，单独使用kalman滤波算法，得到的轨迹更为平滑；在运动补偿过程中，对于不同的飞行状态，通过加速度传感器的检测信息，判定扑翼飞行器是正常飞行状态，还是悬停拍摄状态，若为正常飞行状态时使用相邻帧补偿算法，若为悬停状态时转换为固定帧补偿，如此一来，稳像质量有很大提高；在运动修复的过程中，针对不同环境，对于未定义区域的修复，根据需要选择直接裁剪法、间接裁剪法或者智能拼接法。整套稳像算法，满足速度最低需为30帧/s的实时性要求，并且本方法不需要额外增加处理视频处理芯片，成本低，可移植性高，便于实现。
附图说明
图1为本发明一种扑翼飞行器实时电子稳像方法的半实物仿真流程示意图
图2为本发明实施例所述的半实物仿真的第99帧图像harris角点检测的示意图
图3为本发明实施例所述的半实物仿真的第99帧与第100帧图像klt角点跟踪示意图
图4为本发明实施例所述的半实物仿真的原视频的第99、100和101帧图像
图5为本发明实施例所述的半实物仿真的经过相邻帧补偿后的第99、100和101帧图像
图6为本发明实施例所述的半实物仿真的经过双线性插值算法修复后的第99、100和101帧图像
图7为本发明实施例所述的半实物仿真稳像前后视频帧序列psnr值对比
图8为本发明实施例所述的半实物仿真稳像前后视频帧序列ssim值对比
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，下面结合附图和具体实施对本发明做进一步描述，在此仅用本发明的适宜性实例说明来解释本发明，但并不作为本发明的限定。
在具体的实施例中，利用抖动平台模拟扑翼飞行器飞行产生的随机抖动，该抖动模拟平台的水平移动台模拟扑翼飞行的水平向前飞行的过程，垂直升降台模拟扑翼飞行器翼展上下扑动时导致机身上下波动过程。通过模拟平台上的无线图传模块将其所拍摄的原视频序列传输到地面的工作站(即：计算机)，之后通过计算机的强大计算能力，利用本发明中提到的实时电子稳像算法，对原视频序列进行处理，最终得到稳定的实时画面，并将其输出到监视器。具体过程可参见图1。
本发明的实施例为一种基于扑翼飞行系统的实时电子稳像方法，其特征在于，
所述扑翼飞行系统包括：扑翼飞行器抖动模拟平台、微处理器、无线图传模块、加速传感器、摄像头、锂电池模块、远程服务终端；
所述微处理器分别与所述的摄像头、加速传感器、无线图传模块通过有线方式依次连接；所述无线图传模块与所述远程服务终端通过无线方式依次连接；
所述的微处理器、无线图传模块、加速传感器、摄像头、锂电池模块均放置于所述扑翼飞行器抖动模拟平台；
所述锂电池模块分别与所述的微处理器、无线图传模块、加速传感器、摄像头通过有线方式依次连接；
所述摄像头用于采集视频帧序列，并视频帧序列传输至所述微处理器；所述微处理器将视频帧序列通过无线图传模块无线传输至所述远程服务终端；
所述加速传感器用于采集x轴加速度、y轴加速度、z轴加速度，并将x轴加速度、y轴加速度、z轴加速度传输至所述微处理器；所述微处理器将x轴加速度、y轴加速度、z轴加速度通过无线图传模块无线传输至所述远程服务终端；
所述扑翼飞行器抖动模拟平台选型为：可移动式升降抖动平台；
所述微处理器选型为：stm3f2405；
所述无线图传模块选型为：ts5823l；
所述加速传感器选型为：bmi088；
所述摄像头选型为：1200tvl；
所述锂电池模块选型为：552131锂电池(300mah)；
所述远程服务终端选型为：intel(r)core(tm)i5-4200mcpu@2.50ghz2.49ghz；
所述实时电子稳像方法，其特征在于，包含以下步骤：
步骤1：所述远程服务终端利用harris角点检测算法检测出视频帧序列的参考帧图像中多个特征角点，如图2所示，图中的小圆点为视频帧序列中第99帧图像检测出的特征角点；利用klt算法在视频帧序列的目标帧图像中跟踪视频帧序列的匹配参考帧图像中的特征角点得到跟踪匹配的视频帧序列的参考帧图像中的特征角点，如图3所示，图中的小竖线为第99帧图像检测出的所有特征角点运动到第100帧的运动轨迹；结合跟踪匹配的视频帧序列的参考帧图像中的特征角点构建视频帧序列的目标帧图像与视频帧序列的参考帧图像的刚性变换模型，通过最小二乘法拟合得到最优解，即为视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量；
步骤1所述视频帧序列为：
{f1,f2，…，fm}
其中，m表示视频帧序列中帧图像的数量，fm为视频帧序列的第m帧图像，定义fp为视频帧序列的第p帧参考帧图像，m∈[1,m]，p∈[1,m)，定义fq为视频帧序列的第q帧目标帧图像，q∈(p,m]，该实施例中p＝1，q＝2，m＝369；
步骤1所述，若视频帧序列的第99帧参考帧图像中120个特征角点坐标，表示为：
{(x1,1，y1,1)，(x1,2，y1,2)，…，(x1,120，y1,120)}
其中，n＝120表示视频帧序列的第1帧参考帧图像中特征角点的数量，x1,1表示视频帧序列的第1帧参考帧图像中第1个特征角点的x轴的坐标，y1,1表示视频帧序列的第1帧参考帧图像中第1个特征角点的y轴的坐标；
步骤1所述跟踪匹配的参考帧图像中的特征角点，具体为：
若利用klt算法成功在视频帧序列的第2帧目标帧图像中跟踪匹配到视频帧序列的第1帧参考帧图像中的第1个特征角点(x1,1，y1,1)，则返回布尔值1，以及第1个特征角点在视频帧序列的第2帧目标帧图像中的位置坐标(x'2,1,y'2,1)；
若不成功，则只返回布尔值0；
依次对视频帧序列的第1帧参考帧图像中多个特征角点通过利用klt算法进行跟踪匹配，得到视频帧序列的第1帧参考帧图像中多个特征角点在视频帧序列的第2帧目标帧图像中的对应位置坐标为：
{(x'2，1，y'2，1)，(x'2，2，y'2，2)，…，(x'2,115，y'2,115)}；
步骤1所述目标帧图像与参考帧图像的刚性变换模型，表示为：
其中，(x1,1，y1,1)为视频帧序列的第1帧参考帧图像中的第1个特征角点位置坐标，(x'2,1,y'2,1)为(x1,1，y1,1)在视频帧序列的第2帧目标帧图像中的新位置坐标，δx表示(x'2,1,y'2,1)相对于(x1,1，y1,1)在x轴方向上的位移，δy表示(x'2,1,y'2,1)相对于(x1,1，y1,1)在y轴方向上的位移，α表示(x'2,1,y'2,1)相对于(x1,1，y1,1)沿z轴旋转角度α；
进一步，视频帧序列的第1帧参考帧图像中的第1个特征角点的运动矢量可表示为：
p2,1＝[δx,δy,α]t,其中p2,1又称为视频帧序列的第2帧目标帧图像中的第1个特征角点的局部运动矢量；
步骤1所述通过最小二乘法拟合得到最优解为：
115对特征点的位置坐标{(x1,1，y1,1)，(x'2，1，y'2，1)}…{(x2,115，y2,115)，(x'2,115，y'2,115)}，结合上述刚性模型就会组成一个超定方程组am＝b，具体如下：
其最优解为：
m0＝(ata)-1atb
其中，最优解m0中的元素δx0、δy0和α0一起组成步骤1所述的视频帧序列的第2帧目标帧图像相对视频帧序列的第1帧参考帧图像的全局运动矢量，表示为：
v2＝[xt,2,yt,2,zt,2]t
其中，xt,2＝δx0，xt,2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在x轴方向的位移；
yt,2＝δy0，yt,2表示视频帧序列的第2帧目标帧图像相对视频帧序列的第1帧参考帧图像的在y轴方向上的位移；
zt,2＝α0，zt,2表示视频帧序列的第2帧目标帧图像相对视频帧序列的第1帧参考帧图像的沿z轴方向的旋转角度；
重复v2的计算过程，依次可以得到第1帧、第2帧、…、第369帧相对于视频帧序列的第1帧参考帧图像的全局运动矢量的集合{v1,v2,……,v369}，其中v1＝[0,0,0]t；
步骤2：所述远程服务终端将视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量进行累加处理，得到帧序列间的运动轨迹，分别利用savitzky-golay滤波算法和kalman滤波算法对帧序列间的运动轨迹进行平滑处理，得到只包含主动运动的帧序列间的运动轨迹，通过计算求得目标帧相对于参考帧只包含主动运动的运动矢量；
步骤2所述将视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量进行累加处理，具体为：
其中：l2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像的全局运动矢量和，lx,2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在x轴方向上的矢量和，ly,2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在y轴方向上的矢量和，lz,2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在z轴方向上的矢量和；
重复l2的计算过程可以得到原视频中每一帧相对于第1帧参考帧图像的全局运动矢量和的集合{l1,l2,…,l369}；
步骤2所述帧序列间的运动轨迹，表示为：
{lx，ly，lz}
lx＝{lx,1，lx,2，…，lx,369}
ly＝{ly,1，ly,2，…，ly,369}
lz＝{lz,1，lz,2，…，lz,369}
其中，lx表示369帧序列在x轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，ly表示369帧序列在y轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，lz表示369帧序列在z轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，lx,369表示视频帧序列的第369帧目标帧图像相对于视频帧序列的第1帧参考帧图像在x轴方向上的矢量和，ly,369表示视频帧序列的第369帧目标帧图像相对于视频帧序列的第1帧参考帧图像在y轴方向上的矢量和，lz,369表示视频帧序列的第369帧目标帧图像相对于视频帧序列的第1帧参考帧图像在z轴方向上的矢量和；
步骤2所述的对帧序列间的运动轨迹进行平滑处理，具体过程为：
首先利用savitzky-golay滤波算法分别对lx和lz进行滤波平滑后，得到只包含主动运动的x轴运动轨迹即l'x、只包含主动运动的z轴运动轨迹即l'z；
然后利用kalman滤波算法对ly进行滤波平滑处理，得到只包含主动运动的y轴运动轨迹即l'y；
步骤2所述只包含主动运动的帧序列间的运动轨迹，表示为：
{l'x，l'y，l'z}
l'x＝{l'x,1、l'x,2…l'x,369}
l'y＝{l'y,1、l'y,2…l'y,369}
l'z＝{l'z,1、l'z,2…l'z,369}
其中，l'x表示只包含主动运动的x轴运动轨迹，l'y只包含主动运动的y轴运动轨迹，l'z表示只包含主动运动的z轴运动轨迹，l'x,369表示视频帧序列的第369帧目标帧图像相对于视频帧序列的第1帧参考帧图像在x轴方向上的矢量和，l'y,369表示视频帧序列的第369帧目标帧图像相对于视频帧序列的第1帧参考帧图像在y轴方向上的矢量和，l'z,369表示视频帧序列的第369帧目标帧图像相对于视频帧序列的第1帧参考帧图像在z轴方向上的矢量和；
步骤2所述目标帧相对于参考帧只包含主动运动的运动矢量，表示为：
v'2＝[x't,2,y't,2,z't,2]t
x't,2＝l'x,2-l'x,1
y't,2＝l'y,2-l'y,1
z't,2＝l'z,2-l'z,1
其中：v'2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像只包含主动运动的运动矢量，x't,2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在x轴上只包含主动运动的运动矢量，y't,2表示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在y轴只包含主动运动的运动矢量；z't,2示视频帧序列的第2帧目标帧图像相对于视频帧序列的第1帧参考帧图像在z轴上只包含主动运动的运动矢量；
重复v'2的计算过程，可以得到原视频每一帧图像对于第1帧参考帧图像只包含主动运功的运动矢量集合{v'1,v'2,…,v'369}；
步骤3：通过加速度传感器采集x轴加速度、y轴加速度、z轴加速度，分别计算x轴速度、y轴速度、z轴速度，进一步选择对应补偿算法，对视频帧序列中每幅图像依次进行补偿得到稳定的视频帧序列；稳定的视频帧序列中每幅图像利用双线性插值算法进行图像放大得到放大后的视频帧序列，将放大后的视频帧序列中每幅图像边缘的未定义区域裁剪掉，得到去除了未定义区的稳定帧序列；
步骤3所述加速度传感器采集的x轴加速度、y轴加速度、z轴加速度为：
ax(tn)、ay(tn)、az(tn)；
ax(tn)为tn时刻x轴上加速度，ay(tn)为tn时刻y轴上加速度，az(tn)为tn时刻z轴上加速度；
步骤3所述分别计算x轴速度、y轴速度、z轴速度为：
其中：vx(tn)表示tn时刻飞行器在x轴上的瞬时速度，t0为初始时刻，v(t0)表示初始时刻飞行器在x轴上的瞬时速度，这里有v(t0)＝0；同理有，
步骤3所述进一步选择对应补偿算法为：
本实施例中，分别在x轴、y轴和z轴上设定设定运动阈值βx＝0.3、βy＝1.5和βz＝0.4，条件满足vx(tn)＞0.3，vy(tn)＞1.5且vz(tn)＞0.4时，并且该条件在时间td＝5s内持续满足，由此认定扑翼飞行器的状态为正常飞行运动状态，对应补偿算法选择相邻帧补偿方法；
步骤3所述原视频帧序列{f1,f2，…，f369}经过对应补偿算法得到稳定的帧序列{f1'，f2'，…，f'369}，如图4所示为原视频序列的第99帧、第100帧和第101帧图像，如图5所示为经补偿算法得到的第99帧、第100帧和第101帧稳定的图像。利用相邻帧补偿算法计算过程中，原视频的每一帧图像的补偿向量依次对应集合{v'1,v'2,…,v'369}中的每一个元素；
步骤3所述用利用双线性插值算法将稳定的视频帧序列{f'1，f'2，…，f'369}放大s倍后，得到放大后的视频帧序列{sf'1，sf'2，…，sf'369}，其中，一般s的取值为(1～1.2)，这里s＝1.08；
步骤3所述，利用裁剪的方法依次将放大后的视频帧序列{sf′1，sf′2，…，sf′369}中每帧图像的边缘未定义区裁剪掉，得到去除了未定义区的稳定帧序列{f″1，f″2，…，f″369}，其中帧序列{f″1，f″2，…，f″369}的尺寸与原视频帧序列{f1,f2，…，f369}一致，每帧的大小均为640*360，如图6所示，为去除了未定义区的第99帧、第100帧和第101帧图像，相比于图5未定义区域明显减少，视觉效果更好；
本例中稳像后的视频平均峰值信噪比psnr提高了4.45％，具体如图7所示；平均结构相似性ssim提高了30％，具体如图8所示；该稳像算法的平均处理速度达到32帧/秒，满足实时性要求；
以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，根据本发明所述系统及其实施方法所做的同等变化，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。
1.一种基于扑翼飞行系统的实时电子稳像方法，其特征在于，
所述扑翼飞行系统包括：扑翼飞行器样机、微处理器、无线图传模块、加速传感器、摄像头、锂电池模块、远程服务终端；
所述微处理器分别与所述的摄像头、加速传感器、无线图传模块通过有线方式依次连接；所述无线图传模块与所述远程服务终端通过无线方式依次连接；
所述的微处理器、无线图传模块、加速传感器、摄像头、锂电池模块均放置于所述扑翼飞行器样机上；
所述锂电池模块分别与所述的微处理器、无线图传模块、加速传感器、摄像头通过有线方式依次连接；
所述摄像头用于采集视频帧序列，并视频帧序列传输至所述微处理器；所述微处理器将视频帧序列通过无线图传模块无线传输至所述远程服务终端；
所述加速传感器用于采集x轴加速度、y轴加速度、z轴加速度，并将x轴加速度、y轴加速度、z轴加速度传输至所述微处理器；所述微处理器将x轴加速度、y轴加速度、z轴加速度通过无线图传模块无线传输至所述远程服务终端；
所述实时电子稳像方法包含以下步骤：
步骤1：所述远程服务终端利用harris角点检测算法检测出视频帧序列的参考帧图像中多个特征角点，利用klt算法在视频帧序列的目标帧图像中跟踪视频帧序列的匹配参考帧图像中的特征角点得到跟踪匹配的视频帧序列的参考帧图像中的特征角点，结合跟踪匹配的视频帧序列的参考帧图像中的特征角点构建视频帧序列的目标帧图像与视频帧序列的参考帧图像的刚性变换模型，通过最小二乘法拟合得到最优解，即为视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量；
步骤2：所述远程服务终端将视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量进行累加处理，得到帧序列间的运动轨迹，分别利用savitzky-golay滤波算法和kalman滤波算法对帧序列间的运动轨迹进行平滑处理，得到只包含主动运动的帧序列间的运动轨迹，通过计算求得目标帧相对于参考帧只包含主动运动的运动矢量；
步骤3：通过加速度传感器采集x轴加速度、y轴加速度、z轴加速度，分别计算x轴速度、y轴速度、z轴速度，进一步选择对应补偿算法，对视频帧序列中每幅图像依次进行补偿得到稳定的视频帧序列；稳定的视频帧序列中每幅图像利用双线性插值算法进行图像放大得到放大后的视频帧序列，将放大后的视频帧序列中每幅图像边缘的未定义区域裁剪掉，得到去除了未定义区的稳定帧序列。
2.根据权利要求1所述的基于扑翼飞行系统的实时电子稳像方法，其特征在于，
步骤1所述视频帧序列为：
{f1,f2，…，fm}
其中，m表示视频帧序列中帧图像的数量，fm为视频帧序列的第m帧图像，定义fp为视频帧序列的第p帧参考帧图像，m∈[1,m]，p∈[1,m)，定义fq为视频帧序列的第q帧目标帧图像，q∈(p,m]；
步骤1所述，若视频帧序列的第p帧参考帧图像中多个特征角点，表示为：
{(xp,1，yp,1)，(xp,2，yp,2)，…，(xp,n，yp,n)}
n∈[1,n]
其中，n表示视频帧序列的第p帧参考帧图像中特征角点的数量，xp,n表示视频帧序列的第p帧参考帧图像中第n个特征角点的x轴的坐标，yp,n表示视频帧序列的第p帧参考帧图像中第n个特征角点的y轴的坐标；
步骤1所述跟踪匹配的参考帧图像中的特征角点，具体为：
若利用klt算法成功在视频帧序列的第q帧目标帧图像中跟踪匹配到视频帧序列的第p帧参考帧图像中的第k个特征角点(xp,k，yp,k)，则返回布尔值1，以及第k个特征角点在视频帧序列的第q帧目标帧图像中的位置坐标(x′q,k,y′q,k)；
若不成功，则只返回布尔值0；
依次对视频帧序列的第p帧参考帧图像中多个特征角点通过利用klt算法进行跟踪匹配，得到视频帧序列的第p帧参考帧图像中多个特征角点在视频帧序列的第q帧目标帧图像中的对应位置坐标为：
{(x′q，1，y′q，1)，(x′q，2，y′q，2)，…，(x′q,r，y′q,r)}，且r≤n；
步骤1所述目标帧图像与参考帧图像的刚性变换模型，表示为：
其中，(xp,k，yp,k)为视频帧序列的第p帧参考帧图像中的第k个特征角点位置坐标，(x′q,k,y′q,k)为(xp,k，yp,k)在视频帧序列的第q帧目标帧图像中的新位置坐标，δx表示(x′q,k,y′q,k)相对于(xp,k，yp,k)在x轴方向上的位移，δy表示(x′q,k,y′q,k)相对于(xp,k，yp,k)在y轴方向上的位移，α表示(x′q,k,y′q,k)相对于(xp,k，yp,k)沿z轴旋转角度α；
进一步，第k个特征角点的运动矢量可表示为：
pq,k＝[δx,δy,α]t,其中pq,k又称为视频帧序列的第q帧目标帧图像中的第k个特征角点的局部运动矢量；
步骤1所述通过最小二乘法拟合得到最优解为：
r对特征点的位置坐标{(xp,1，yp,1)，(x′q，1，y′q，1)}…{(xq,r，yq,r)，(x′q,r，y′q,r)}，结合上述刚性模型就会组成一个超定方程组am＝b，具体如下：
其最优解为：
m0＝(ata)-1atb
其中，最优解m0中的元素δx0、δy0和α0一起组成步骤1所述的视频帧序列的第q帧目标帧图像相对视频帧序列的第p帧参考帧图像的全局运动矢量，表示为：
vq＝[xt,q,yt,q,zt,q]t；
其中，xt,q＝δx0，xt,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向的位移；
yt,q＝δy0，yt,q表示视频帧序列的第q帧目标帧图像相对视频帧序列的第p帧参考帧图像的在y轴方向上的位移；
zt,q＝α0，zt,q表示视频帧序列的第q帧目标帧图像相对视频帧序列的第p帧参考帧图像的沿z轴方向的旋转角度；
重复vq的计算过程，依次可以得到第(p 1)帧、第(p 2)帧、…、第q帧相对于视频帧序列的第p帧参考帧图像的全局运动矢量的集合{v(p 1),v(p 2),……,vq}。
3.根据权利要求1所述的基于扑翼飞行系统的实时电子稳像方法，其特征在于，
步骤2所述将视频帧序列的参考帧图像与视频帧序列的目标帧图像之间的全局运动矢量进行累加处理，具体为：
其中：lq表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像的全局运动矢量和，lx,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向上的矢量和，ly,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴方向上的矢量和，lz,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴方向上的矢量和，q∈(p,m]；
步骤2所述帧序列间的运动轨迹，表示为：
{lx，ly，lz}
lx＝{lx,1，lx,2，…，lx,m}
ly＝{ly,1，ly,2，…，ly,m}
lz＝{lz,1，lz,2，…，lz,m}
其中，lx表示m帧序列在x轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，ly表示m帧序列在y轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，lz表示m帧序列在z轴上的一条既包含主动运动又包含随机抖动运动的运动轨迹，lx,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向上的矢量和，ly,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴方向上的矢量和，lz,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴方向上的矢量和，q∈(p,m]；
步骤2所述的对帧序列间的运动轨迹进行平滑处理，具体过程为：
首先利用savitzky-golay滤波算法分别对lx和lz进行滤波平滑后，得到只包含主动运动的x轴运动轨迹即l′x、只包含主动运动的z轴运动轨迹即l′z；
然后利用kalman滤波算法对ly进行滤波平滑处理，得到只包含主动运动的y轴运动轨迹即l′y；
步骤2所述只包含主动运动的帧序列间的运动轨迹，表示为：
{l′x，l′y，l′z}
l′x＝{l′x,1、l′x,2…l′x,m}
l′y＝{l′y,1、l′y,2…l′y,m}
l′z＝{l′z,1、l′z,2…l′z,m}
其中，l′x表示只包含主动运动的x轴运动轨迹，l′z表示只包含主动运动的z轴运动轨迹，l′y只包含主动运动的y轴运动轨迹，l′x,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴方向上的矢量和，l′y,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴方向上的矢量和，l′z,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴方向上的矢量和，q∈(p,m]；
步骤2所述目标帧相对于参考帧只包含主动运动的运动矢量，表示为：
v′q＝[x′t,q,y′t,q,z′t,q]t
x′t,q＝l′x,q-l′x,p
y′t,q＝l′y,q-l′y,p
z′t,q＝l′z,q-l′z,p
其中：v′q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像只包含主动运动的运动矢量，x′t,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在x轴上只包含主动运动的运动矢量，y′t,q表示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在y轴只包含主动运动的运动矢量；z′t,q示视频帧序列的第q帧目标帧图像相对于视频帧序列的第p帧参考帧图像在z轴上只包含主动运动的运动矢量。
4.根据权利要求1所述的基于扑翼飞行系统的实时电子稳像方法，其特征在于，
步骤3所述加速度传感器采集的x轴加速度、y轴加速度、z轴加速度为：
ax(tn)、ay(tn)、az(tn)；
ax(tn)为tn时刻x轴上加速度，ay(tn)为tn时刻y轴上加速度，az(tn)为tn时刻z轴上加速度；
步骤3所述分别计算x轴速度、y轴速度、z轴速度为：
其中：vx(tn)表示tn时刻飞行器在x轴上的瞬时速度，t0为初始时刻，v(t0)表示初始时刻飞行器在x轴上的瞬时速度；同理有，
步骤3所述进一步选择对应补偿算法为：
分别在x轴、y轴和z轴上设定设定运动阈值βx、βy和βz，当条件满足vx(tn)＞βx，vy(tn)＞βy且vz(tn)＞βz时，并且该条件一定时间td内持续满足，认定当前扑翼飞行器的状态为正常飞行运动状态，对应补偿算法选择相邻帧补偿方法；
否则为悬停状态，对应补偿算法为选择固定帧补偿方法；
步骤3所述原视频帧序列{f1,f2，…，fm}经过对应补偿算法得到稳定的帧序列{f′1，f′2，…，f′m}，相邻帧补偿算法中的补偿向量为(x′t,i，y′t,i，z′t,i)，固定帧补偿中的补偿向量为(l′x,i，l′y,i，l′z,i)；
步骤3所述用利用双线性插值算法将稳定的视频帧序列{f′1，f′2，…，f′m}放大s倍后，得到放大后的视频帧序列{sf′1，sf′2，…，sf′m}，其中，一般s的取值为(1～1.2)；
步骤3所述，利用裁剪的方法依次将放大后的视频帧序列{sf′1，sf′2，…，sf′m}中每帧图像的边缘未定义区裁剪掉，得到去除了未定义区的稳定帧序列{f″1，f″2，…，f″m}，其中有帧序列{f″1，f″2，…，f″m}的尺寸与原视频帧序列{f1,f2，…，fm}一致。
技术总结