CN116132819B - 一种航空光学相机阵列的实时拍照规划方法 - Google Patents

Abstract

本发明公开了一种航空光学相机阵列的实时拍照规划方法，包括：步骤S1、对参与计算的图像边界像素对进行稀疏化；步骤S2、采用二分法和数字高程模型迭代求解图像边界像素在辅助坐标系中的图像边界物点集合；步骤S3、根据要求的重叠率和图像边界物点集合，计算重叠区域边界物点集合；步骤S4、在重叠区域边界物点集合中，选择与视场后平面的航向距离最小的物点作为重叠区域临界物点；步骤S5、根据重叠区域临界物点，计算下一帧曝光平面的方程。本发明有效考虑了相机阵列的安装方式和成像模式，能更好的适应地形起伏的地貌、变化的载机速度，提高了航向重叠率的完成精度，有效提高了侦察、测绘任务对外部条件的适应能力。

Description

一种航空光学相机阵列的实时拍照规划方法

技术领域

本发明涉及种航空光学相机拍照技术领域，具体为一种航空光学相机阵列的实时拍照规划方法。

背景技术

航空光学相机是航空遥感中获得地物信息的一种重要载荷。为了提高侦察效率，可以将N个相机进行结构上的刚性连接，组成航空相机阵列系统，以视场拼接的方式获得更大的视场角。

航空相机阵列的拍照规划有两种方式，一种是提前进行任务规划，即使用地面站的任务规划软件，在载机的飞行路径上设定若干曝光点，当载机飞到曝光点位置即执行一次拍照。第二种是实时拍照任务规划，即在载机飞行的过程中，实时规划曝光点，这种方法更灵活、适应性更强，不需要和飞行计划绑定，并且可以避免飞机位置偏离预设曝光点导致的漏拍问题。

现阶段常用的实时拍照规划方法为根据拍照间隔执行拍照。该方法无法适应相机阵列中各相机的安装特点和成像方式、无法适应变化的地形和变化的速度，容易导致图像重叠率不准确、航拍图像拉缝的问题，为此，我们提出一种航空光学相机阵列的实时拍照规划方法。

发明内容

本发明的目的在于提供一种航空光学相机阵列的实时拍照规划方法，解决了现有的问题。

为实现上述目的，本发明提供如下技术方案：一种航空光学相机阵列的实时拍照规划方法，包括：

步骤S1、对参与计算的图像边界像素对进行稀疏化；

步骤S2、采用二分法和数字高程模型迭代求解图像边界像素在辅助坐标系中的图像边界物点集合；

步骤S3、根据要求的重叠率和图像边界物点集合，计算重叠区域边界物点集合；

步骤S4、在重叠区域边界物点集合中，选择与视场后平面的航向距离最小的点作为重叠区域临界物点；

步骤S5、根据重叠区域临界物点，计算下一帧曝光平面的方程；

步骤S6、将POS新的地理坐标转换到辅助坐标系，代入到曝光平面方程中，判断是否执行曝光。

优选的，步骤S1中具体包括以下步骤:

步骤S11、针对相机阵列中的每个相机，计算其第i个边缘像素点对应的横向GSD：

其中，θ_fov是每个像素的瞬时视场角：

θ_cross为该相机在横向上偏离导航坐标系中Z轴矢量的角度：

θ_cross＝θ_roll+θ_dev

其中θ_dev为当前相机偏离相机阵列坐标系Z轴矢量的角度，θ_roll为相机阵列坐标系的横滚角；

步骤S12、根据每个像素的横向GSD，计算出稀疏的图像边界像素对集合P^c。

优选的，步骤S2中具体包括以下步骤:

步骤S21、将数字高程模型中的坐标从地理坐标系转换为辅助坐标系，其中，辅助坐标系选择以当前曝光位置为原点的导航坐标系；

步骤S211、将数字高程模型中的坐标从地理坐标系转换到地心地固坐标系；

从地理坐标系到地心地固坐标系的转换算子为：

X^e＝(R_E(B)+H)cos(B)cos(L)

Y^e＝(R_E(B)+H)cos(B)sin(L)

Z^e＝[(1-e²)R_E(B)+H]sinB

式中，

L表示经度，B表示纬度，H表示高度。

R₀＝6378137m

e＝0.0818191908425

将DEM中所有点的地理坐标(L_DEM,B_DEM,H_DEM)和辅助坐标系原点的地理坐标(L_S,B_S,H_S)通过转换，分别得到和

步骤S212、将从地心地固坐标系转换到辅助坐标系。从地心地固坐标系到辅助坐标系的转换公式为：

式中，为从地心地固坐标系到辅助坐标系的旋转矩阵，具体如下：

步骤S22、根据DEM中的高度边界得到待求高度的初值

其中，

和分别为DEM所有中的Z^DEM的最大值和最小值；

步骤S23、通过视轴方程，计算P^c中每个元素(x,y)对应的地物的平面坐标(X^j,Y^j)；

其中，输入的参数包括：P^c中的坐标(x,y)，定位定姿系统获得的姿态数据(航向角ψ，俯仰角θ，横滚角φ)，相机阵列坐标系到导航坐标系的旋转矩阵本次迭代(第j轮)中的高度值Z^j；

计算式如下：

其中，(X_I,Y_I,Z_I)通过下式得到：

其中，为POS本体坐标系到导航坐标系的旋转矩阵：

x为(x,y)在相机坐标系中的坐标：

x＝(x-x₀,y-y₀,f)^T

其中(x₀,y₀)是主点、f是主距；

(x,y)为集合P中的元素；

步骤S24、根据计算得到的在DEM中找到对应的栅格，得到对应的Z^DEM；

步骤S25、用二分法迭代求解P^c对应的图像边界物点集合P^aux。

优选的，步骤S3具体包括以下步骤：

步骤S31：根据要求的重叠率和图像边界物点集合，得到重叠区域边界点的平面坐标：

其中，ρ是重叠率。

步骤S32：将重叠区域边界点的平面坐标代入到DEM中，得到第j个相机对应的重叠区域边界点集合C^j。

优选的，步骤S4具体包括以下步骤：

步骤S41、对于第j个相机，计算C^j中每个点到该相机视场后平面Π的航向交点D；

交点D的坐标是下面方程的解：

其中，(n_Πx,n_Πy,n_Πz)^T为相机视场后平面Π的法向量

其中，n是飞行方向上的像素个数，b是像素尺寸，f是焦距。

步骤S42、计算所有相机中每个C_iD_i向量的模，并从中选出模最小对应的C点作为相机阵列的重叠区域临界物点。

优选的，步骤S5具体包括以下步骤：

定义曝光平面H为经过重叠区域临界物点的平面，其法向量为重叠区域临界物点所属相机的视场后平面法向量

H平面的表达式如下：

n_Πx(X-X_C)+n_Πy(Y-Y_C)+n_Πz(Z-Z_C)＝0

其中，(X_C,Y_C,Z_C)^T为重叠区域临界物点在辅助坐标系中的坐标。

优选的，步骤S6具体包括以下步骤：

将POS获得新的地理坐标(L_t,B_t,H_t)转换到辅助坐标系：

其中，(L₀,B₀,H₀)为辅助坐标系原点的地理坐标；

将X_t代入到曝光平面H的方程中，如果满足

|n_Πx(X_t-X_C)+n_Πy(Y_t-Y_C)+n_Πz(Z_t-Z_C)|＜δ，则执行曝光，其中δ为阈值。

与现有技术相比，本发明的有益效果如下：

本发明引入了数字高程模型(DEM)参与实时拍照规划的计算，避免了传统方法无法适应地表高程变化的问题；

本发明通过计算下一帧的曝光平面，避免了传统方法无法适应速度矢量变化的问题。

本发明通过计算相机阵列中与视场后平面航向距离最小的重叠区域临界点，避免了传统方法无法适应相机的倾斜安装方式、倾斜成像方式的问题。

本发明通过稀疏化相机阵列的边界像素点，加速了计算过程；

本发明通过二分搜索，加速了求解边界像素点地物坐标的计算过程。

附图说明

图1为本发明流程示意图；

图2为本发明边界像素示意图；

图3为本发明相机坐标系示意图；

图4为本发明坐标计算过程示意图；

图5为本发明地理坐标系转换为辅助坐标过程示意图；

图6为本发明数字高程模型示意图；

图7为本发明二分法求解物点坐标系1示意图；

图8为本发明二分法求解物点坐标系2示意图；

图9为本发明相机照片实物图1；

图10为本发明重叠区域边界点坐标示意图；

图11为本发明相机照片实物图2；

图12为本发明航向覆盖区域示意图；

图13为本发明最短距离计算步骤示意图；

图14为本发明相机照片实物图3；

图15为本发明H平面示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。

如图1所示，一种航空光学相机阵列的实时拍照规划方法，包括：

步骤S1、因为数字高程模型的分辨率通常小于光电传感器的地面采样距离(GSD)，因此首先对参与计算的边界像素对进行稀疏化，其中，“边界像素点”、“边界像素对”和“航向”、“横向”的概念见附图2所示；

步骤S11、以中相机坐标系作为相机阵列坐标系，原点为中相机的投影中心，轴系朝向和导航坐标系相同；

见附图3所示，针对相机阵列中的每个相机，计算其第i个边缘像素点对应的横向GSD：

其中，θ_fov是每个像素的瞬时视场角：

θ_cross为该相机在横向上偏离导航坐标系中Z轴矢量的角度：

θ_cross＝θ_roll+θ_dev

其中θ_dev为当前相机偏离相机阵列坐标系Z轴矢量的角度，θ_roll为相机阵列坐标系的横滚角；

步骤S12、计算出稀疏的边缘像素对在横向上的索引集合；

具体算法过程如下：

通过索引集合B，可以得到稀疏的边缘像素对的集合

P^c＝{(x_i,y₁)，(x_i,y_m)},i∈Β。

步骤S2、遍历集合P^c中每个像素，用二分法计算该像素对应的物点在辅助坐标系中的坐标，得到集合其中，具体计算过程如附图4所示；

步骤S21、见附图5所示，将数字高程模型中的坐标从地理坐标系转换为辅助坐标系，其中，辅助坐标系选择以当前曝光位置为原点的导航坐标系；

步骤S211、从地理坐标系到地心地固坐标系的转换算子为：

X^e＝(R_E(B)+H)cos(B)cos(L)

Y^e＝(R_E(B)+H)cos(B)sin(L)

Z^e＝[(1-e²)R_E(B)+H]sinB

式中，

L表示经度，B表示纬度，H表示高度。

R₀＝6378137m

e＝0.0818191908425

将DEM某点的地理坐标(L_DEM,B_DEM,H_DEM)和辅助坐标系原点的地理坐标(L_S,B_S,H_S)通过转换，分别得到和

步骤S212、从地心地固坐标系到辅助坐标系的转换公式为：

式中，为从地心地固坐标系到辅助坐标系的旋转矩阵，具体如下：

步骤S22、见附图6所示，根据DEM中的高度边界得到待求高度的初值

其中，

和分别为X^DEM中Z^DEM的最大值和最小值；

步骤S23、通过视轴方程，计算P^c中每个元素(x,y)对应的地物的平面坐标(X^j,Y^j)；

其中，输入的参数包括：P^c中的坐标(x,y)，定位定姿系统获得的姿态数据(航向角ψ，俯仰角θ，横滚角φ)，相机阵列坐标系到导航坐标系的旋转矩阵本次迭代(第j轮)中的高度值Zj；

计算式如下：

其中，(X_I,Y_I,Z_I)通过下式得到：

其中，为POS本体坐标系到导航坐标系的旋转矩阵：

x为(x,y)在相机坐标系中的坐标：

x＝(x-x₀,y-y₀,f)^T

其中(x₀,y₀)是主点、f是主距；

(x,y)为集合P中的元素；

步骤S24、根据计算得到的在DEM中找到对应的栅格，得到对应的Z^DEM；

步骤S25、用二分法迭代求解物点坐标，见附图7所示，设阈值为δ，如果|Z^j-Z^DEM|＜δ，则计算结束，当前像素(x,y)对应的物点坐标为(X^j,Y^j,Z^j)；如果则

见附图8所示，如果则

遍历P^c中的所有像素后，最终得到P^c对应的所有物点在辅助坐标系中的集合

以中相机为例，P^aux见附图9所示。

步骤S3、根据要求的重叠率和集合P^aux，计算重叠区域边界点坐标集合

见附图10所示，其中BC曲线表示的地表是时间上相邻两张图像的重叠区域，而C是重叠区域的边界点；

步骤S31、根据P^e，得到P^aux在辅助坐标系下的Z＝0平面的投影坐标集合

步骤S32、根据重叠率，得到重叠区域边界点的平面坐标：

步骤S33、见附图11所示，将代入到DEM中，得到重叠区域边界点的三维坐标：对于相机阵列中第j个相机，它的所有C_i组成集合C^j。

步骤S4、选择与视场后平面的航向距离最小的点；

只要下一张曝光的图像能够以重叠率ρ覆盖该点，那么该图像在横向上的任意位置都可以满足重叠率ρ；

见附图12所示，C_m与视场后平面的航向距离最小，这样只要下一张图像覆盖了C_m、就一定能覆盖C_n等其他点；

见附图13所示，具体步骤如下：首先，对每个相机，计算其点到视场后平面的航向交点D，然后，在所有相机中，选择最短CD对应的C；

步骤S41、对于第j个相机，计算C^j中所有点到该相机视场后平面Π的航向交点D；

交点D的坐标是下面方程的解：

其中，(n_Πx,n_Πy,n_Πz)^T为相机视场后平面Π的法向量

步骤S42、见附图14所示，计算所有相机的每个C_iD_i向量的模，并从中选出模最小对应的C点。

步骤S5、计算下一帧曝光平面H的方程；

见附图15所示，定义曝光平面H为经过C点的平面，其法向量为C点所属相机的视场后平面法向量

设C点在辅助坐标系中的坐标为(X_C,Y_C,Z_C)^T，则H平面的表达式如下：

n_Πx(X-X_C)+n_Πy(Y-Y_C)+n_Πz(Z-Z_C)＝0。

步骤S6、将POS新的地理坐标转换到辅助坐标系；

当POS坐标位置在H平面上时，执行曝光；

设POS获得新的地理坐标为(L_t,B_t,H_t)，将其转换到辅助坐标系：

其中，(L₀,B₀,H₀)为辅助坐标系原点的地理坐标；

将X_t代入到曝光平面H的方程中，如果满足

|n_Πx(X_t-X_C)+n_Πy(Y_t-Y_C)+n_Πz(Z_t-Z_C)|＜δ，则执行曝光，其中δ为阈值。

尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同物限定。

Claims (4)

1.一种航空光学相机阵列的实时拍照规划方法，其特征在于，包括：

步骤S1、对参与计算的图像边界像素对进行稀疏化；

步骤S2、采用二分法和数字高程模型迭代求解图像边界像素在辅助坐标系中的图像边界物点集合；

步骤S21、将数字高程模型中的坐标从地理坐标系转换为辅助坐标系，其中，辅助坐标系选择以当前曝光位置为原点的导航坐标系；

步骤S211、将数字高程模型中的坐标从地理坐标系转换到地心地固坐标系；

从地理坐标系到地心地固坐标系的转换算子为：

X^e＝(R_E(B)+H)cos(B)cos(L)

Y^e＝(R_E(B)+H)cos(B)sin(L)

Z^e＝[(1-e²)R_E(B)+H]sin B

式中，

L表示经度，B表示纬度，H表示高度；

R₀＝6378137m

e＝0.0818191908425

将DEM中所有点的地理坐标(L_DEM,B_DEM,H_DEM)和辅助坐标系原点的地理坐标(L_S,B_S,H_S)通过转换，分别得到和

步骤S212、将从地心地固坐标系转换到辅助坐标系，从地心地固坐标系到辅助坐标系的转换公式为：

式中，为从地心地固坐标系到辅助坐标系的旋转矩阵，具体如下：

步骤S22、根据DEM中的高度边界得到待求高度的初值；

其中，

和分别为DEM所有中的Z^DEM的最大值和最小值；

步骤S23、通过视轴方程，计算P^c中每个元素(x,y)对应的地物的平面坐标(X^j,Y^j)；

其中，输入的参数包括：P^c中的坐标(x,y)，定位定姿系统获得的姿态数据(航向角ψ，俯仰角θ，横滚角φ)，相机阵列坐标系到导航坐标系的旋转矩阵本次迭代(第j轮)中的高度值Z^j；

计算式如下：

其中，(X_I,Y_I,Z_I)通过下式得到：

其中，为POS本体坐标系到导航坐标系的旋转矩阵：

x为(x,y)在相机坐标系中的坐标：

x＝(x-x₀,y-y₀,f)^T

其中(x₀,y₀)是主点、f是主距；

(x,y)为集合P中的元素；

步骤S24、根据计算得到的在DEM中找到对应的栅格，得到对应的Z^DEM；

步骤S25、用二分法迭代求解P^c对应的图像边界物点集合P^aux；

步骤S3、根据要求的重叠率和图像边界物点集合，计算重叠区域边界物点集合；

步骤S4、在重叠区域边界物点集合中，选择与视场后平面的航向距离最小的点作为重叠区域临界物点；

步骤S41、对于第j个相机，计算C^j中每个点到该相机视场后平面Π的航向交点D；

交点D的坐标是下面方程的解：

其中，(n_Πx,n_Πy,n_Πz)^T为相机视场后平面Π的法向量

其中，n是飞行方向上的像素个数，b是像素尺寸，f是焦距；

步骤S42、计算所有相机中每个C_iD_i向量的模，并从中选出模最小对应的C点作为相机阵列的重叠区域临界物点；

步骤S5、根据重叠区域临界物点，计算下一帧曝光平面的方程；

定义曝光平面Η为经过重叠区域临界物点的平面，其法向量为重叠区域临界物点所属相机的视场后平面法向量

Η平面的表达式如下：

n_Πx(X-X_C)+n_Πy(Y-Y_C)+n_Πz(Z-Z_C)＝0

其中，(X_C,Y_C,Z_C)^T为重叠区域临界物点在辅助坐标系中的坐标；

步骤S6、将POS所获取的新的地理坐标转换到辅助坐标系，代入到曝光平面方程中，判断是否执行曝光。

2.根据权利要求1所述的一种航空光学相机阵列的实时拍照规划方法，其特征在于，步骤S1中具体包括以下步骤:

步骤S11、针对相机阵列中的每个相机，计算其第i个边缘像素点对应的横向GSD：

其中，θ_fov是每个像素的瞬时视场角：

θ_cross为该相机在横向上偏离导航坐标系中Z轴矢量的角度：

θ_cross＝θ_roll+θ_dev

其中θ_dev为当前相机偏离相机阵列坐标系Z轴矢量的角度，θ_roll为相机阵列坐标系的横滚角；

步骤S12、根据每个像素的横向GSD，计算出稀疏的图像边界像素对集合P^c。

3.根据权利要求1所述的一种航空光学相机阵列的实时拍照规划方法，其特征在于，步骤S3具体包括以下步骤：

步骤S31：根据要求的重叠率和图像边界物点集合，得到重叠区域边界点的平面坐标：

其中，ρ是重叠率；

步骤S32：将重叠区域边界点的平面坐标代入到DEM中，得到第j个相机对应的重叠区域边界点集合C^j。

4.根据权利要求1所述的一种航空光学相机阵列的实时拍照规划方法，其特征在于，步骤S6具体包括以下步骤：

将POS获得新的地理坐标(L_t,B_t,H_t)转换到辅助坐标系：

其中，(L₀,B₀,H₀)为辅助坐标系原点的地理坐标；

将X_t代入到曝光平面Η的方程中，如果满足|n_Πx(X_t-X_C)+n_Πy(Y_t-Y_C)+n_Πz(Z_t-Z_C)|<δ，则执行曝光，其中δ为阈值。