光学学报, 2018, 38 (12): 1228001, 网络出版: 2019-05-10   

基于微通道板的单光子激光测高技术研究 下载: 1055次封面文章

Photon-Counting Laser Altimetry Based on Microchannel Plate
作者单位
1 中国科学院西安光学精密机械研究所, 陕西 西安 710119
2 中国科学院大学, 北京 100049
摘要
与传统的激光测高技术相比,单光子激光测高技术具有数据量大、质量轻、测距精度高等优势,是激光测高技术的发展趋势。建立数学模型对单光子激光测高特性进行研究,用数值计算估计了单光子激光测高的性能,并建立了地物模型,用蒙特卡罗方法进行了仿真,提出了测高数据的滤波方法和一种利用遥感图像优化高程信息的算法。结果表明:在正午阳光背景最强烈的条件下,典型地物模型单光子激光测高的均方根误差为6.1 cm,用算法优化后的误差为2.6 cm。
Abstract
Compared with the traditional laser altimetry technology, the photon-counting laser altimetry technology has the advantages of large data, light weight, and high ranging precision, which is the development trend of laser altimeter technolgoy. In this paper, we establish a mathematical model to study the characteristics of the photon-counting laser altimetry. The performance of the photon-counting laser altimeter is estimated by numerical calculation. The ground object model is established, and the simulation is carried out with Monte Carlo method. A filtering method for the altimetry data and an algorithm for optimizing the elevation information using the remote sensing images are proposed. The results show that the root-mean-square error of the photon-counting laser altimeter is 6.1 cm under the condition of the noonday background with the most intense sun for the typical ground model, and the error after optimization by the algorithm is 2.6 cm.

1 引言

激光测高技术是一种先进的测量手段,在激光器出现之初,其作为传统雷达技术与激光技术的结合产物就已被提出来[1]。卫星激光测高技术以卫星为平台,搭载激光器、探测器和计时器,通过测量激光从激光器发射到地面反射之后再被探测器接收的时间,就可以知道激光传输路径的精确距离。再结合精确的卫星轨道和姿态信息,就可以获取地面三维模型。激光测高技术可以被应用于气候研究、调查森林覆盖率、城市规划等。1994年,美国第一次在探月卫星Clementine号上搭载激光测高仪[2-4]来获取高精度的月球表面特征信息。目前,激光测高技术已得到快速发展和广泛应用[5]

美国的Clementine、ICESat-1、LOLA等[3-4,6-7],以及我国的嫦娥一号、嫦娥二号、资源三号02星[8-10]都搭载有激光测高仪,日本、印度等也有成功发射星载激光测高仪的记录[11-13],这表明激光测高仪在遥感观测中具有重要意义。目前发射的激光测高仪都采用工作在线性模式的雪崩光电二极管(APD),可探测的峰值功率在10-9~10-8 W量级,需要激光器的单脉冲能量在10-2~10-1 J量级[14],大脉冲能量制约了激光器脉冲的重复频率,同时也会影响激光器的寿命。激光测高仪的发展趋势是采用单光子探测器,单光子探测器的灵敏度可达10-15~10-14 W量级[14],可以提高脉冲的重复频率,减小足迹间隔,获得更高的空间分辨率。目前,基于卫星平台的单光子激光测高技术还在研究阶段,NASA计划中的ICESat-2和LIST卫星都将搭载单光子激光测高仪[15-16]

单光子激光测高仪是激光测高仪的发展趋势,单光子激光测高技术与以前的技术有很大区别,点云数据的处理是一个难题,目前国内对这方面的研究还很少。本文对单光子激光测高技术进行了数学建模,对误差进行了理论分析和数值计算,对单光子器件死时间带来的系统误差进行了理论分析和数值计算,并建立了地物模型进行仿真计算,提出一种优化高程信息的方法。

2 单光子激光测高技术

2.1 单光子激光测高原理

星载激光测高仪类似于激光测距仪,用飞行时间法(TOF)对地面成像,可以获得地面目标的高程信息[17]

单光子激光测高的原理是用窄脉冲、高重复频率的激光照射地面,然后用探测器读出返回的光子,每个光子可以精确定位三维坐标。

2.2 探测器

单光子探测器主要有两种:微通道板(MCP)和APD[16,18-19]。APD根据其淬火电路的不同又可分为被动淬火APD(PQ-APD)和主动淬火APD(AQ-APD)。

APD是一种高性能的光电传感器,普遍应用在各种场景中。当APD偏置电压超过其雪崩击穿电压时,称为盖革模式APD,可以响应单光子[20]。APD由于特殊的工作原理,在雪崩发生之后不能自行恢复到探测状态,需要配合淬火电路来重置[21]。重置过程中,探测器不能响应到来的光子,导致其死时间较长,由于每个APD单元都有独立的淬火电路,导致APD阵列的集成度不够,APD阵列像素尺寸在0.5~2 mm之间,像素间隙为0.3~0.4 mm,典型的分辨率有4 pixel×4 pixel、8 pixel×8 pixel、16 pixel×16 pixel,像素总数只有几百[22-23]。集成度不够高导致APD阵列体积大且很重,极大地增加了测高系统的设计难度。

MCP是一种二维平面型电子倍增器,由大量平行排布的单通道电子倍增器组成。通道内壁有较高的倍增系数,光子从光阴极入射转化为电子,电子在电场作用下连续碰撞内壁实现多次倍增,从而增强光子图像[24]。通道的孔径通常为6~12 μm,器件像素总数可以达到百万量级[25]

探测器的死时间会影响计数率。在强烈的阳光背景下,背景光子可能会触发探测器,导致探测器进入死时间,无法被信号光子触发,所以信号光子的计数率会下降,用信号光子计数率下降因子Fscr来描述[26]

Fscr=exp(-nsτd),(1)

式中:ns为太阳光子的平均计数率;τd为探测器死时间。取三种探测器的典型性能参数,PQ-APD的死时间为1600 ns,AQ-APD的死时间为50 ns,MCP的死时间为2 ns[26],将它们代入(1)式,可得三种探测器的Fscr,如图1所示。

图 1. 探测器的计数率下降因子

Fig. 1. Signal count reduction factor of different detectors

下载图片 查看所有图片

图1中可以看出:MCP的计数率最高;在背景比较暗的情况下,AQ-APD和MCP的表现比较接近,而MCP即使是在背景亮度很高的情况下也有很好的计数率表现;随着背景光子增多,两种APD的性能急剧下降。

与APD相比,MCP有一个缺点,即提取过量的电荷会降低其寿命。此外,MCP的电流会饱和,在强光刺激下,MCP输出电流接近饱和时,增益会降低,这种饱和机制能在强光环境下保护探测器。

2.3 信号预滤波

对于单光子激光测高仪,有三种预滤波技术可以使用[26],它们分别是光谱滤波、空间滤波和时间滤波。光谱滤波是指选用与激光器波长匹配的窄带滤光片,在保证信号光子透过率的情况下,尽量缩小带宽,从而有效滤除背景噪声。空间滤波是指根据激光器投射在地面的光斑尺寸来设计探测器的光学系统,使视场与光斑匹配,获取尽可能多的信号光子。时间滤波是指探测器受选通门控制,只关心地面高程范围内返回的光子,而不响应云层、大气散射的光子。

空间滤波能增大信号的能量,增加探测到的事件数;时间滤波能有效避免云层和大气散射的干扰;空间滤波和时间滤波能提高信噪比,但不能提高探测器的Fscr。窄带滤光片能有效滤除太阳光中的大部分成分,能提高信噪比(SNR)和Fscr

2.4 探测器的选择

与APD相比,MCP有诸多优势:MCP是二维元件,很容易达到百万量级的像素数(APD想要获取二维信息就必须拼接成APD阵列,其像素数也仅有几百量级);MCP具有高增益、增益可变、电流自饱和、防强光等优点;MCP搭配交叉延迟线(CDL)可以获取光子的(x,y)坐标;计数率能达到106 s-1量级[19],极高的计数率能带来极高的事件数、密集的点云,搭配合理的点云滤波算法就能准确还原地物的高程信息。本研究将以MCP为例,对激光测高仪的性能进行分析及仿真计算。

3 激光测距的性能

3.1 随机误差

3.1.1 随机误差模型

若要获得三维的地面高程信息,就需要建立地面坐标和高程信息之间的映射。垂直于测距方向的误差,即被测距点在水平坐标上的定位误差,会影响测距精度。在测距过程中,对卫星精确定位、定轨以及姿态控制是星载激光测高仪精确工作的基本要求。将一组特定的激光测高信息考虑为成千上万激光测高点云中的一个点,则每一组测高信息应该包括水平方向上的坐标(x,y)以及高程坐标z。由于地面有粗糙度和坡度,所以水平方向的误差会影响高程的精度。水平方向的坐标是通过对卫星精确定位、定轨和定姿,加上激光器光斑模式的先验知识来确定的。误差Δ用半峰全宽(FWHM)来表示,水平方向的误差来源主要有4个[18],可以表示为

Δx2=Δd2+Δa2+Δp2+Δl2,(2)

式中:Δd为衍射误差;Δa为大气湍流引起的横向误差;Δp为卫星姿态控制引起的横向误差;Δl为卫星定位误差。

取光斑的中心为理想(x,y)坐标,衍射误差为

Δd=RλD,(3)

式中:R为卫星轨道高度;D为扩束望远镜口径;λ为激光的中心波长。激光器发射的脉冲经望远镜扩束后投射到地面形成光斑,取光斑中心为理想(x,y)坐标,衍射误差就是激光投射在地面的光斑的半峰全宽。

大气湍流引起的横向误差为

Δa=hturbΔζ,(4)

式中:hturb为对流层顶的高度;Δζ为大气湍流参数,用于表示湍流在激光传播途径中引起的偏移角度。

卫星指向不准确引起的水平方向误差可以表示为

Δp=RΔθ,(5)

式中:R为卫星轨道高度;Δθ为指向误差。

沿测距方向的误差Δz主要有三部分组成,表达式为

Δz2=Δt2+Δs2+Δl2,(6)

式中:Δt为飞行时间误差引起的测距误差;Δs为地表斜面引起的测距误差;Δl为卫星定位误差。

飞行时间误差引起的测距误差可以表示为

Δt=ΔTOFc2,(7)ΔTOF2=Δpulse2+Δdetermination2,(8)

式中:ΔTOF为飞行时间误差,飞行时间误差主要有两个来源:激光脉冲的半峰全宽脉宽Δpulse和判决电路延迟误差Δdetermination

地表斜面引起的测距误差可以表示为

Δs=tanβΔx,(9)

式中:tan β为被测地表的坡度;Δx为水平方向上的坐标误差。

随机误差的模型为加性高斯噪声,概率密度分布符合均值为0,半峰全宽为Δz,高斯分布的半峰全宽(VFWHM)和方差σ的转换关系为VFWHM=2.355σ

3.1.2 随机误差数值计算

为了计算星载激光测高仪的随机误差,以ICESat-2上搭载的激光测高仪ATLAS参数作为基础[27]设置仿真参数。ATLAS的性能参数如表1所示。

表 1. ATLAS的性能参数

Table 1. Performance parameters of ATLAS

ParameterValue
Laser wavelength /nm532
Spot diameter /m10
Pulse duration /ns1
Repetition rate /kHz10
Repetition period (30 km) /ms0.1
Peak laser power /mJ250-900
Orbital height /km600
Step /m0.688

查看所有表

光斑直径即(2)和(3)式中的Δd取值为10 m。(4)式中的hturb为对流层顶高度,典型值取中纬度地区的对流层顶高度10 km。Δζ为大气湍流参数,在非恶劣天气条件下,取值10-5 rad[18]。将hturbΔζ代入(4)式可得Δa=0.1 m。(5)式中R为激光测高仪距离地面的高度,对于星载激光测高仪,R就是卫星轨道高度,取600 km。(5)式中的Δθ是指向角度误差,由先验知识可知,Δθ≪10-5 rad[18]。以最坏的情况来计算,Δp=6 m。由先验知识可知,卫星定位误差Δl≪10 cm。以最坏的情况计算,Δl=10 cm[18]。将ΔdΔaΔpΔl代入(2)式,综合4个误差来源,可得到Δx=11.66 m。

在(8)式中,Δpulse为激光脉冲的半峰全宽脉宽,取值为1 ns,Δdetermination为判决电路延迟带来的误差,一般低于100 ps[18],取最坏情况下的100 ps。将其代入(8)式,得到ΔTOF=1.005 ns。将(8)式的结果代入(7)式,把时间误差转换为测距误差,得到Δt=0.1507 cm。(9)式中tan β为被测表面的坡度,取高速公路的典型坡度0.03代入(9)式,得Δs=0.35 m。将ΔtΔsΔl代入(6)式,综合3个误差来源,得到Δz=0.394 m。

3.2 回波信号分析

3.2.1 回波光子数

太阳光的成分比较复杂,是激光探测的主要干扰光源,太阳光的一部分通过大气散射进入探测器的镜头,一部分被地表反射进入镜头,通过探测器的窄带滤波器后被单光子器件响应。由于太阳发光在时间上是均匀的,所以这部分光子均匀地散布在整个空间,不带有特定的时间特征,这部分噪声就是背景噪声。

对于单光子计数,入瞳处接收到的背景噪声光子数可以表示为

Nnoise=Lhv·AtelescopeπR2·Aground·ΔλF·ΔtW,(10)

式中:L为地面辐亮度;h为普朗克常数;v为激的中心光频率;hv为单光子能量;Atelescope为望远镜入瞳的面积;R为卫星轨道高度;Aground为地元面积;ΔλF为滤波器带宽;ΔtW为探测器的时间窗口。

激光器发射脉冲后经地面反射回到望远镜入瞳处的信号光子数可以表示为

Nsig=Qlaserhv·α·T2·AtelescopeπR2,(11)

式中:Qlaser为脉冲能量;α为地面反射率;T为大气透过率。

激光的带宽可以表示为

ΔtL·ΔvL=0.44,(12)ΔλL=λΔvLv,(13)

式中:ΔtL为激光的脉宽;ΔvL为激光的带宽,用频率来描述;ΔλL为激光的波长带宽;λ为激光的中心波长。(12)式描述的是高斯型脉冲激光的时间带宽积[28],式(13)描述波长带宽和频率带宽的转换关系。因为一般的滤光片用波长来描述通带中心和通带宽度,所以有必要把ΔvL转换为ΔλL

3.2.2 信噪比的数值计算

对于单光子探测的激光脉冲,其单脉冲的能量低于传统线性探测巨脉冲的能量,但脉冲的重复频率高,ATLAS脉冲的重复频率为10 kHz[27],两个相邻脉冲的时间间隔为0.1 ms。为了防止混淆,用不同的时间窗口来区分不同的出射脉冲,ICESat-2的轨道高度R=600 km[27],考虑地表高度为-424~8848 m(死海和珠穆朗玛峰的高度差),如果只关心地表的高程信息,则只关心发射脉冲后的3.941~4.0028 ms的时间窗口,时间窗口宽为0.0618 ms。结果表明,对于每一个回波脉冲,都能对应一个特定的脉冲发射时间,考虑到轨道可能存在误差,以及为了方便判决电路,只记录每个到达探测器的光子的时间戳,在后续的电路中将光子到达时间对应到一个特定的发射脉冲的时间,这意味着时间窗口是覆盖全部时间段的,即时间窗口为0.1 ms。

ATLAS激光器的中心波长λ=532 nm[27],中心频率v=5.6391×1011 Hz,将它们代入(12)~(13)式可得ΔλL=0.0016 nm。对于高斯型脉冲,带宽大于6σ就可以保证99%的能量通过,常见的1 nm窄带滤波器即可满足需求。

使用ModTran来计算探测器接收到的背景光强度,正午太阳直射地表时对探测器的影响最大,本研究比较关心太阳光中532 nm的分量,所以选用6种不同的大气模型。在太阳直射地表和探测器天顶角为180°的条件下,计算经大气散射和地表反射后探测器观察到的地表方向的辐亮度,以及大气的透过率,结果如表2所示。

保守估计下,选择最坏的参数,滤波器ΔλF=1 nm,时间窗口ΔtW=0.1 ms,大气模型选择透过率最低的亚极地夏季大气,T=0.621,望远镜入瞳口径取0.2 m。将这些参数代入(10)式,可得:地面反射率为0时,背景噪声的功率为2.7827×10-13 W,光子数为74.4;地面反射率为1时,背景噪声功率为4.3556×10-12 W,光子数为1165.7。

表 2. 典型大气模型的影响

Table 2. Influence of typical atmosphere model

Atmosphere modelRadiance /(W·cm-2·μm-1·sr-1)Transmittance
Albedo is 0Albedo is 1
Tropical103.935×10-36.400×10-20.625
Mid-latitude summer4.057×10-36.365×10-20.622
Mid-latitude winter3.911×10-36.320×10-20.640
Sub-arctic summer4.060×10-36.355×10-20.621
Sub-arctic winter3.915×10-36.336×10-20.640
US standard 19764.061×10-36.358×10-20.622

查看所有表

将同样的参数代入(11)式可得:当地面反射率为1时,探测器接收到的信号的单脉冲能量为9.641×10-15 J,光子数为25803.1。

不同地表反射率下信噪比的变化如图2所示。

图 2. 反射率和信噪比的关系

Fig. 2. Relationship between signal-to-noise ratio and albedo

下载图片 查看所有图片

图2可以看出,在单光子探测器的激光测高仪中,有必要使用能量大的激光器,因为激光器的能量关系到信噪比,这将直接影响测高的精度。取地表反射率为0.6,信号光子数为15481.8,背景噪声光子数为729.2,信噪比为21.2。

3.2.3 回波光子分布

在激光测距中,时间t和高程有对应关系,且与光速c有关,可以表示为

2×(R-z)=c·t,(14)

式中:z为被测地面的海拔高度。为了方便计算误差以及后续仿真,建立回波模型时将距离坐标作为横轴,这样做与将时间作为横轴是等价的。

假设被测地面的高度为100 m,测高误差为σ=0.167 m的高斯误差,回波信号的光子数为15481.8,回波脉冲波形如图3所示。

图 3. 回波脉冲波形

Fig. 3. Echo pulse waveform

下载图片 查看所有图片

背景噪声的光子总数为729.2,背景噪声分布如图4所示。

图 4. 背景噪声分布

Fig. 4. Background noise distribution

下载图片 查看所有图片

3.3 探测概率分析

3.3.1 探测概率模型

MCP的死时间约为2 ns,CDL也会引入纳秒级别的死时间,根据以上的计算可知回波脉冲宽度为2.6 ns,MCP/CDL的系统死时间大于回波脉冲宽度,整个回波脉冲最多只能触发一次计数,所以设定不合理的探测率可能会引入系统误差。

每个独立的光子触发探测器是独立事件,探测器探测到的光子的概率为

P=1-(1-Pphoton)N,(15)

式中:Pphoton为单个光子被探测到的概率;N为单脉冲的光子数;P为探测到至少一个光子的概率。

探测器被光子触发的概率分布满足几何分布,可以表示为

Pindex(x=k)=(1-Pphoton)k-1·Pphoton,(16)Pphoton=1-(1-P)1N,(17)

式中:Pindex(x=k)为探测器被第k个光子触发的概率。在已知回波信号模型时,可以得到与高程z相关的探测概率模型,可以表示为

Pr(z)=-(1-Pphoton)F(-z)·f(-z),(18)F(z)=-zf(t)dt,(19)

式中:f(z)为光子密度函数,表示单位高程内接收到的信号光子数;F(z)为f(z)的累积分布函数。

3.3.2 探测概率的数值计算

在整个回波脉冲中,若干个光子进入望远镜,探测率过低会导致探测的光子数不足,信息量不足,很难重建地面的高程模型。但探测率也不是越高越好,若探测率过高,探测器容易被回波脉冲的前沿触发,引入系统误差,甚至可能被噪声光子触发,导致探测器一直处于死时间,回波脉冲就不能以正确的时间触发探测器。

考虑探测器第一次被触发的时间以及引入的误差的分布,以几何分布的模型为基础,计算得到测距误差的概率分布,结果如图5所示。

图 5. 不同探测率下的误差概率分布

Fig. 5. Probability distribution of errors at different detection rates

下载图片 查看所有图片

不同探测率下峰值探测率和半峰全宽如表3所示。

表 3. 不同探测率下的峰值探测率和半峰全宽

Table 3. Peak detection rate and FWHM at different detection rates

Detection rate /%Peak detection rate /%FWHM /cm
100.2439.32
501.2139.30
802.0638.45
902.4735.56
993.2233.19
99.93.6425.07

查看所有表

测距误差的分布不满足高斯分布,但呈类似于高斯分布的形状,过高的探测率会导致探测器响应回波脉冲的前沿,使测量得到的飞行时间偏小,测量的高程偏高。在图6中,随着探测率增加,误差分布的峰值增大,半峰全宽变窄,波形整体向右偏移,这意味着增大探测率会导致系统误差变大。

图 6. 系统误差和探测率的关系

Fig. 6. Relationship between system error and detection rate

下载图片 查看所有图片

由于探测器响应回波脉冲的前沿,所以会导致测量距离偏小,其带来的误差是正值,而且随着探测率增加,误差增大。在一般的激光测距系统中,单光子探测器的探测率选用80%就能兼顾误差和测距效率,如图7所示。

随着探测率增加,峰值高度增高,半峰全宽降低。峰值半峰全宽积定义为峰值高度与半峰全宽的乘积,是一个无量纲的量,从图7(c)中可以看出峰值半峰全宽积约等于探测率。

4 仿真

4.1 提取高程信息

为了验证系统的性能,以及验证点云信息处理算法,对激光测高系统进行建模仿真。首先建立地物模型,地物模型如图8所示。

图 7. 不同参数与探测率的关系。(a)半峰全宽与探测率的关系;(b)峰值探测率与探测率的关系;(c)峰值半峰全宽积与探测率的关系

Fig. 7. Relationship between different parameters and detection rate. (a) Relationship between FWHM and detection rate; (b) relationship between peak detection rate and detection rate; (c) relationship between detection rate and the product of peak height and FWHM

下载图片 查看所有图片

图 8. 输入的高程模型

Fig. 8. Input topography model

下载图片 查看所有图片

单地元尺寸为10 m,整个地物模型包含64×64个地元,分为地面背景和前景,地面背景高度为0,前景有4个部分:B1是横截面尺寸为300 m×400 m的方盒,高度为1 m;B2是横截面尺寸为150 m×150 m的方盒,高度为1.75 m;B3是横截面尺寸为300 m×20 m的窄条形方盒,高度为1.4 m;C1是直径为80 m的圆柱体,高度为1 m。

使用蒙特卡罗模型进行仿真,光斑直径为10 m,则匹配的地元尺寸为10 m。根据600 km的近地轨道以及10 kHz的重复频率,可以计算得到两个光斑之间的步长为0.688 m。对每个像元都进行14次脉冲发射,因为单脉冲触发探测器的概率为80%,所以平均有11个事件被记录。

点云图是一系列带有(x,y)坐标的事件在三维坐标系中的散点图,根据前面的计算可知Δz为0.39 m。图9分别是Δz取0.1,0.2,0.3,0.39 m的点云图。

显然,误差越大,从点云中提取高程信息就会越困难。为了滤除背景噪声,以及抵消横向误差Δx,综合邻近的3×3范围内像素的所有点云数据,对高程坐标z进行中值滤波。滤波的结果如图10所示。

图 9. 事件点云图。(a) Δz=0.1 m; (b)Δz=0.2 m; (c) Δz=0.3 m; (d) Δz=0.39 m

Fig. 9. Point clouds of event. (a) Δz=0.1 m; (b) Δz=0.2 m; (c) Δz=0.3 m; (d) Δz=0.39 m

下载图片 查看所有图片

图 10. 输出的高程模型

Fig. 10. Output topography model

下载图片 查看所有图片

滤波之后基本能还原地物模型的形状,但因为误差和噪声的影响,平台处不平整。高程信息的直方图如图11所示。

图 11. 输出高程的直方图

Fig. 11. Histogram of output topography

下载图片 查看所有图片

图11中可以看出,受误差和噪声影响,4个不同高程的尖峰都有一定宽度,B1和B3的尖峰甚至出现混叠。在更恶劣的条件下,可能无法区分两个尖峰。

评价测高精度一般用均方根误差,其定义为

eRMS=1MNm,n[hout(m,n)-hin(m,n)]2,(20)

式中:eRMS为均方根误差;hin为高程的真实值;hout为高程的测量值;MN为被测表面的像元数;(m,n)为特定像元的坐标。计算得到高程均方根误差为6.1 cm。

4.2 优化的高程信息

可以注意到,激光测高与图像分割有相似的形式,高程信息的突变处往往伴随光影的明显变化,比如森林的边界、楼房的阴影、海洋沿岸。对于大部分卫星,激光测高仪是卫星上众多载荷中的一个,一般具有遥感观测能力的卫星,都配有遥感相机。典型的遥感相机有R、G、B三个通道,或加一个近红外NI通道共4个通道。可以使用多光谱相机所成的像来提取图像的边界,并对带有误差的高程信息进行优化。优化的目标函数表示为[29]

E(h)=hTLh+λ'(h-h˙)T(h-h˙),(21)

式中: h˙为利用上述方法直接测量得到的粗测高程信息;h为优化的高程信息;L为多光谱图像的拉普拉斯矩阵,用来描述多光谱图像邻近像素之间的差别,用两个像素在RGB坐标系中的马氏距离来描述;λ'为权重系数,用来调整两个条件的权重。h优化时遵循两个条件:1)与直接测量结果接近;2)在多光谱图像的边界处变化。目标函数E(h)等号右边有两项,第一项约束高程信息的变化,使高程信息变化与多光谱图像的边界相匹配。第二项用直接测量结果 h˙来约束高程信息h。拉普拉斯矩阵L的定义为[30]

Li,j=k|(i,j)wkδi,j-1|wk|1+(Ii-μk)TΣk+ε|wk|U3-1(Ij-μk)(22)

式中:ijk为像素的序号;δi,j为克罗内克符号函数;wk为以像素k为中心的窗口中所有像素的集合,在这里选择3×3的窗口;IiIj分别为像素ij的RGB三通道值;μkwk中所有像素的RGB平均值;Σkwk中像素RGB值的协方差矩阵;ε为正则化系数;U3为3×3的单位矩阵。在计算拉普拉斯矩阵时,用一维坐标表示图像中像素的序号,在64×64的图像中,序号取值范围为1~4096。

使目标函数极小的h可以用稀疏线性方程组表示为

(L+λ'U)h=λ'h˙,(23)

式中:U为与L同样大小的单位矩阵。h的解可以表示为

h=(L+λ'U)-1λ'h˙(24)

把地物模型转化为RGB三通道的遥感图像,并给地面、B1、B2、B3、C1涂上不同的颜色,如图12所示。

图 12. 地物模型的遥感图像

Fig. 12. Remote sensing image of ground object model

下载图片 查看所有图片

ε=10-7,将遥感图像代入(22)式得到拉普拉斯矩阵。将λ'=10-4代入(24)式优化高程模型,得到优化后的高程模型如图13所示。

图 13. 优化后的高程模型

Fig. 13. Topography model after optimizing

下载图片 查看所有图片

优化后的高程直方图如图14所示。

图 14. 优化后的高程直方图

Fig. 14. Topography histogram after optimizing

下载图片 查看所有图片

这样的修正可以消除大部分噪声,并可以还原地物模型的真实形状。优化后,高程的均方根误差为2.6 cm,说明用遥感图像优化高程信息的方法是有效的。

5 结论

目前,基于单光子探测器的激光测高仪尚未应用在卫星平台上,理论研究也不完善,本研究对单光子激光测高技术进行梳理,从激光的数学模型出发,到大气模型、地物模型,对激光测高的各个环节都进行了讨论。代入典型值进行数值计算,得到激光回波的半峰全宽为0.39 m,标准差为0.17 m。为研究激光测高点云的处理方法,搭建仿真地物模型,用蒙特卡罗方法进行仿真计算,采用中值滤波处理点云,计算得到高程的均方根误差为6.1 cm,低于激光的标准差,这是因为高程信息是由多个光子事件联合计算得到的,相当于多次采样取平均值。对于单光子激光测高仪,如果能提高激光脉冲的重复频率,就可以进一步提高测高精度。最后提出一种利用遥感的多光谱图像来优化高程信息的算法,先验知识表明高程信息变化的地方往往是遥感图像RGB通道发生变化的地方,所以用遥感图像的边界来约束高程信息进行优化。仿真结果表明优化后的高程的均方根误差为2.6 cm,优于中值滤波的结果。

参考文献

[1] 李然, 王成, 苏国中, 等. 星载激光雷达的发展与应用[J]. 科技导报, 2007, 25(14): 58-63.

    Li R, Wang C, Su G Z, et al. Development and applications of spaceborne LiDAR[J]. Science & Technology Review, 2007, 25(14): 58-63.

[2] Nozette S, Rustan P, Pleasance L P, et al. The clementine mission to the moon: scientific overview[J]. Science, 1994, 266(5192): 1835-1839.

[3] Nozette S, Lichtenberg C L, Spudis P, et al. The clementine bistatic radar experiment[J]. Science, 1996, 274(5292): 1495-1498.

[4] 唐新明, 李国元. 激光测高卫星的发展与展望[J]. 国际太空, 2017( 11): 13- 18.

    Tang XM, Li GY. Development and prospect of laser altimetry satellite[J]. Space International, 2017( 11): 13- 18.

[5] 杨馥, 贺岩, 周田华, 等. 基于伪随机码调制和单光子计数的星载测高计仿真[J]. 光学学报, 2009, 29(1): 21-26.

    Yang F, He Y, Zhou T H, et al. Simulation of space-borne altimeter based on pseudorandom modulation and single-photon counting[J]. Acta Optica Sinica, 2009, 29(1): 21-26.

[6] Schutz BE. Laser altimetry and lidar from ICESat/GLAS[C]. 2001 International Geoscience and Remote Sensing Symposium, Sydney, 2001: 1016- 1019.

[7] Smith D E, Zuber M T, Jackson G B, et al. The lunar orbiter laser altimeter investigation on the lunar reconnaissance orbiter mission[J]. Space Science Reviews, 2010, 150(1/2/3/4): 209-241.

[8] 王建宇, 舒嵘, 陈卫标, 等. 嫦娥一号卫星载激光高度计[J]. 中国科学: 物理学力学天文学, 2010, 40(8): 1063-1070.

    Wang J Y, Shu R, Chen W B, et al. Lather altimeters of Chang'e-1[J]. Scientia Sinica Physica, Mechanica & Astronomica, 2010, 40(8): 1063-1070.

[9] . 嫦娥一号卫星的初步科学成果与嫦娥二号卫星的使命[J]. 航天器工程, 2010, 19(5): 1-6.

    Ouyang Z Y. Science results of Chang'e-1 lunar orbiter and mission goals of Chang'e-2[J]. Spacecraft Engineering, 2010, 19(5): 1-6.

[10] 宋博, 李旭, 郑伟, 等. 02)星载高精度激光测距技术的实现[J]. 光电子技术, 2017, 37(1): 61-65.

    Song B, Li X, Zheng W, et al. The implementation of high precision space-borne laser ranging technology in ZY-3(02) satellite[J]. Optoelectronic Technology, 2017, 37(1): 61-65.

[11] Santovito M R, Tommasi L, Sgarzi G, et al. A laser altimeter for BepiColombo mission: instrument design and performance model[J]. Planetary and Space Science, 2006, 54(7): 645-660.

[12] Araki H, Tazawa S, Noda H, et al. Lunar global shape and polar topography derived from Kaguya-LALT laser altimetry[J]. Science, 2009, 323(5916): 897-900.

[13] Kamalakar J A. Bhaskar K V S, Prasad A S L, et al. Lunar ranging instrument for Chandrayaan-1[J]. Journal of Earth System Science, 2005, 114(6): 725-731.

[14] 季云飞, 耿林, 冯国旭, 等. 激光测高技术的发展趋势[J]. 激光与红外, 2011, 41(8): 830-833.

    Ji Y F, Geng L, Feng G X, et al. Progress and prospect of laser altimeter technology[J]. Laser & Infrared, 2011, 41(8): 830-833.

[15] Abdalati W, Zwally H J, Bindschadler R, et al. The ICESat-2 laser altimetry mission[J]. Proceedings of the IEEE, 2010, 98(5): 735-751.

[16] Yu AW, Harding DJ, KrainakM, et al. Development of an airborne lidar surface topography simulator[C]. 2011 Laser Applications to Photonic Applications, CLEO: Applications and Technology, Baltimore, 2011: 1- 2.

[17] 蔡厚智, 刘进元, 付文勇, 等. 基于微通道板选通的飞行时间测量技术[J]. 光学学报, 2018, 38(2): 0204002.

    Cai H Z, Liu J Y, Fu W Y, et al. Measurement technology of time of flight based on gated microchannel plates[J]. Acta Optica Sinica, 2018, 38(2): 0204002.

[18] Priedhorsky WC, Smith RC, ChengH. Laser ranging and mapping with a photon-counting detector[C].Proceedings of SPIE, 1995: 441- 52.

[19] Baron MH, Priedhorsky WC. Crossed-delay line detector for ground- and space-based applications[C]. Proceedings of SPIE, 1993: 188- 198.

[20] Aull BF, Marino RM. Three-dimensional imaging with arrays of Geiger-mode avalanche photodiodes[C]. Proceedings of SPIE, 2005, 6014: 60140D.

[21] 张国青, 张英堂, 翟学军, 等. 多像素光子计数器的信噪比特性[J]. 光学学报, 2013, 33(3): 0304001.

    Zhang G Q, Zhang Y T, Zhai X J, et al. Signal-to-noise ratio properties of multi-pixel photon counter[J]. Acta Optica Sinica, 2013, 33(3): 0304001.

[22] Shao Y, Silverman R W, Farrell R, et al. Design studies of a high resolution PET detector using APD arrays[J]. IEEE Transactions on Nuclear Science, 2000, 47(3): 1051-1057.

[23] Kataoka J, Koizumi M, Tanaka S, et al. Development of large-area, reverse-type APD-arrays for high-resolution medical imaging[J]. Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, 2009, 604(1/2): 323-326.

[24] 白廷柱, 金伟其. 光电成像原理[M]. 北京: 北京理工大学出版社, 2013: 127- 128.

    Bai TZ, Jin WQ. Photoelectric imaging principle[M]. Beijing: Beijing Institute of Technology Press, 2013: 127- 128.

[25] Vallerga J V. Siegmund O H W. 2 K×2 K resolution element photon counting MCP sensor with >200 kHz event rate capability[J]. Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, 2000, 442(1/2/3): 159-163.

[26] DegnanJ. Impact of receiver deadtime on photon-counting SLR and altimetry during daylight operations[C]. 16th International Workshop on Laser Ranging, 2008: 339- 346.

[27] 许艺腾. 单光子激光测高数据处理技术研究[D]. 西安: 西安科技大学, 2017: 5- 6.

    Xu YT. Research on data processing technology of single photon laser altimetry[D]. Xi'an: Xi'an University of Science andTechnology, 2017: 5- 6.

[28] 陈聪, 陈海燕. 洛伦兹与超高斯矩形超短脉冲的时间带宽积[J]. 长江大学学报(自然科学版), 2013, 10(1): 10-11.

    Chen C, Chen H Y. Time-band width products for Lorentz and super-Gaussian rectangular line-shape ultrashort pulse lasers[J]. Journal of Yangtze University(Natural Science Edition), 2013, 10(1): 10-11.

[29] He KM, SunJ, Tang XO. Single image haze removal using dark channel prior[C]. Miami: 2009 IEEE Conference on Computer Vision and Pattern Recognition, 2009: 1956- 1963.

[30] Levin A, Lischinski D, Weiss Y. A closed-form solution to natural image matting[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2008, 30(2): 228-242.

米小什, 赵惠, 樊学武, 盛立志. 基于微通道板的单光子激光测高技术研究[J]. 光学学报, 2018, 38(12): 1228001. Xiaoshi Mi, Hui Zhao, Xuewu Fan, Lizhi Sheng. Photon-Counting Laser Altimetry Based on Microchannel Plate[J]. Acta Optica Sinica, 2018, 38(12): 1228001.

本文已被 3 篇论文引用
被引统计数据来源于中国光学期刊网
引用该论文: TXT   |   EndNote

相关论文

加载中...

关于本站 Cookie 的使用提示

中国光学期刊网使用基于 cookie 的技术来更好地为您提供各项服务,点击此处了解我们的隐私策略。 如您需继续使用本网站,请您授权我们使用本地 cookie 来保存部分信息。
全站搜索
您最值得信赖的光电行业旗舰网络服务平台!