。
[0077] S4014,得到相机光学节点的运行位置为：
[0079] 其中，[《，K，马f为卫星位置归一化向量，[乂乂，剖r为相机光学节点在卫星本体 坐标系的安装位置向量在地球固连坐标系的形式。
[0080] S4015,设相机成像的积分时间值为Ato, t+A to时刻的相机光学节点的运行位置 马；^"f >*时刻相机光学节点的速度矢量为：
[0082]该速度矢量在地面的投影为：
[0084]其中Hobt为轨道高度。
[00化]设CCD尺寸为a相，a单位为米，相机主距为f，则地面分辨率GR = a/巧化bt。
[0086]实际积分时间为：
[0088]重新计算t时刻相机光学节点的速度矢量为：
[0090] 优选的：
[0091] 每片TDI CCD的安装误差用一个WCCD列坐标为变量的二阶多项式表示，根据相机 光学系统的崎变模型将不同视场的崎变量换算为像面位置误差，与安装误差多项式的常数 项相加，计算每个CCD在垂直TDI积分方向的崎变量OCSi和每个CCD在沿TDI积分方向的崎变 量OASi，两个CCD列之间的亚像元点的崎变量为：
[0092] 〇wk) = 〇Asr0.5<c-i<0.5 [OOW] 〇cs(c) = 〇cs(i+i)0.5<c-i。
[0094] 其中，i表示CCD序列号，C表示列坐标，i < C < i+1。
[00巧]优选的：
[0096] 首先将某一时刻的TDI CCD成像面细分为m*m个更小的子CCD,运时整个TDI CCD像 面可W看作一幅分辨率为模拟影像分辨率m倍的框幅式影像。设TDI CCD积分级数为M，TDI CCD单行探元数为N，则TDI CCD像面构成的框幅式影像（W下称为TDI影像)行数为M*m，列数 为
[0097] 设TDI影像四角像面坐标为^1，71]。= 1，2,3,4)，相机主距为'。则101影像四角 与光学系统主点构成向量11〇 = ^;[，7；[，-門（1 = 1，2,3,4);
[0098] 根据相机安装角、卫星平台S轴姿态角、相机主点位置P[X，Y，Z]、速度向量V[vx， vy，VZ]计算的相机安装矩阵Mset、姿态矩阵Matt、轨道-物方转换矩阵Mobt将像空间向量uO 转换为物方向量U;
[0099] 设U对应地面点高程为H，根据P、H、u可W计算地面点平面坐标化，B)。将化，B，H)通 过Wor 1 d V i ew-3对应视角的准核线影像RPC计算像面坐标，并获取该点的高程化rue;
[0100] 将H替换为化rue,重复前面的计算过程，直到Ih-H化ue I小于一定的限差化，优选 的，所述限差为0.1m;
[0101] 运样就得到TDI影像四角坐标在Worldview-3对应视角的准核线影像上的四个坐 标，也就确定了所需数据的概略范围。将该范围在TDI积分方向对应的Worldview-3准核线 影像方向外扩一定距离确定参与计算的数据范围。
[0102] 优选的：
[0103] 将确定数据范围内的Worldview-3准核线影像像点对应地面点通过TDI影像的瞬 时外方位元素投射到TDI影像上，加上相应的TDI CCD崎变改正量后，按子CCD大小取整后查 找坐标对应子CCD的行、列号，每个CCD都有一个存储单元记录投射地面元在Worldview-S准 核线影像的坐标和地面元与TDI影像投影中屯、的距离。所有地面点都完成投射后，有投射记 录的子CCD将挑选出距离最近的地面元作为可见地面元，并将对应的Worldview-3准核线影 像DN值转换为像面福照度记录下来，没有投射记录的子CCD福照度通过周围有记录的子CCD 福照度距离倒数加权法插值获得。
[0104] 优选的：
[0105] S801，计算各细分时间段中屯、时刻的TDI CCD离散平均焦面福照度；
[0106] S802,计算一个积分时间段内的TDI CCD时间平均焦面福照度;具体为:将通过步 骤S801得到的TDI CCD的每个CCD在各离散时刻的TDI CCD离散平均焦面福照度取平均，得 到一个积分时间段内的TDI CCD时间平均焦面福照度。
[0107] 本发明还提供了一种计算机可读介质，所述介质中存储有计算机指令，所述指令 执行前述任意一种方法。
[0108] 本发明W30厘米分辨率的Worldview-3全色立体影像作为仿真数据源，对50-80厘 米分辨率的立体测绘卫星影像进行成像仿真，在仿真链路中把平台参数、相机参数有机地 串联起来，实现了几何仿真和福射仿真的结合。
[0109] 本发明的有益效果是:本发明实现了基于甚高分辨率卫星Worldview-3立体影像 的分米级星载TDI CCD立体测绘相机成像仿真，模拟的两个相机的影像源于视角接近的输 入源立体影像。最大程度地避免了同一数据源条件下地面几何模型过度简化引起的立体影 像逼真度不高的问题。通过输入甚高分辨率卫星影像的像点=维坐标向TDI像面投射的方 法计算TDI影像子CCD的观测点，可W实现较快速的视点检测，避免了传统迭代方法在高程 突变位置不收敛的问题。对TDI CCD相机多级动态积分的过程进行严格模拟，实现了卫星轨 道、姿态数据，相机几何、福射参数，TDI CCD多级动态积分运一融合了严格几何转换和福射 转换过程的高精度模拟，保证了模拟结果的可信度。
【附图说明】
[0110] 图1为根据本发明实施例的、分米级星载TDI CCD立体测绘相机成像仿真方法的流 程图；
[0111] 图2为TDI CCD与细分CCD示意图。
【具体实施方式】
[0112] 下面参照附图对本发明进行更全面的描述，其中说明本发明的示例性实施例。
[0113] 如图1所示，一种分米级星载TDI CCD立体测绘相机成像仿真方法，其特征在于，包 括W下步骤：
[0114] 步骤Sl，获取基础数据，并对所述基础数据进行预处理，生成卫星立体影像每一点 的影像-几何模型。
[0115] 根据本发明的实施例，通过导入数据的方式获取基础数据，导入的基础数据为 Worldview-S卫星立体影像和辅助RPC(Rational化Iynomial Coefficient)参数。本发明 采用Worldview-S卫星立体影像数据，是因为Worldview-S卫星立体影像是目前商业应用立 体测绘卫星影像中分辨率最高的。其他分米级立体测绘卫星的分辨率一般0.5米-0.8米左 右。因此可W将其他设计分辨率低于Worldview-3的分米级立体测绘卫星影像仿真作为 Wor 1 d V i ew-3影像的退化处理。
[0116] 对基础数据进行预处理包括：
[0117] (Sl 1 )，通过多次前方交会和后方交会迭代计算，消除RPC参数的相对误差；
[0118] 具体为：首先通过SIFT匹配获取Worldview-3全色立体影像间的9个W上(9、16或 25,即3*3、4*4或5巧均匀分布)分布均匀的同名像点；通过同名像点前方交会获取物方地面 点坐标;执行后方交会，所述后方交会为W地面点和对应像点作为控制点校正RPC参数的分 子常数项和一次项。执行S次W上前方交会、后方交会迭代计算，消除RPC参数的相对误差。
[0119] (S12)，通过准核线立体影像制作和半全局匹配算法计算左右影像的视差图，根据 所述视差图计算同名像点坐标，进而前方交会计算准核线影像每一点的物方大地坐标化， B，H)，其中，L表示经度、B表示缔度、H表示楠球高。
[0120] 步骤S12具体可W为：通过投影轨迹法计算Worldview-3准核线立体影像和各自 RPC参数;通过半全局匹配算法计算左右影像的视差图，通过左右一致性检测剔除误匹配点 并通过距离倒数加权法拟合误匹配点的视差;根据视差图计算同名像点坐标，前方交会计 算所述准核线立体影像每一点的物方大地坐标化，B，H)，其中，L表示经度、B表示缔度、H表 示楠球高。
[0121] 通过W上预处理操作，Worldview-S准核线影像上的每个像点不仅包括灰度信息， 还包括用于跟仿真影像建立映射关系的几何位置信息。
[0122] 在立体影像的仿真过程中，每个视角的仿真相机对应的数据为观测角度接近的 Worldview-3准核线立体影像、RPC参数和准核线立体影像每一点的物方大地坐标。运是保 证立体影像仿真精度、体现立体仿真影像间几何、福射差异性的数据基础，因此，本发明采 用立体影像上每一点的"影像-几何模型"（即，每一点的物方大地坐标化，B，H)和影像灰度） 作为模拟数据源，W提高仿真精度，运是本发明的改进点之一。
[0123] 步骤S2,计算卫星轨道数据，所述卫星轨道数据包括地球固连坐标系下卫星的扫 描时刻、位置、速度数据。
[0124] 本发明通过步骤201-210所述的计算流程来计算卫星轨道数据，步骤201-210所述 的计算流程不是严格的动力学计算方法，而是一种简化的数学模型，运样可W更有效的计 算卫星轨道数据。
[0125] 轨道数据计算的目的是建立卫星本体坐标系的原点(恥，打，2。)在卫星运行的每一 时刻，在WGS84地球固连坐标系中的运行轨迹模型，卫星原点轨迹数据用卫星相机成像时间 段内每一时刻卫星本体坐标系原点的位置速度矢量(X。，Yd，Zd，Vxd，Vyd，Vzd )来表示。模拟卫 星轨道数据的频率是卫星下传轨道数据频率的10倍，模拟成像任一时刻的卫星位置、速度 矢量由该时刻前后各两组位置、速度矢量均通过拉格朗日插值获得。
[0126] 本发明用规则圆轨道来描述卫星轨道，卫星轨道数据的计算具体包括W下步骤：
[0127] 步骤S201，计算惯性系轨道数据。
[012引在惯性参考系的XOZ平面上，取一组卫星测量坐标系原点坐标，各原点坐标与惯性 参考系原点的距离相同，均为R，R等于平均地球半径加上平均轨道高度。假设第一个原点对 应时刻为0，相邻两个质屯、对应时刻相差0.1秒，由R可求得轨道周期T:
[0130] 其中G为万有引力常数，等于6.67*10-11米V(千克*秒2)，M为地球质量，等于5.98* 1〇24千克。
[0131] 降交点时刻Tdsd(秒)减去T/4就是第一个质屯、对应的实际时间To,即
[0132] To = Tdsd-T/4
[0133] W下取样时刻依次加0.1秒为相应的实际时间。d0 = 36O/(T*lO)为相邻两个原点 的地屯、角。第一个质屯、在Z轴上（Z = R),第i(i = 0，l，一)个取样时刻的卫星原点位置为： 义= i?xcos(0) 7 = 0 Z-Rx sin(0) (I)
[0134] 其中，目= (90-d 目 Xi)/180.0XJi。
[0135] 设轨道倾角为1，将取样坐标绕X轴顺时针旋转1-90度，然后绕巧由逆时针旋转L度。
[0136] 步骤S202,计算第i个时刻卫星原点(Xi，Yi，Zi)对应的初始大地坐标化i ,Bi，Hi)。 步骤S202具体包括：
[0137] 对第i个取样时刻，将卫星原点^1，￥1，2。（第1个取样时刻的=轴坐标)位置绕2 轴旋转i X化度。化=(1 /2400)°，其为0.1秒内地球自转角度：
[0139] 其中，
[0140] 计算第i个时刻卫星原点(Xi，Yi，Zi)对应的初始大地坐标化i ,Bi，Hi)，计算公式 为：
[0142] L = atan(Y/X)
[0143] H=PZcos(B)-N
[0144] 其中，（L，B，H)为化i，Bi，出)序列的任一组数据，（X，Y，Z)为对应(Xi，Yi，Zi)序列的 数据，（Xi，Yi，Zi)为第i个时刻卫星原点=维坐标。其中：
[01 加]0=atan(Z/P*a/b)
[0151] 步骤S203,计算经度