本發(fā)明屬于圖像仿真領(lǐng)域，涉及一種基于輻射幾何一體化設(shè)計(jì)的光學(xué)遙感圖像仿真方法。

背景技術(shù)：

遙感圖像仿真能夠?qū)b感器在空間環(huán)境條件運(yùn)行中遇到的各種狀態(tài)進(jìn)行模擬。通過仿真分析，設(shè)計(jì)人員可以在遙感器研制初期可靠地掌握有關(guān)信息，為提高遙感器性能，獲得用戶滿意的圖像而服務(wù)。

目前由于計(jì)算機(jī)仿真技術(shù)的快速發(fā)展，可以比較逼真地模仿各種目標(biāo)與背景的各種輻射特性，模擬各種照明條件，對(duì)不同目標(biāo)進(jìn)行各種數(shù)學(xué)仿真。美國(guó)柯達(dá)公司設(shè)計(jì)了遙感器仿真軟件，從1979年開始應(yīng)用共經(jīng)過了超過20000幅(次)成像驗(yàn)證。軟件考慮了成像過程中的15個(gè)環(huán)節(jié)，主要關(guān)注輻射傳輸鏈路的模擬，成像幾何關(guān)系語焉不詳。北京理工大學(xué)開發(fā)了圖像仿真系統(tǒng)，由實(shí)驗(yàn)室內(nèi)可見光—近紅外運(yùn)動(dòng)目標(biāo)生成裝置、動(dòng)態(tài)特性測(cè)試系統(tǒng)、光量子噪聲生成裝置等組成，但不是針對(duì)真實(shí)地面三維信息進(jìn)行遙感成像仿真。中科院長(zhǎng)春光機(jī)所開發(fā)了圖像仿真軟件系統(tǒng)，將選定的輸入圖像目標(biāo)經(jīng)過處理，考察在假定的條件下，圖像目標(biāo)的質(zhì)量下降(退化)程度，用于圖像質(zhì)量的預(yù)估。該軟件也僅考慮輻射傳輸?shù)那闆r。

隨著衛(wèi)星平臺(tái)機(jī)動(dòng)能力的增強(qiáng)，多角度全方位高時(shí)間分辨率的要求日益增加，成像幾何關(guān)系對(duì)圖像仿真的影響也越來越高，所以急需將遙感成像幾何輻射聯(lián)合仿真，才能更真實(shí)的反映成像條件，提高圖像仿真精度。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明所要解決的問題是：顧及成像幾何對(duì)圖像仿真的影響，提供一種基于輻射幾何一體化設(shè)計(jì)的光學(xué)遙感圖像仿真方法，提高光學(xué)遙感圖像仿真精度。

本發(fā)明的技術(shù)方案是：一種基于輻射幾何一體化設(shè)計(jì)的光學(xué)遙感圖像仿真方法，包括以下步驟：

1)確定成像起始時(shí)間t₀，成像時(shí)間t，時(shí)標(biāo)間隔Δt，采樣點(diǎn)數(shù)n＝t/Δt，n取整為自然數(shù)，第p個(gè)采樣點(diǎn)處對(duì)應(yīng)的時(shí)間為t₀+p×Δt；以時(shí)間為順序，在每個(gè)采樣點(diǎn)處仿真得到包括成像時(shí)側(cè)擺角、俯仰角和偏航角的姿態(tài)參數(shù)和成像時(shí)以攝影中心三維坐標(biāo)表示的軌道參數(shù)(X_i,Y_i,Z_i)，其中i＝1,2…n，為采樣點(diǎn)數(shù)；將每個(gè)姿態(tài)參數(shù)擬合為以時(shí)間為自變量，姿態(tài)參數(shù)為因變量的姿態(tài)仿真三次多項(xiàng)式；將每個(gè)軌道參數(shù)擬合為以時(shí)間為自變量，軌道參數(shù)為因變量的軌道仿真三次多項(xiàng)式；

2)確定光學(xué)系統(tǒng)焦距f，積分級(jí)數(shù)integral_level，積分時(shí)間integral_time,每一行探元數(shù)量為num_CCD，仿真圖像行數(shù)num_line＝t/(integral_level×integral_time)，num_line取整數(shù)，第l行成像時(shí)刻為t_l＝t₀+l×(integral_level×integral_time)；l行對(duì)應(yīng)的姿態(tài)數(shù)據(jù)和軌道數(shù)據(jù)(X_l,Y_l,Z_l)由t_l帶入步驟1)中的姿態(tài)仿真三次多項(xiàng)式和軌道仿真三次多項(xiàng)式獲得；

3)設(shè)定探元尺寸p，軌道高度H，地面采樣間隔GSD＝pH/f，成像范圍為num_line×num_CCD×GSD；

4)利用起始行軌道數(shù)據(jù)(X₁,Y₁,Z₁)，得到仿真場(chǎng)景范圍為X方向X₁-num_CCD×GSD/2到X₁+num_CCD×GSD/2，Y方向?yàn)閅₁到Y(jié)₁+num_line×GSD；對(duì)已有的地面數(shù)字高程模型在仿真場(chǎng)景范圍內(nèi)進(jìn)行分辨率為GSD的重采樣，得到仿真場(chǎng)景高程數(shù)據(jù)；

5)對(duì)第l行建立任一探元的成像幾何模型：

其中s＝1,2,……num_CCD，為探元編號(hào)，(ψ_x(s),ψ_y(s))為探元s在相機(jī)坐標(biāo)系下的指向角，λ＝H/f為攝影比例因子，

為第l行相方坐標(biāo)系到物方坐標(biāo)系的旋轉(zhuǎn)矩陣，(X_g,Y_g,Z_g)為探元s物方坐標(biāo)系下的坐標(biāo)，(x₀,x₁,x₂,x₃,y₀,y₁,y₂,y₃)為內(nèi)定標(biāo)參數(shù)模型的系數(shù)；

第l行任一探元s，將仿真場(chǎng)景高程數(shù)據(jù)平均值為Z₀代入成像幾何模型計(jì)算得到該探元初始地面點(diǎn)物方平面坐標(biāo)(X_g',Y_g')，在仿真場(chǎng)景高程數(shù)據(jù)中讀取(X_g',Y_g')對(duì)應(yīng)的高程Z_g'，再由成像幾何模型計(jì)算獲得新的平面坐標(biāo)(X_g",Y_g")，在仿真場(chǎng)景高程數(shù)據(jù)中讀取(X_g",Y_g")對(duì)應(yīng)的高程Z_g"，迭代計(jì)算直到本次高程計(jì)算結(jié)果與上一次高程計(jì)算結(jié)果兩次計(jì)算的到的高程之差小于0.1米，計(jì)算結(jié)束，將本次計(jì)算結(jié)果作為得到探元s所對(duì)應(yīng)的物方坐標(biāo)(X_g,Y_g,Z_g)；

遍歷所有探元得到每個(gè)探元對(duì)應(yīng)的大小為num_line×num_CCD×3的物方坐標(biāo)三維矩陣RG；

6)以GSD大小為采樣間隔實(shí)地測(cè)得仿真場(chǎng)景中各點(diǎn)三維點(diǎn)的反射率，得到仿真場(chǎng)景反射率數(shù)據(jù)，按照RG中每個(gè)點(diǎn)的三維坐標(biāo)重采樣得到其對(duì)應(yīng)的反射率值，根據(jù)地理位置獲取太陽福照度及大氣透過率，太陽福照度乘以大氣透過率再乘以反射率，得到各點(diǎn)入瞳福亮度Ri，輸出仿真場(chǎng)景入瞳輻亮度數(shù)據(jù)；

7)給定探測(cè)器最大輻亮度Rmax和最小輻亮度Rmin，雜光系數(shù)SLC，漸暈系數(shù)VC，面遮攔因子BC，光學(xué)系統(tǒng)平均透過率PR，大氣MTF矩陣aMTF，光學(xué)系統(tǒng)MTF矩陣oMTF，仿真場(chǎng)景入瞳輻亮度數(shù)據(jù)中各點(diǎn)值Ri到探測(cè)器的能量QP＝(Ri×PR×SLC×VC×BC)*aMTF*oMTF，其中*為卷積符號(hào)，若QP≥Rmax，則QP＝Rmax，若QP≤Rmin，則QP＝Rmin；

8)給定光電轉(zhuǎn)換效率LER，平均量子效率QE，探測(cè)器MTF矩陣dMTF；探測(cè)器獲得的電子數(shù)量NE＝(QP×LER×QE)*dMTF，設(shè)定飽和電子數(shù)fNE，若NE≥fNE，則NE＝fNE，若NE＜fNE，則NE＝(QP×LER×QE)*dMTF；設(shè)定光譜響應(yīng)函數(shù)fs(NE)，量化噪聲矩陣QNoise，探測(cè)器獲得的DN值DN＝fs(NE)*QNoise；

9)輸出結(jié)果圖像。

本發(fā)明與現(xiàn)有技術(shù)相比的優(yōu)點(diǎn)在于：首先建立探測(cè)器各探元成像幾何模型，由幾何模型獲取成像位置三維坐標(biāo)，再由三維坐標(biāo)確定該探元對(duì)應(yīng)的仿真場(chǎng)景反射率，考慮太陽福照度和大氣透過率，在經(jīng)過光學(xué)仿真系統(tǒng)和探測(cè)器仿真系統(tǒng)得到仿真圖像。同時(shí)考慮幾何成像和輻射成像關(guān)系，提高空間遙感光學(xué)成像仿真精度。

附圖說明

圖1為本發(fā)明的操作流程圖；

具體實(shí)施方式

下面結(jié)合附圖對(duì)本發(fā)明做進(jìn)一步詳細(xì)描述，結(jié)合圖1本發(fā)明流程圖，概括起來，本方法的實(shí)施可以分為三個(gè)步驟：

1.成像幾何建模

1.1實(shí)施例中確定成像起始時(shí)間t₀＝0，成像時(shí)間t＝144s，時(shí)標(biāo)間隔Δt＝1s，采樣點(diǎn)數(shù)n＝t/Δt＝144，則生成的姿態(tài)文件包含144組姿態(tài)參數(shù)生成的軌道文件包含144組軌道參數(shù)(X_i,Y_i,Z_i)，其中i＝1,2…144；將144個(gè)姿態(tài)參數(shù)擬合為以時(shí)間為自變量，姿態(tài)參數(shù)為因變量的姿態(tài)仿真三次函數(shù)；將144個(gè)軌道參數(shù)擬合為以時(shí)間為自變量，軌道參數(shù)為因變量的軌道仿真三次函數(shù)；

1.2探元尺寸p＝10μm，軌道高度H＝500km，焦距f＝5m，則地面采樣間隔GSD＝pH/f＝1m。圖像積分級(jí)數(shù)為integral_level＝48，integral_time＝0.3ms，每一行探元數(shù)量為10000個(gè)，則仿真圖像行數(shù)num_line＝t/(integral_level×integral_time)＝10000。

第l行成像時(shí)刻為t_l＝l×(48×0.3ms)；l行對(duì)應(yīng)的姿態(tài)數(shù)據(jù)和軌道數(shù)據(jù)(X_l,Y_l,Z_l)由t_l帶入步驟1.1中的姿態(tài)仿真三次函數(shù)和軌道仿真三次函數(shù)獲得；

成像范圍為num_line×num_CCD×GSD＝10000m×10000m。

1.3利用起始行軌道數(shù)據(jù)(X₁,Y₁,Z₁)中的沿軌方向值Y₁和Y₁+10000×GSD作為仿真場(chǎng)景范圍Y方向的上下限，X₁-5000×GSD和X₁+5000×GSD作為仿真場(chǎng)景范圍X方向的上下限，得到仿真場(chǎng)景范圍。對(duì)已有的地面數(shù)字高程模型在仿真場(chǎng)景范圍內(nèi)進(jìn)行分辨率為GSD的重采樣，得到仿真場(chǎng)景高程數(shù)據(jù)。

1.4由式建立第l行每一個(gè)探元的成像幾何模型。其中s為探元號(hào)，為探元s在相機(jī)坐標(biāo)系下的指向角，λ＝H/f為攝影比例因子，

為第l行相方坐標(biāo)系到物方坐標(biāo)系的旋轉(zhuǎn)矩陣，由第l行姿態(tài)數(shù)據(jù)計(jì)算得到，(X_g,Y_g,Z_g)為探元s物方坐標(biāo)系下的坐標(biāo)，(X_l,Y_l,Z_l)第l行軌道數(shù)據(jù)，(x₀,x₁,x₂,x₃,y₀,y₁,y₂,y₃)為內(nèi)定標(biāo)參數(shù)模型的系數(shù)。

第l行任一探元s，將仿真場(chǎng)景高程數(shù)據(jù)平均值為Z₀代入成像幾何模型計(jì)算得到該探元初始地面點(diǎn)物方平面坐標(biāo)(X_g',Y_g')，在仿真場(chǎng)景高程數(shù)據(jù)中讀取(X_g',Y_g')對(duì)應(yīng)的高程Z_g'，再由成像幾何模型計(jì)算獲得新的平面坐標(biāo)(X_g",Y_g")，在仿真場(chǎng)景高程數(shù)據(jù)中讀取(X_g",Y_g")對(duì)應(yīng)的高程Z_g"，迭代計(jì)算直到本次高程計(jì)算結(jié)果與上一次高程計(jì)算結(jié)果兩次計(jì)算的到的高程之差小于0.1米，計(jì)算結(jié)束，將本次計(jì)算結(jié)果作為得到探元s所對(duì)應(yīng)的物方坐標(biāo)(X_g,Y_g,Z_g)；

遍歷所有探元得到每個(gè)探元對(duì)應(yīng)的物方坐標(biāo)三維矩陣RG(10000×10000×3)。

2.成像輻射建模

2.1以GSD大小為采樣間隔實(shí)地測(cè)得仿真場(chǎng)景中各點(diǎn)三維點(diǎn)的反射率，得到仿真場(chǎng)景反射率數(shù)據(jù)，按照RG中每個(gè)點(diǎn)的三維坐標(biāo)雙三次重采樣得到其對(duì)應(yīng)的反射率值，根據(jù)地理位置獲取太陽福照度及大氣透過率，太陽福照度乘以大氣透過率再乘以反射率，得到各點(diǎn)入瞳福亮度Ri，輸出仿真場(chǎng)景入瞳輻亮度數(shù)據(jù)。

2.2給定探測(cè)器最大輻亮度Rmax＝75和最小輻亮度Rmin＝8，雜光系數(shù)SLC＝0.04，漸暈系數(shù)VC＝0.2，面遮攔因子BC＝0，光學(xué)系統(tǒng)平均透過率PR＝0.35，大氣MTF矩陣aMTF，光學(xué)系統(tǒng)MTF矩陣oMTF，仿真場(chǎng)景入瞳輻亮度數(shù)據(jù)中各點(diǎn)值Ri到探測(cè)器的能量QP＝(Ri×PR×SLC×VC×BC)*aMTF*oMTF，其中*為卷積符號(hào)，若QP≥Rmax，則QP＝Rmax，若QP≤Rmin，則QP＝Rmin。

2.3給定光電轉(zhuǎn)換效率LER＝11×10^-6，平均量子效率QE＝0.25，探測(cè)器MTF矩陣dMTF。探測(cè)器獲得的電子數(shù)量NE＝(QP×LER×QE)*dMTF，設(shè)定飽和電子數(shù)fNE＝86，若NE≥fNE，則NE＝fNE，若NE＜fNE，則NE＝(QP×LER×QE)*dMTF，設(shè)定光譜響應(yīng)函數(shù)fs(NE)，量化噪聲矩陣QNoise，探測(cè)器獲得的DN值DN＝fs(NE)*QNoise。

2.4按照10000×10000每個(gè)像素的排列順序輸出DN值文件，即為仿真圖像。

本發(fā)明中未作詳細(xì)描述的內(nèi)容屬于本領(lǐng)域?qū)I(yè)技術(shù)人員的公知技術(shù)。