楊云飛，李金平,2

(1.云南師范大學(xué)地理學(xué)部，昆明 650500；2.西部資源環(huán)境地理信息技術(shù)教育部工程研究中心，昆明 650500)

近年來(lái)城市爆炸時(shí)常發(fā)生，造成了巨大的人員傷亡和財(cái)產(chǎn)損失，如天津港爆炸[1]、響水爆炸[2-4]和貝魯特爆炸[5]等。差分干涉測(cè)量(Differential InSAR，D-InSAR)技術(shù)是目前國(guó)內(nèi)外在InSAR應(yīng)用上最為成熟的技術(shù)，這項(xiàng)技術(shù)通過(guò)星載合成孔徑雷達(dá)監(jiān)測(cè)地球表面形變精度在厘米級(jí)甚至毫米級(jí)。Gabriel等[6]在1989年首次驗(yàn)證了D-InSAR技術(shù)具有監(jiān)測(cè)地表微小形變的能力；Goldstein等[7]在1993年把D-InSAR技術(shù)應(yīng)用于南極的冰川移動(dòng)監(jiān)測(cè)并獲得成功，成果發(fā)表于Science。同年Massonnet等[8]利用D-InSAR技術(shù)獲取了1992年加利福尼亞Landers地震同震形變場(chǎng)。此后，D-InSAR技術(shù)被廣泛應(yīng)用于由各種原因引發(fā)的地表形變的監(jiān)測(cè)中，包括地震、滑坡、火山運(yùn)動(dòng)、冰川移動(dòng)等自然現(xiàn)象[9-14]，還包括礦山開(kāi)采、地下水抽取等人類活動(dòng)[15-18]。

本文采用D-InSAR技術(shù)對(duì)爆炸產(chǎn)生地表形變場(chǎng)進(jìn)行監(jiān)測(cè)，同時(shí)對(duì)比天然地震產(chǎn)生形變場(chǎng)，得到爆炸產(chǎn)生的形變場(chǎng)特點(diǎn)，為爆炸后的重建提供形變信息支持；并基于形變值與薩道夫斯基公式構(gòu)建城市爆炸形變模型。

1 研究區(qū)概況

天津港爆炸中心位于東經(jīng)117.74°、北緯39.04°，從衛(wèi)星影像(見(jiàn)圖1)看，爆炸約554 m范圍內(nèi)設(shè)施嚴(yán)重?fù)p毀。響水爆炸位于東經(jīng)119.76°、北緯34.34°，爆炸340 m范圍內(nèi)建筑損毀嚴(yán)重(見(jiàn)圖2)。貝魯特爆炸位于東經(jīng)35.52°、北緯33.90°貝魯特港口內(nèi)，爆炸范圍約553 m內(nèi)港口設(shè)施摧毀殆盡，在約943 m范圍內(nèi)樓房有明顯震碎物品墜落(見(jiàn)圖3)。作為對(duì)比的天然地震選取2017年伊朗Mw7.3，其震中位于伊朗-伊拉克邊界歐亞板塊與阿拉伯板塊碰撞交匯處的扎格羅斯褶皺沖斷帶，構(gòu)造活動(dòng)強(qiáng)烈地震多發(fā)(根據(jù)地球物理中板塊運(yùn)動(dòng)模型NUVEL -1A計(jì)算結(jié)果顯示，阿拉伯板塊以30 mm/yr撞向歐亞板塊[19])；從地貌特征來(lái)看，震中區(qū)域植被稀疏，利于合成孔徑雷達(dá)進(jìn)行干涉測(cè)量；以上兩個(gè)因素也造就了此區(qū)域作為InSAR技術(shù)監(jiān)測(cè)地震引發(fā)形變較為理想的試驗(yàn)場(chǎng)(見(jiàn)圖4)。

圖1 天津港爆炸衛(wèi)星影像Fig.1 Satellite image of Tianjin Port explosion

圖2 響水爆炸衛(wèi)星影像Fig.2 Satellite image of Xiangshui explosion

圖3 貝魯特爆炸衛(wèi)星影像Fig.3 Satellite image of Beirut explosion

2 數(shù)據(jù)及數(shù)據(jù)處理

Sentinel-1A數(shù)據(jù)來(lái)自于歐州航天局于2014-04-03發(fā)射的全新C波段SAR衛(wèi)星， 它可以提供具有周期性且精確的大范圍形變監(jiān)測(cè)，性能比ERS-1/2、ENVISAT ASAR衛(wèi)星獲取的SAR影像更為完善，具有噪聲小、影像干涉能力更強(qiáng)、監(jiān)測(cè)精度更高、周期更短和免費(fèi)申請(qǐng)獲取等特點(diǎn)[20]。本文分別選取覆蓋該爆炸、地震事件發(fā)生前后區(qū)域的一組 Sentinel-1A影像數(shù)據(jù)(SLC格式)，通過(guò)D-InSAR技術(shù)來(lái)獲取爆炸產(chǎn)生形變場(chǎng)和地震產(chǎn)生形變場(chǎng)，影像參數(shù)如表1所示。

表1 研究區(qū)影像參數(shù)

D-InSAR技術(shù)其數(shù)據(jù)處理流程可分為：數(shù)據(jù)輸入、干涉圖生成、濾波和相干性計(jì)算、去平地效應(yīng)及相位解纏、軌道精煉和重去平、相位轉(zhuǎn)形變及地理編碼[20]。鑒于軌道存在的系統(tǒng)誤差、平地效應(yīng)及不同地物具有不同相干性的實(shí)際特點(diǎn)，本次數(shù)據(jù)處理過(guò)程中通過(guò)引入外部數(shù)據(jù)和多次調(diào)整參數(shù)，獲得了較為可靠結(jié)果，具體如下：數(shù)據(jù)輸入中通過(guò)添加衛(wèi)星精密軌道數(shù)據(jù)AUX_POEORB降低系統(tǒng)誤差產(chǎn)生的影響提高數(shù)據(jù)處理精度；干涉噪聲先采用距離向比方位向?yàn)?∶4的多視處理，之后選擇高通濾波(Goldstein)去除相位噪聲；利用STRM數(shù)據(jù)(30 m×30 m)外部數(shù)字高程模型(DEM)去除地形相位產(chǎn)生的影響；采用最小費(fèi)用流法進(jìn)行相位解纏，將纏繞的相位校正到高相干性參考點(diǎn)上，更好地還原地表形變引起的真實(shí)相位變化；在估算和去除殘余的恒定相位后，軌道誤差依然存在，由此選取GCP控制點(diǎn)，對(duì)干涉圖進(jìn)行重新估算，減小誤差；經(jīng)過(guò)以上處理獲取的純凈相位信息轉(zhuǎn)化為形變位移信息(包括生成雷達(dá)視線方向形變、垂直方向形變)并進(jìn)行地理編碼得到實(shí)際應(yīng)用中對(duì)應(yīng)的地表形變。

3 爆炸形變特征分析

天津港爆炸中心附近分別出現(xiàn)了一個(gè)較大沉降和隆升區(qū)域，其東向出現(xiàn)兩個(gè)較大的沉降漏斗區(qū)域，整體監(jiān)測(cè)區(qū)域垂直向形變范圍在-3.7 ～ 4.8 cm(見(jiàn)圖5)；響水爆炸中心附近出現(xiàn)了一個(gè)較小沉降區(qū)域但整體主要表現(xiàn)為隆升，其東向形變整體較西向劇烈，整體監(jiān)測(cè)區(qū)域垂直向形變范圍在-4.5～3.8 cm(見(jiàn)圖6)；貝魯特港爆炸中心附近主要表現(xiàn)為沉降但有小片區(qū)域表現(xiàn)為隆升，其東向主要表現(xiàn)為隆升西向主要表現(xiàn)為沉降，整體監(jiān)測(cè)區(qū)域垂直向形變范圍在-4.9～1.8 cm(見(jiàn)圖7)。從天津港爆炸、響水爆炸及貝魯特港爆炸雷達(dá)視線向形變場(chǎng)和垂直形變場(chǎng)綜合來(lái)看：爆炸中心周邊出現(xiàn)有密集且不連續(xù)(隆升與沉降交錯(cuò))的形變，結(jié)合爆炸前后衛(wèi)星影像地表隆升與沉降交錯(cuò)區(qū)域多發(fā)于建筑物密集且損毀嚴(yán)重區(qū)域，符合實(shí)際中建筑倒塌現(xiàn)象。由形變判斷發(fā)生建筑損毀區(qū)域，天津港爆炸損毀范圍半徑約在2 034 m；響水爆炸損毀半徑約在830 m；貝魯特港爆炸損毀半徑約在2 100 m。通過(guò)對(duì)比天然地震產(chǎn)生形變場(chǎng)(見(jiàn)圖8)，地表爆炸所引發(fā)的地表形變更多的表現(xiàn)出一種雜亂的不連續(xù)波動(dòng)。

圖6 響水爆炸形變場(chǎng)Fig.6 Deformation field of Xiangshui explosion

圖7 貝魯特港口爆炸形變場(chǎng)Fig.7 Deformation field of Beirut Port explosion

圖8 地震形變場(chǎng)Fig.8 Earthquake deformation field

4 城市爆炸形變模型構(gòu)建

綜合爆炸前后遙感影像目視解譯及爆炸形變場(chǎng)(InSAR)獲取信息，如表2所示。由表2發(fā)現(xiàn)從遙感影像目視解譯能較為準(zhǔn)確地獲取爆炸源周邊建筑物倒塌的范圍(范圍內(nèi)損壞等級(jí)為完全破壞)；綜合InSAR技術(shù)監(jiān)測(cè)到的爆炸形變范圍，可以獲取形變邊界至建筑物倒塌邊界范圍內(nèi)的建筑物受損情況(范圍內(nèi)損壞等級(jí)為次嚴(yán)重破壞至基本無(wú)破壞)。

表2 遙感方式獲取爆炸影響信息

由于爆炸發(fā)生在日常生產(chǎn)運(yùn)輸?shù)牡胤?化工園區(qū)、港口)可認(rèn)為爆炸是在普通地面發(fā)生，在此前提下爆炸產(chǎn)生的空氣沖擊波超壓可以運(yùn)用薩道夫斯基公式[21]計(jì)算:

(1)

式中：Δp為空氣沖擊波的超壓值,kPa；m為T(mén)NT炸藥的質(zhì)量，kg；R為爆炸半徑，m。

由公式(1)結(jié)合表2遙感手段獲取的邊界范圍信息計(jì)算空氣沖擊波的超壓值(見(jiàn)表3)。

表3 爆炸在遙感獲取邊界范圍下超壓值

在已知爆炸當(dāng)量情況下，通過(guò)遙感手段可快速獲取損毀范圍和形變范圍，即在較短的時(shí)間里獲得損毀嚴(yán)重范圍和次損毀嚴(yán)重范圍內(nèi)空氣沖擊波的超壓值。

由于D-InSAR形變監(jiān)測(cè)具有臨界形變梯度[22]，由于臨近爆炸中心，破壞等級(jí)在完全破壞范圍內(nèi)存在大量區(qū)域真實(shí)形變超過(guò)D-InSAR形變監(jiān)測(cè)能力范圍，因此本文在城市爆炸形變模型構(gòu)建中不考慮完全破壞區(qū)域；通過(guò)獲取的次嚴(yán)重區(qū)域至基本不受影響區(qū)域內(nèi)空氣沖擊波的超壓值及范圍內(nèi)形變值，以形變值作為自變量，超壓值作為因變量，通過(guò)分析3次爆炸形變值、超壓值曲線發(fā)現(xiàn)3次爆炸具有一致的曲線特征(見(jiàn)圖9)；通過(guò)對(duì)3次爆炸形變曲線擬合得到基于3次爆炸的爆炸形變模型：y=-0.160x3+0.273x2-1.854x+13.150。

圖9 次嚴(yán)重區(qū)域至基本不受影響區(qū)域形變—空氣沖擊波超壓模型Fig. 9 Deformation from severe area to basically unaffected area-air shock wave overpressure model

5 結(jié)論

1)通過(guò)獲取3次爆炸事件所引發(fā)的地表形變場(chǎng)，發(fā)現(xiàn)爆炸所引發(fā)的形變相較于天然地震具有不連續(xù)波動(dòng)的特點(diǎn)。

2)結(jié)合爆炸前后衛(wèi)星影像發(fā)現(xiàn)形變隆升區(qū)域?yàn)榻ㄖ芗瘏^(qū)，以此為依據(jù)可通過(guò)形變隆升的區(qū)域評(píng)估爆炸帶來(lái)的毀傷范圍。

3)對(duì)于較大沉降區(qū)域，可能處于不穩(wěn)定的地質(zhì)區(qū)，災(zāi)害后重建時(shí)需重點(diǎn)評(píng)估。

4)通過(guò)總結(jié)3次爆炸形變與空氣沖擊波超壓值曲線特點(diǎn)，構(gòu)建了次嚴(yán)重區(qū)域至基本不受影響區(qū)域形變—空氣沖擊波超壓模型。