受災區域衛星遙感監測的條帶分解方法研究

許友善，林廣宏，趙 晗，謝亞恩

（1.中國船舶集團有限公司第七二二研究所，湖北 武漢 430200；2.中國船舶及海洋工程設計研究院，上海200011；3.哈爾濱工程大學航天與建筑工程學院，黑龍江 哈爾濱 150001）

引言

我國是自然災害多發的國家之一，其災害種類多、分布廣、發生頻率高、受災損失嚴重。在軌衛星具有觀測頻次高、遙感范圍廣等特點，在災害監測評估中發揮著重要的作用［6］。衛星遙感一直是森林草原火災、雪災、洪澇、地震等自然災害的災后監測和評估的重要手段。

由于受災區域在災后的監測條件復雜，傳統的地面監測、氣象衛星監測等手段很難實現對受災區域的全面、高精度的監測，使得災后救援難以有針對性地進行。因此，可以充分利用遙感衛星對受災區域進行監測［1-5］。由于受災區域面積較大，衛星幅寬有限，因此，用遙感衛星對受災區域進行監測需要首先對受災區域進行區域分解，使其能被衛星多次推掃覆蓋。傳統的衛星在軌觀測任務規劃通常為點目標規劃，隨著敏捷衛星的發展，衛星在軌機動能力的提升為大區域的條帶拼幅成像提供了技術支持［7］。因此文中針對受災區域災后監測，給出了基于條帶分解的敏捷衛星災后快速監測方法，能夠有效分析單顆敏捷衛星對受災區域的覆蓋觀測時長，為各種自然災害的災后救援工作提供了數據支撐。

1 數學模型

1.1 集合和參數設計

定義T={v1,v2,...,vN}表示待觀測的受災區域，vi表示目標區域的第i個頂點，N表示目標區域的頂點數量；S={Si|Si=(a,e,i,Ω,ω,tp)}，表示衛星的集合，每個衛星對應一組軌道六要素，a,e,i,Ω,ω,tp分別表示衛星的軌道半長軸、軌道偏心率、軌道傾角、升交點赤經、近地點幅角、近地點時刻；LLA={LLAi|LLAi=(loni,lati),1 ≤i≤Norbit}表示預報的衛星軌道星下點經緯度數據集合，Norbit表示預報的衛星軌道點數；PV={PV1,PV2,...,PVNorbit}表示預報的衛星軌道位置和速度數據。Width表示衛星的視場幅寬；Cover表示相鄰條帶邊緣所需的重疊寬度。

1.2 變量設計

m表示進行條帶劃分時，將目標區域沿垂直于星下點的方向均等分成m個條帶；Strip={strip1,strip2,...,stripm}表示劃分后的條帶集合；stripj={vs,ve}表示第j個條帶劃分后得到的觀測起始點和結束點；W={W1,W2,...,Wm} 表示計算得到的所有條帶的成像時間窗口；Wj={tjs,tje} 表示衛星對條帶stripj的實際成像時間窗口，tjs、tje分別表示實際成像時間窗口的開始時刻、結束時刻；TM表示2次條帶成像任務間機動時長；γj,k為二進制變量，表示條帶stripj的執行情況，γj,k= 1表示條帶stripj在第k軌被執行，否則γj,k= 0［8］。

2 條帶分解

2.1 計算條帶分解的參考直線

計算目標區域每個頂點與同緯度下每一軌衛星星下點之間的經度差，第k軌的經度差之和計算如下。

則可以得到所有軌道對應的經度差之和的集合{sum1,sum2,...,sumk,...,summ} ，sumk的最小值對應的軌道K即為參考軌道。

確定參考軌道后，可以根據目標區域所有頂點與參考軌道的經度差{Δlon1,Δlon2,...,ΔlonN} 判斷區域與參考軌道之間的位置關系。若所有Δlon均大于零，則區域在參考軌道右側且不與參考軌道相交，若所有Δlon均小于零，則區域在參考軌道左側且不與參考軌道相交，若Δlon既有大于零的值也有小于零的值，則參考軌道穿過受災區域，若有Δlon等于零，則表明有頂點位于參考軌道之上。

在區域所有頂點中尋找緯度最大值latmax和最小值latmin，則可以得到中間緯度。

在參考軌道第K軌中尋找latmiddle對應的星下點LLAmiddle作為條帶分解的參考點，如圖1中的點LLAmiddle。

圖1 參考直線示意圖Fig.1 Reference line of aera

找到LLAmiddle兩側的星下點LLA1、LLA2，將其轉換到高斯坐標系得到P1=(x1,y1)、P2=(x2,y2（)后文中提到的x,y坐標均為高斯坐標系下坐標）則可以得到參考直線的斜率。

因此條帶分解的參考直線可以寫為：

2.2 確定條帶邊界線

確定條帶邊界線首先需要計算邊界的平移量Δx，再根據邊界平移量確定最終的條帶數量m。圖2 給出了條帶邊界平移量示意圖，根據給定的衛星幅寬Width計算得到條帶邊界的平移量如下。

圖2 條帶邊界平移量示意圖Fig.2 Displacement of strip edge

計算所有頂點到參考直線的距離，找到距離最小值dmin和最大值dmax，從而計算得到條帶的數量如下

確定條帶邊界線首先需要檢查頂點順序并重新排列頂點。若參考直線不穿過受災區域，則判斷v1是否為距離參考直線最近的點，若不是則找到最近的頂點vi作為第一個頂點并將原來的順序進行如下調換（原頂點是按照順時針順序給出的）：

若參考直線穿過受災區域，則參考圖3，首先找到參考直線兩側距離最遠的2個頂點v1,v4，在兩點之中找到一個距離相對較近的點v4作為第一個點，同樣按照上述方式對頂點順序進行調整。

圖3 參考直線穿過受災區域示意圖Fig.3 Reference line crossing disaster aera

在得到邊界平移量和第一個頂點之后，即可得到每個條帶邊界的直線，此處v1指的是調整順序后得到的第一個頂點。Δx的正負按照左加右減的原則確定，例如受災區域在參考直線左側，或參考直線穿過受災區域且第一個頂點在直線右側，條帶邊界線從第一個頂點開始依次向左平移Δx，則Δx為正。同理，若受災區域在參考直線右側，或參考直線穿過受災區域且第一個頂點在直線左側，條帶邊界線從第一個頂點開始依次向右平移Δx，則Δx為負。第一個條帶邊界線穿過第一個頂點，平行于參考直線，直線方程如下：

第2個條帶邊界線方程如下：

第j個條帶邊界線方程如下：

由于劃分共得到m個條帶，則條帶邊界共有m+ 1個。

2.3 確定條帶的觀測起始點和結束點

條帶的觀測起始點和結束點確定的垂直于參考直線的2 條邊界線我們稱其為起始邊線和結束邊線，條帶的起始邊線和結束邊線可能經過條帶中存在的受災區域頂點，也可能經過條帶邊界線與受災區域邊界線的某一個交點，因此要確定條帶的起始邊線和結束邊線，首先需要找到條帶中兩端邊界線所有可能的經過的點。

首先找到條帶中存在的受災區域的頂點。依次判斷每個條帶與每個頂點之間的位置關系，得到每個條帶中的頂點，若頂點位于某一條帶邊界線上，則兩側條帶所考慮的受災區域頂點中均包括該頂點，最終可以得到每個條帶中的頂點集合Vertex={vertex1,vertex2,...,vertexm} 。其中，vertexj={v1,...,vi} ，vertexj也可能為空集，即條帶中不存在頂點。

在確定條帶中存在的受災區域的頂點后，再找到所有區域邊界線與條帶邊界線之間的交點。根據集合Vertex={vertex1,vertex2,...,vertexm} 可以確定每個頂點所在的條帶，假設2個相鄰的頂點vi和vi+1分別位于條帶s1和s2中，則頂點vi和vi+1之間的條帶邊界線為

式中(s1,s2]表示條帶邊界線s1+ 1,s1+ 2,...,s2，2 個相鄰頂點vi和vi+1之間共s2-s1個條帶邊界線。因此與受災區域邊界edgei有交點的條帶邊界線集合={strips1+1,...,strips2} 。

計算受災區域邊界edgei與條帶邊界線={strips1+1,...,strips2} 的交點，記錄到每個條帶的頂點集合Vertex={vertex1,vertex2,...,vertexm} 中。則頂點集合中包括每個條帶中存在的受災區域頂點和條帶邊界線與受災區域邊界的交點。

作一條過頂點v1并垂直于參考直線l0的直線L，直線方程如下：

計算每個條帶中的所有頂點qp到直線L的距離，如圖4所示：

式中，d的值可能小于零。

從大到小（考慮d的正負）對計算得到的某一條帶中的所有d進行排序，找到d的最大值和最小值對應的2 個點，即為條帶起始和結束兩端直線需要經過的點。以圖4 中的條帶strip2為例，計算條帶中的4個頂點與直線L的距離，由于頂點位于L兩側，因此d有正負，此處正負只表示方位。對4個頂點的d進行排序，則可以確定條帶最上方和最下方的2個頂點為v5、q3。分別過v5、q3作垂直于l0的直線，并計算直線與條帶邊界線的交點，得到條帶的4個頂點{point1,point2,point3,point4} ，取條帶的起始邊線和結束邊線的中點為條帶兩端的觀測起始點和結束點，則條帶的觀測起始點和結束點可以表示為

圖4 條帶兩端頂點計算示意圖Fig.4 Vertex calculation figure of strip

最后判斷星下點的走向，以圖1 中的LLAmiddle為t0時刻，則LLA1={lon1,lat1} 、LLA2={lon2,lat2} 分別對應t0- 1和t0+ 1時刻。

若lat1＜lat2，則衛星上行，則選取point01、point02中緯度較小的點為觀測起始點，另一點為觀測結束點，反之衛星下行，選取point01、point02中緯度較大的點為起始點，另一點為觀測結束點。由于上述計算是在高斯坐標系中進行的，因此需要將得到的觀測起始點和結束點轉換為經緯度，得到條帶stripj的觀測起始點和結束點stripj={vs,ve}。最終計算m個條帶的觀測起始點和結束點，得到Strip={strip1,strip2,...,stripm}用作成像任務規劃的輸入。

3 成像規劃

計算條帶stripj的觀測起始點的開始可見時刻，并作為條帶stripj的實際觀測開始時刻tjs。根據衛星的地速vg計算條帶stripj的觀測時長，具體計算如下：

式中P1,P2分別表示條帶stripj的起始觀測點和結束觀測點在地固坐標系中的坐標。因此，當條帶的實際觀測開始時刻確定后，其實際觀測的結束時刻為：

衛星的可見性如圖5所示，本文假設衛星是敏捷衛星，軌道偏心率為0，衛星最近可對星下點進行觀測，最遠可觀測到與衛星地球表面的切點。衛星可見性約束可以表述為0 ＜θ＜θmax，其中RE、RS分別表示地球半徑和衛星半長軸。

圖5 衛星可見性示意圖Fig.5 Satellite visibility figure

根據上述可見性約束，依次計算條帶Strip={strip1,strip2,...,stripm}的可見性。

遍歷衛星的LLA數據，計算衛星與目標間的地心角θ，若0 ＜θ＜θmax，則當前時刻衛星對目標可見，并以當前時刻作為條帶stripj的實際觀測開始時刻tjs，并可以通過式（17）計算得到條帶stripj成像結束時刻。重復上述方法依次計算每個條帶的實際成像時間。但是在計算每個條帶的實際成像時間時，還需要考慮任務間的沖突，即條帶stripl的成像時段對于已規劃的任務序列中的任意其他條帶stripj應滿足tls-TM≥tje或tle+TM≤tjs。

最終，遍歷所有條帶，計算得到整個受災區域的所有條帶的實際成像時間窗口W={W1,W2,...,Wm} 。

4 仿真驗證

4.1 仿真場景設計

表1 給出了仿真用到的衛星軌道參數。衛星的軌道預報時長為24 h，步長為1s，圖6 中的曲線即為衛星的所有星下點軌跡。衛星的視場幅寬Width= 60 km，相鄰條帶邊緣所需的重疊寬度Cover= 2 km。

表1 衛星軌道參數Table 1 Orbital parameters of satellite

（1）仿真場景1

考慮到四川地區自然災害頻發，因此在四川附近地區選取一個五邊形區域作為受災區域1，表2給出了目標區域1的頂點，圖6給出了受災區域1與軌道的示意圖。

表2 目標區域1頂點Table 2 Vertex of disaster aera 1

圖6 軌道與受災區域1示意圖Fig.6 Orbit and disaster aera 1

圖7 軌道與受災區域2示意圖Fig.7 Orbit and disaster aera 2

（2）仿真場景2

為驗證算法在處理受災區域與衛星下點不同位置關系下的準確性，設計一種受災區域位于星下點兩側的仿真場景，表3給出了目標區域的頂點，圖6給出了受災區域與軌道的示意圖。

表3 目標區域2頂點Table 3 Vertex of disaster aera 2

4.2 仿真

（1）仿真場景1計算結果

圖8 給出了條帶計算結果，圖9 給出了條帶實際成像時段規劃結果示意圖。在不考慮軌道機動的情況下，以本文算例給出的軌道初始條件，得到的所有條帶的最早成像時刻為19842 s（5.511 h），成像結束時刻為21 003 s（5.834 h），從開始成像到成像結束共耗時1 161 s（0.323 h），其中成像時長共1 004 s（0.279 h），任務間的機動及等待成像共耗時157 s。

圖8 受災區域1條帶分解結果Fig.8 Strips splitting result of disaster aera 1

圖9 受災區域1成像時段Fig.9 Imaging time result of disaster aera 1

（2）仿真場景2結果

圖10 給出了條帶計算結果，圖11 給出了條帶實際成像時段規劃結果示意圖。在不考慮軌道機動的情況下，以本文算例給出的軌道初始條件，得到的所有條帶的最早成像時刻為19 991 s（5.553 h），成像結束時刻為20 951 s（5.820 h），從開始成像到成像結束共耗時960 s（0.267 h），其中成像時長共820 s（0.228 h），任務間的機動及等待成像共耗時140 s。

圖10 受災區域2條帶分解結果Fig.10 Strips splitting result of disaster aera 2

圖11 受災區域2成像時段Fig.11 Imaging time result of disaster aera 2

5 結論

針對災后小范圍區域觀測，氣象衛星則難以滿足精度要求的問題，文中給出了基于條帶分解的受災區域衛星遙感觀測方法，提供了災后小范圍區域的觀測思路。通過將受災區域進行條帶分解，使敏捷遙感衛星能夠對受災區域進行單軌連續觀測或多軌連續觀測，從而達到對受災區域的災后應急觀測的目的。本文通過仿真算例，驗證了本文給出的基于條帶分解的區域目標的成像時段的規劃方法，算法可行有效，為災后小范圍區域的災后監視提供了參考。