李海超 满益云
(钱学森空间技术实验室,北京 100094)
线阵CCD倾斜布置敏捷卫星辐射定标的姿态计算方法
李海超 满益云
(钱学森空间技术实验室,北京 100094)
针对搭载倾斜布置的多片线阵CCD相机的敏捷卫星上所有探测器对同一区域成像的相对辐射定标,提出了基于地面成像轨迹重合的偏航角速度计算模型与方法,充分利用敏捷卫星对偏航角速度的控制潜力。建立了地面物点到相机焦平面的空间坐标变换关系;推导了斜视成像时地心角的计算公式;给出了对同一区域成像的敏捷卫星偏航角速度计算方法。以5片线阵CCD共30 000个探测器、相邻两片倾斜0.3°的敏捷卫星作为仿真对象,给出了所有探测器在约30 s的定标成像过程中的地面轨迹间隔宽度,结果表明:为了实现所有探测器对同一区域成像,需要间隔一定周期调整偏航角速度,当调整周期为0.36 s时,对应的地面轨迹宽度在赤道附近优于1个像元、在纬度50°附近优于6个像元。文章对敏捷卫星基于同一区域成像的相对辐射定标可提供一定的参考价值。
敏捷卫星;线阵CCD相机;相对辐射定标;地面轨迹;偏航角速度
遥感图像的质量由获取的图像数据质量与图像处理的质量两方面决定,而图像数据的获取质量由遥感器的综合性能决定。由于光学遥感器中各探测器之间的响应特性不一致,对于同样输入的辐射能量其输出响应存在一定差异,在获取的图像中表现为条带现象[1]。因此,为获取高质量的遥感图像,遥感器有必要进行相对辐射定标,给出每个探测器的相对辐射定标系数,以消除或降低各个探测器之间的响应不一致性带来的影响。目前,在轨卫星的相对辐射定标主要通过地面均匀场地、在轨遥感图像统计等手段实现。然而,地面均匀场地法需要寻找或铺设足够大的均匀场地,如自然湖泊、草原、冰川、人工定标场等,其特点是寻找难度大、建设成本高;目前应用较广的在轨遥感图像统计方法[2-3],需要多轨图像数据,因此,该类方法的前提是确保各轨图像数据的稳定性。
随着遥感卫星技术的发展,具备敏捷机动能力的卫星已经成为遥感卫星的一个重要发展方向[4]。国外近几年广泛应用的Side-Slither定标方法[5],正是利用了卫星敏捷的姿态控制将卫星偏航角旋转90°,使得多片线阵CCD对同一区域成像成为可能,实现了所有的探测器对同一区域成像的相对辐射定标。近年来,以法国昴宿星(Pleiades)为代表的敏捷卫星,为减小光学畸变的影响,该卫星的全色谱段探测器采取将5片线阵CCD阵列倾斜布置在焦平面上,用“弧线”形状代替常用的“直线”形状布置方式[6]。因此,对于这类倾斜布置的线阵CCD实现Side-Slither定标,需要借助卫星较强的敏捷能力。
本文为实现对同一区域成像的相对辐射定标,提供一种多片倾斜布置的敏捷卫星偏航角速度计算方法:建立了卫星斜视成像的空间坐标变换关系;在此基础上,给出了基于空间坐标变换关系的对同一区域成像的卫星偏航角速度计算方法;定标成像过程中根据敏捷卫星姿态控制能力不断调整卫星偏航角速度。
2.1 坐标系定义
从地面景物到光学遥感器的坐标变化过程中,需要定义如下7个坐标系[7],本文中坐标系一律采用右手坐标系:
(1)卫星轨道坐标系B(B0-B1B2B3):原点B0在轨道上,B1轴在轨道面上指向卫星飞行方向,B3轴过地心并指向天顶,B2轴与轨道面垂直,卫星以角速度Ω在轨道上运动。
(2)卫星本体坐标系S(S0-S1S2S3):原点S0与B0重合,卫星无姿态运动时S与B系重合,滚动角φ、俯仰角θ、偏航角ψ即指S系在B系内的三轴姿态(顺序φ→θ→ψ)。
(3)相机坐标系C(C0-C1C2C3):无安装误差时相机坐标系与本体坐标系重合。
(4)焦平面坐标系P(P0-P1P2P3):多片线阵CCD的中点P0为焦平面坐标系的原点;P1,P2,P3三个坐标轴过P0点,且分别平行于本体坐标系的滚动轴、俯仰轴和偏航轴。
(5)地面坐标系G(G0-G1G2G3):原点G0为P0在地面的成像点,G1轴平行于卫星飞行方向,G3轴过地心指向天顶,G2轴垂直于G1G3构成的平面。
(6)参考坐标系I(O-I1I2I3):原点O在地心,I2轴指向北极,I3轴指向卫星轨道升交点,I1轴垂直I2和I3构成的平面。
(7)固定坐标系E(O-E1E2E3):该坐标系固联于地球,原点O与I系原点重合(在地心),E2轴指向北极,与I2轴重合,该坐标系在I坐标系内绕E2(即I2)、逆时针方向以地球角速度ω自转,该坐标系在起始时刻与I系重合。
2.2 卫星斜视成像的坐标变换
为适应卫星斜视成像的需求,本文给出了非星下点斜视成像的坐标变换模型。图1为斜视成像示意图,其对应的从地面坐标系G到焦平面坐标系P的坐标转换顺序见表1。
根据表1和图1,可得到地面坐标系G中的点(G1,G2,G3)到焦平面坐标系P中的点(P1,P2,P3)的变换关系,为表示两个坐标系间的平移、旋转和缩放,采用齐次坐标变换矩阵表示见式(1)。
图1 从地面坐标系G到焦平面坐标系PFig.1 Transformation from ground coordinate G to focal plane coordinate P
其中,R为地球半径;h为地面点G0的地形高度;H为航天器的轨道高度;f为相机焦距;γ=γ0+Ωt,γ0为摄影时刻升交点与卫星之间所对应的地心角,Ω为卫星轨道运行相对地心的角速度;i0为卫星轨道倾角;ω为地球自转角速度;ψ=ψ0+dψt,θ=θ0+dθt,φ=φ0+dφt;ψ0、θ0、φ0分别表示卫星的初始偏航角、俯仰角、滚动角,dψ,dθ,dφ分别表示卫星的偏航、俯仰和滚动角速度。
图1(b)所示,∠K0OK1和∠K1OG0是卫星斜视成像对应的地心角,计算方法如下。
在ΔB0OK1中:
图2是卫星从常规成像模式进入在轨定标成像模式的示意图,其中,轨道S位置为常规成像模式,此时首先计算卫星对应的初始偏流角:
其中,VP1和VP2分别为相机像面的前向像移速度和横向像移速度。设定某一时刻地面上一物点经过相机成像,对应的像点的坐标为(P11,P21),经过足够短的时间Δt之后,该地面上的物点在像平面上的像点运动到(P12,P22)位置处,利用差分法计算像移速度[8]:
然后,经过卫星偏航角旋转90°―β0的姿态调整,在S0位置卫星进入相对辐射定标成像模式。考虑偏流角β0的主要目的是消除地球自转对成像的影响[9-10]。
图2 从常规成像模式进入在轨定标成像模式Fig.2 From conventional imaging to on-orbit calibration imaging
3.1 倾斜布置的线阵CCD相对辐射定标原理
根据Side-Slither相对辐射定标原理[5],线阵CCD的所有探测器应该依次对地面同一区域成像(即要求所有探测器的地面成像轨迹理论上应该重合),获取同一区域的成像数据。对于倾斜布置的多片线阵CCD的相对辐射定标原理描述如下:
如图3(a)所示,设定位置S0为卫星进行相对辐射定标成像的起始位置(对于时刻t0=0);如图3(b)所示,该位置所对应的焦平面坐标系作为原焦平面坐标系P0-P1P2P3,对应的地面坐标系作为原地面坐标系G0-G1G2G3。则此时线阵CCD对应的地面轨迹称为初始地面轨迹。
随着卫星在轨道的运动,原焦平面坐标系以线阵CCD的像元尺寸de作为步长沿着卫星飞行方向进行平移,同时原焦平面坐标系绕着坐标原点P0以偏航角速度进行旋转。经过时间t1=n1×te后,原焦平面坐标系的原点P0运动到新焦平面坐标系P0'-P1'P2'P3'的原点P0',卫星到达位置S1,如图3(c)所示,其中te为曝光时间,n1为原焦平面坐标系移动的像元个数。此时刻线阵CCD对应的地面轨迹理论上应该要求与初始地面轨迹重合。
由于本文中多片线阵CCD的倾斜布置,为确保所有探测器先后依次对地面同一区域成像,因此,成像过程中必须对卫星偏航角速度进行不断调整。
图3 其地面轨迹作为条件计算偏航角速度的示意图Fig.3 Diagram of satellite YAV computation for all CCD detectors passed through the same trajectories
3.2 偏航角速度计算步骤
为实现对同一区域成像的偏航角速度计算方法分为以下6个步骤:
(1)在t0时刻根据式(1)计算倾斜布置的线阵CCD对应的初始地面轨迹,如图3(b)所示。
(2)经过时间t1到达位置S1之后,如图3(c)所示,P0'在原焦平面坐标系中的坐标P0'_old(p10'_old,p20'_old,0)计算为:
(3)设定从t0时刻(位置S0)开始,卫星以偏航角速度运动,则从t0到t1时间段内,卫星偏航角旋转的大小为Δβ=t1×,即P0'-P1'P2'P3'坐标系相对P0-P1P2P3坐标系旋转的角度。
(4)计算t1时刻倾斜布置的线阵CCD对应的地面轨迹:①设线阵CCD的任意点m在新焦平面坐标系中的坐标表示为(p1m_new,p2m_new,0),该坐标由线阵CCD在焦平面上的数学模型计算得到;②根据P0-P1P2P3和P0'-P1'P2'P3'之间的旋转与平移关系,点m在原焦平面坐标系中的坐标(p1m_old,p2m_old,0)计算公式见式(13);③根据式(1)计算(p1m_old,p2m_old,0)在原地面坐标系G0-G1G2G3中对应的坐标(g1M,g2M,0);④重复步骤①~③计算其它点的坐标,得到线阵CCD对应的地面轨迹。
(6)计算(5)中与初始地面轨迹距离最近的地面轨迹,其所对应的角速度即所求的偏航角速度。
本文仿真数据:5片线阵CCD倾斜布置,相邻两片夹角0.3°,每片CCD有L=6000个探测器,一共30 000个探测器,焦距2.6 m,积分时间362.59μs;卫星轨道高度H=645.338 km,地球半径R=6 378.137 km,卫星轨道周期98min,轨道倾角98°。
焦平面坐标系如图4所示,其中实际在轨运行时,5片倾斜布置的线阵CCD的数学模型可以通过在轨定标得到,本仿真实验中用圆拟合来建立线阵CCD的数学模型,将5片CCD的首末端点和中点共11个探测器,即图中的点A~J,作为分析对象,根据它们的地面轨迹来计算偏航角速度。
图4 焦平面坐标系Fig.4 Focal plane coordinate
4.1 赤道附近、偏航角速度调整周期
在t0=0时刻成像对应的纬度是0°、滚动角0°、控制偏航角速度的调整周期为0.36 s条件下,得到在大约32 s时间内成像的地面轨迹如图5所示,5片CCD中11个探测器的地面轨迹用不同的曲线表示。其中,图5(a)为大约32 s时间内的整个轨迹,图中红色虚线框区域是11个探测器都经过的区域,统计得到整个红色虚线框区域内的地面轨迹间距均值为2.157 m(其中,地面轨迹间距定义为11个探测器地面轨迹的最大距离,如图5(b)所示);为了清楚表示多个探测器之间的地面轨迹情况,如图5(b)给出了图5(a)中的区域1的放大图;图5(c)是调整周期0.36 s下的卫星偏航角速度变化曲线,可以看出,为获得图5(a)的地面轨迹,偏航角速度变化趋势先增大后减小,从―0.067 61(°)/s增大到―0.066 17(°)/s,再降低到―0.067 47(°)/s。
4.2 纬度40°附近、偏航角速度调整周期
在t0=0时刻成像对应的纬度是北纬40°、滚动角0°、控制偏航角速度的调整周期为3.6 s条件下,得到线阵CCD中的11个探测器对应的地面轨迹如图6所示。其中,图6(a)中的红色虚线框区域是11个探测器都经过的区域,统计得到整个红色虚线框区域内的轨迹距离平均值为31.05 m;图6(b)给出了图6(a)中区域2的放大图,在该区域内探测器点E与点J对应的地面轨迹间距最大,约为21 m;图6(c)是偏航角速度变化曲线。
图5 五片阵列中11个探测器的地面轨迹(赤道附近,偏航角速度调整周期为0.36 s)Fig.5 Trajectories of 11 detectors in 5 arrays(latitude 0°,YAV adjusting period is 0.36s)
图6 五片阵列中11个探测器的地面轨迹(北纬40°,偏航角速度调整周期为3.6 s)Fig.6 Trajectories of 11 detectors in 5 arrays(latitude 40°N,YAV adjusting period is 3.6s)
4.3 不同纬度不同滚动角下的地面轨迹间距分析
表2给出偏航角速度调整周期0.36 s条件下,给出了在不同纬度、不同滚动角成像时所对应的地面轨迹间距。可以看出,赤道附近对应的地面轨迹间距非常窄,例如纬度0°对应的地面轨迹间距2.157 m,5°对应的地面轨迹间距为1.378 m,因此低纬度区域的地面轨迹间距优于1个像元(对于2.5 m分辨率来说)。随着纬度的增大,地面轨迹间距越来越宽,达到50°时,间距达到了15 m(对于2.5 m分辨率来说为6个像元)。低纬度赤道附近的轨迹间距随着滚动角的增大变化明显,例如,滚动角从0°到40°其对应的轨迹间距从2.157 m变化到22.287 m。
表2 不同纬度不同滚动角下的地面轨迹间距Table 2 Trajectory interval at different latitudes and roll angles m
4.4 地面轨迹间距误差影响因素及不确定性分析
由于轨迹计算的不确定性,计算轨迹与实际轨迹引起了误差,引起偏航角速度计算误差的主要影响因素有姿态角(设定最大测量误差0.01°)、姿态角速度(设定最大测量误差0.001°)、轨道高度(设定最大测量误差100 m)、地心角(设定最大测量误差0.1°)等。
本部分对卫星在纬度0°、滚动角0°执行定标成像时,各个误差因素引起的地面轨迹间距的不确定性进行定量分析。由于各个误差因素的最大测量误差不一致,为便于在图中表示,用横坐标表示最大误差的倍数(0.1、0.2、……1.0),如图7所示。
图7 不同参数误差带来的地面轨迹间距的误差Fig.7 Trajectory interval error caused by various error factors
从图7中可以看出,姿态角速率测量误差对地面轨迹间距的影响程度最大,随着姿态角速率测量精度的下降,引起的地面轨迹间距急剧增大,轨迹间距从2.157 m(姿态角速率无误差)下降到41.837 m(姿态角速率误差0.001(°)/s),因此,为确保所有探测器对同一区域成像,姿态角速率测量精度应该严格要求;轨道高度测量精度对地面轨迹间距也有一定的影响,轨迹间距从2.157 m(轨道高度无测量误差)下降到8.347 m(轨道高度误差100 m);姿态角和地心角测量精度对轨迹间距几乎不产生影响。因此,为满足相对辐射定标需求的同一区域成像这一条件,应该重点减小姿态角速率的测量误差。
在定标成像模式下执行倾斜布置线阵CCD的Side-Slither相对辐射定标时,需要所有探测器对同一个区域进行成像,较目前传统的相对辐射定标方法,建立了相同的输入数据源。本文提出了基于严密数学模型的倾斜布置多片线阵CCD相机的相对辐射定标方法,建立了斜视成像条件下从地理坐标系到焦平面坐标系的转换关系,据此计算卫星偏航角速度,适应于斜视成像的相对辐射定标。
从仿真实验可以看出,理论上通过控制卫星姿态可以得到较小的地面轨迹间距,使得所有探测器近似地依次对地面同一区域成像,但是,对卫星姿态控制要求较高,这对于在轨卫星实现难度大。因此,在轨进行Side-Slither相对辐射定标时,为降低对卫星姿态控制的要求,可以选择具有一定宽度的均匀地面场景,从而可以适当放宽对地面轨迹间距的限制约束,根本上降低对卫星姿态控制的要求。
(References)
[1]Corsini G,Diani M,Walzel T.Striping removal in MOS-B data[J].IEEE Transactions on Geoscience and Remote Sensing,2000,38(3):1439-1446
[2]Horn B K P,Woodham R J.Destriping Landsat MSS images by histogram modification[J].Computer Graphics and Image Processing,1979,10:69-83
[3]Gadallah F L,Csillag F,Smith E J M.Destriping multi sensor imagery with moment matching[J].International Journal of Remote Sensing,2000,21:2505-2511
[4]张新伟,戴君,刘付强.敏捷遥感卫星工作模式研究[J].航天器工程,2011,20(4):32-38 Zhang Xinwei,Dai Jun,Liu Fuqiang.Satellite with agile attitude control[J].Spacecraft Engineering,2011,20(4):32-38(in Chinese)
[5]Cody A,Denis N,Andreas B,et al.Radiometric correction of Rapid Eye Imagery using the on-orbit sideslither method[C]//SPIE 8180,Image and Signal Processing for Remote Sensing XVII.Bellingham:SPIE,2011
[6]Baillarin S,Lebegue L,Kubik P.PLEIADES-HR system qualification:a focus on ground processing and image products performances,a few months before launch[C].Proceedings IGARSS.New York:IEEE,2009
[7]王家骐,于平,颜昌翔,等.航天光学遥感器像移速度矢计算数学模型[J].光学学报,2004,24(12):1585-1588 Wang Jiaqi,Yu Ping,Yan Changxiang,et al.Space optical remote sensor image motion velocity vector computational modeling[J].Acta Optica Sinica,2004,24(12):1585-1588(in Chinese)
[8]王运,颜昌翔.基于差分法的空间相机像移速度矢量计算[J].光学精密工程,2011,19(5):1054-1060 Wang Yun,Yan Changxiang.Computation of image motion velocity vector for space camera based on difference method[J].Optics and Precision Engineering,2011,19(5):1054-1060(in Chinese)
[9]景泉.敏捷卫星偏流角计算模型研究[J].航天器工程,2012,21(4):16-20 Jing Quan.Research on computation model of bias angle for agile satellite[J].Spacecraft Engineering,2012,21(4):16-20(in Chinese)
[10]朱兴鸿,邸国栋,陆春玲.一种遥感卫星偏流角修正的仿真分析方法[J].航天器工程,2013,22(1):39-43 Zhu Xinghong,Di Guodong,Lu Chunling.Simulation analysis on drift angle adjustment for remote sensing satellite[J].Spacecraft Engineering,2013,22(1):39-43(in Chinese)
(编辑:李多)
Attitude Computation of Radiometric Calibration for Agile Satellite with Tilted Linear CCD
LI Haichao MAN Yiyun
(Qian Xuesen Laboratory of Space Technology,Beijing 100094,China)
In the light of the agile satellite relative radiometric calibration method with all detectors of tilted linear CCD viewing the same scene,this paper provides the modeling and computation of satellite yaw angular velocity(YAV)based on the same trajectories of all detectors,which makes full use of the ability of yaw angle control.The space coordinate transform process from an object point on the earth to the image point in the focal plane is found.The formulas of geocentric angles are derived to be used for squint imaging mode.Furthermore,the calculation method of YAV for satellite CCD camera imaging the same scene is given.Taking an agile satellite equipped with 30,000 detectors in five CCD arrays and 0.3°tilted angle each other as an example,Simulation experiments are carried out,and the ground trajectory interval width is given during calibration process more than 30s imaging.Computation results indicate that YAV should be adjusted with certain interval period for all the detectors imaging the same scene.With YAV adjusting period of 0.36 s,the ground trajectory interval width is less than one pixel near the equator,and about six pixels at 50°latitude.This paper has important reference value for relative radiometriccalibration of the agile satellite based on all detectors imaging the same scene.
agile satellite;linear CCD camera;relative radiometric calibration;ground trajectory;yaw angular velocity(YAV)
TP722;V19
:ADOI:10.3969/j.issn.1673-8748.2015.05.007
2015-07-23;
:2015-08-30
李海超,男,博士,高级工程师,从事遥感卫星定标、遥感图像处理等研究工作。Email:lihaichao2000@163. com。