基于曲波變換的航空電磁數據調平方法研究

高玲琦， 殷長春* ， 王寧， 劉云鶴， 蘇揚， 熊彬

1 吉林大學地球探測科學與技術學院， 長春 130026 2 桂林理工大學地球科學學院， 桂林 541006

0 引言

航空電磁法是基于直升機、固定翼飛機等移動搭載平臺的地球物理方法.它利用不接地回線進行電磁發射并通過磁傳感器接收電磁信號，實現對地下目體的探測.航空電磁具有效率高、成本低的特點，特別適合高山、沙漠、湖泊沼澤、森林覆蓋等地形地質條件復雜地區的資源探測.目前航空電磁法已在礦產資源、地下水和地熱、環境工程、城市地下空間、地質災害監測等領域發揮了重要的作用(Reed，1981；Wolfgram and Golden，2001；Haas，2004；Smith et al.，2006；Viezzoli et al.，2010；Wiederhold et al.，2010；Minsley et al.，2012).

由于采用飛機平臺搭載電磁設備，航空電磁數據中會產生很多地面電磁法遇不到的噪聲和誤差.最典型包括由飛機飛行狀態不同導致電磁系統發生變化，從而產生的調平誤差.調平誤差與隨機噪聲存在本質區別，調平誤差存在特定的物理成因，不具隨機性，因此需要采取特殊方法對其進行處理.Huang和Fraser(1999)將調平誤差分為測區調平誤差、測線調平誤差以及短波長調平誤差.測區調平誤差主要由不同架次之間電磁系統物理條件變化引起，在數據中表現為區塊狀的異常抬升.測線調平誤差是由不同測線間由于溫度、壓力、風速和風向等條件發生變化，導致儀器狀態發生變化而產生的(Huang，2008；Beiki et al.，2010；Siemon,2009).航空電磁(特別是直升機航空電磁)通常采用相鄰測線來回飛行的“S”形觀測方式，因此某一條測線上電磁系統面向太陽(或者迎風)，相鄰測線上則背對太陽(或者避風)，造成相鄰測線間系統狀態發生變化，進而在電磁信號中產生條帶狀調平誤差(稱為窗簾效應).另外，由儀器噪聲、飛行速度突然改變及局部熱變化等因素還可引起短波長調平誤差，在數據中表現為局部干擾(Huang，2008).

航空電磁數據調平旨在消除調平誤差，改善數據質量.傳統的位場數據調平大多基于位場隨高度變化緩慢的特點，采用切割線調平(Yarger et al.，1978；Bandy et al.，1990；Minty，1991；Luyendyk，1997；Ferraccioli et al.，1998).然而，航空電磁信號對飛行高度非常敏感，在切割線與測線交叉處，難以保證飛行高度一致，因此切割線調平無法進行電磁數據調平.Huang和Fraser(1999)提出了一種針對視電阻率數據的自動調平方法.該方法通過一維與二維數據非線性相關濾波器的組合實現數據調平.Fedi和Florio(2003)提出了一種基于小波變換的調平技術，然而該方法受到小波變換僅有三個方向特征的限制.Mauring和Kihle(2006)提出了基于移動微分中值濾波器的調平方法，該方法可以應用于非平行分布的測線數據調平.Huang(2008)提出了基于測線數據相關性的調平方法.該方法基于相鄰測線間信號具有較大的相關性，而調平誤差不具有相關性的特點，通過最小二乘擬合相鄰測線的數據差，實現航空電磁數據調平.Beiki等(2010)提出了基于一維和二維滑動窗口多項式擬合的航空電磁和航空磁測數據調平方法.Zhang等(2018)提出了基于主成分分析的調平方法并將其應用于航空電磁數據與航磁數據的調平.然而，上述方法均存在一定的局限性.其中，基于一維及二維濾波的調平方法本質上對數據整體進行了平滑處理，容易導致高頻細節的丟失，而基于誤差多項式擬合的方法由于其階數限制，常常不能完美地擬合誤差，造成調平過度或殘余.

為改善上述方法的局限性，本文提出一種基于曲波變換的調平方法，通過曲波變換的多尺度和多方向特征，調平誤差得以在特定尺度和方向內被精確擬合并去除，從而避免多項式在擬合誤差時的不準確性.另一方面，由于大多數尺度及方向未經調平處理，故數據畸變小，且不容易導致高頻細節的丟失.在曲波變換等多尺度分析方法中，閾值法通過對曲波系數設置閾值實現數據去噪和數據重構等.Donoho和Johnstone(1994)首先提出了用于小波變換重構的硬閾值方法，之后，Donoho(1995)又提出了一種用于噪聲去除的軟閾值方法.目前兩種方法均獲得廣泛應用.本文實現航空電磁數據調平，我們首先對飛行高度進行歸一化處理(Huang，2008)，在此基礎上對數據進行多尺度分析，通過設置閾值進行方向濾波有效去除數據中的調平誤差.這種新調平方法與基于曲波變換數據去噪的不同之處在于：數據去噪僅需對信號進行多尺度分析，通過對高階尺度設置閾值即可去除隨機噪聲(王寧等，2020)；然而，調平誤差具有方向性，但不具有隨機性，因此需要在特定方向和尺度上對其進行識別，進而通過設置閾值法去除調平誤差.下面我們首先介紹曲波變換及其多尺度和多方向性，進而通過對調平誤差進行分析，討論如何通過在不同尺度和不同方向上設置閾值對調平誤差進行去除.最后，我們通過對合成數據以及愛爾蘭Waterford地區頻率域航空電磁數據的調平驗證本文調平方法的有效性.

1 基于曲波變換的數據調平方法

1.1 曲波變換

曲波變換是在小波變換基礎上發展而來，同時繼承了小波變換的多尺度性質.除此之外，曲波變換還具有多方向性，這使得我們能夠利用曲波變換更好地表征二維或更高維的信號.本文將利用這些曲波變換特性實現航空電磁數據調平.為了方便理解，我們將首先介紹曲波變換的基本理論和離散算法.由于航空電磁數據是離散的，我們主要描述離散形式的曲波變換.

1.1.1 波數域中的笛卡爾窗函數

在二維波數域中，離散曲波變換通過笛卡爾窗函數平滑地提取不同尺度和不同方向上的數據信息(Candès et al.，2006).笛卡爾窗函數由徑向窗與角度窗一起共同構建.其中，徑向窗是一系列同心的正方形，可定義為

(1)

(2)

其中，ω=(ω1,ω2)為二維波數域中的變量，j為尺度參數，Φ由兩個一維低通濾波窗口的乘積定義，即：

Φj(ω1,ω2)=φ(2-jω1)φ(2-jω2)，

(3)

二維波數域是由所有徑向窗支撐而成，因此：

(4)

離散曲波變換中的角度窗Vj(ω)呈放射狀，可定義為

Vj(ω)=V(2?j/2_|ω2/ω1),

(5)

其中，|_j/2_|表示對j/2取整.在對徑向窗與角度窗分別進行定義后，離散曲波變換的笛卡爾窗函數可表示為

(6)

進而，我們引入等間隔的斜率：tanθl∶=l·2-?j/2_|,l=-2-?j/2_|,…,2?j/2_|-1，則(6)式可表示為

(7)

其中，剪切矩陣Sθ可表示為

(8)

圖1 離散曲波變換定義的楔形窗(參考Ying et al.，2005)Fig.1 Wedge-shaped window defined by discrete curvelet transform (Refer to Ying et al., 2005)

圖1展示了二維波數域中的窗函數.其中，最內層代表粗尺度.該尺度內的曲波不具有方向性(相當于小波).在粗尺度之外，徑向窗呈同心的正方形狀，而角度窗呈圖中虛線所示的放射狀.圖中灰色的楔形區域由徑向窗與角度窗共同支撐，它包含了特定尺度和方向上的數據信息.

1.1.2 空間域中的曲波函數

離散曲波函數為連續曲波函數的近似，它保留了連續曲波函數的所有性質(Ying et al.，2005).為了介紹離散曲波函數，我們首先定義連續曲波函數.在連續曲波變換理論中，母曲波φj(x)由楔形窗的逆傅里葉變換得到，進而所有尺度為2-j的曲波均可由母曲波經位移和旋轉得到(Candès et al.，2006).

(9)

其中，x=(x1,x2)為空間域自變量，旋轉矩陣定義為

(10)

由上面的討論，我們可以通過二維信號與曲波函數的內積計算曲波系數，即：

(11)

通過二維傅里葉變換，我們也可以在波數域中計算曲波系數，即：

(12)

其中，Uj為連續曲波變換中的窗函數，可由φj的傅里葉變換得到.

在定義了連續曲波函數后，我們可以利用與(9)式相似的形式定義離散曲波函數：

(13)

其中，b為離散變量，可表示為

(14)

圖2展示了四個不同的曲波函數與其對應的楔形窗.圖2a中，空間域內的曲波函數由不同顏色的橢圓標記.其中，綠色代表第2尺度內的曲波函數，藍色代表第3尺度的曲波函數，而紅色和黃色代表第4尺度內不同方向的曲波函數.圖2b中展示了波數域中對應的楔形窗.由圖可以看出，隨著尺度參數j增加，對應的曲波函數變小，其各向異性增強.值得注意的是，由于二維傅里葉變換存在90°的旋轉關系，空間域內的曲波函數與波數域內的窗函數在方向上互相正交.

圖3展示曲波函數在兩個互相垂直的方向上的分布特征.沿主軸方向，其波形從中心向兩側緩慢衰減，而在垂直于主軸方向，其波形從中心向兩側呈現快速震蕩衰減.此外，曲波函數波形的長度與寬度還滿足如下各向異性尺度關系：

(15)

曲波變換的各向異性特征使得其能更好地表示不同方向上的特征.這種特征非常適合于航空電磁數據中測線誤差調平.

1.1.3 基于Wrapping算法的離散曲波變換

Candès等(2006)提出了用于快速離散曲波變換的卷繞(Wrapping)和USFFT算法.本文中我們采用卷繞算法.(13)式中我們已經定義了離散曲波函數，我們采取與(12)式類似的形式在波數域中計算曲波系數，即：

(16)

式(16)經離散化可得：

(17)

必須指出，由于楔形區域無法直接進行傅里葉變換，利用式(17)進行離散曲波變換存在難度.為了解決此問題，我們在圖4給出“Wrapping”算法的實施步驟.圖4中高度為L1j、寬度為L2j黑色平行四邊形為波數域中楔形的支撐空間，該空間內包含的信息由灰色橢圓形區域表示.為進行二維傅里葉變換，我們定義一個圍繞坐標原點、邊長分別為L1j和L2j的正四邊形，并將包含在黑色平行四邊形中的數據信息Wrapping到圍繞原點的正四邊形中，從而可進行二維傅里葉變換.

圖2 空間域與波數域中的曲波(a) 空間域中四個不同的曲波函數； (b) 波數域中對應的楔形支撐空間(參考Herrmann and Hennenfent, 2008).Fig.2 Curvelets in spatial and wavenumber domains(a) Four different curvelets in spatial domain; (b) Wedge-shaped support in wavenumber domain(Refer to Herrmann and Hennenfent, 2008).

圖3 不同方向曲波函數波形(a) 沿主軸方向； (b) 垂直于主軸方向.Fig.3 Curvelet waveforms in varied directions(a) Along the major axis; (b) Perpendicular to the minor axis.

Wrapping算法可歸納為如下4個步驟：

圖4 曲波變換中的Wrapping算法(參考Candès et al.，2006)Fig.4 Wrapping algorithm in curvelet transform(Refer to Candès et al.，2006)

1.2 基于閾值法的電磁數據調平

本節我們首先解釋曲波變換可用于航空電磁數據調平的原因，然后討論具體調平方法.我們這里不考慮測區調平誤差，而將注意力集中于測線調平誤差和短波長調平誤差.這是由于前者可采用較為簡單的加/減背景場實現，而測線調平相對來說要困難的多.然而，曲波變換特別適合于對測線誤差進行調平.事實上，如前所述，測線調平誤差主要是由飛行方向改變導致的系統狀態變化引起的，它們在垂直于測線的方向上具有相對穩定的空間尺度.另外，調平誤差沿測線方向比較穩定(或者伴隨緩慢的非線性漂移).所有這些導致該誤差具有沿測線極強的方向性.因此，我們可以對信號進行多尺度和多方向分析，進而在曲波域特定尺度和角度方向上找出與調平誤差相關的信息，進行有用信號和調平誤差的有效分離，實現調平誤差去除.

圖5展示了理論模型的重構結果.圖5a給出理論數據，我們分別將振幅為最大異常的12.5%及振幅為0的調平誤差交替加入到圖5a相鄰測線的數據中，得到了圖5b中包含調平誤差的數據.考慮到理論模型共有60×60單元，我們在曲波域設置了4個尺度，分別在粗尺度上設置了小波系數，在第2尺度上設置了16個角度方向，在第3和第4尺度上設置了32個角度方向.進而，我們提取不同尺度和不同角度方向中的曲波系數，將信號重構回空間域，得到特定尺度和角度方向中的信號成分.從圖5給出的結果可以看出，隨著尺度參數j的增大，空間尺度逐漸變得精細.從圖5f可以看出，調平誤差集中在第4尺度，而從圖5g—k可以看出，調平誤差具有明顯的方向性.對比圖5h和圖5k可以看出，第4與第20個角度具有相同的方向，這是由于它們在波數域中的支撐楔處于同一方向上的兩個頂端.然而，由于這兩個支撐楔存在位置、形態和對稱性差異，它們包含的信息完全不同.最后，考慮到調平誤差僅存在于第4尺度的第4、5、20、21角度方向，我們僅對這些尺度和方向進行數據調平，得到的結果由圖5l給出.對比圖5l和圖5a可以看出，利用本文算法數據中的調平誤差被去除，信號得到有效恢復.

圖5 理論模型及不同尺度、不同角度方向的重構結果(a) 原始數據； (b) 加入調平誤差的數據； (c) 粗尺度重構結果； (d) 第2尺度重構結果； (e) 第3尺度重構結果； (f) 第4尺度重構結果； (g) 第4尺度第1角度重構結果； (h) 第4尺度第4角度重構結果； (i) 第4尺度第5角度重構結果； (j) 第4尺度第12角度重構結果； (k) 第4尺度第20角度重構結果； (l) 曲波變換調平結果.Fig.5 Theoretical model and reconstruction results from curvelet transform at different scales and angles(a) Original data；(b) Data with leveling errors added; (c) Coarse-scale reconstruction; (d) 2nd scale reconstruction result; (e) 3rd scale reconstruction; (f) 4th scale reconstruction; (g) Reconstruction of 1st angle from 4th scale; (h) Reconstruction of 4th angle from 4th scale; (i) Reconstruction of 5th angle from 4th scale; (j) Reconstruction of 12th angle from 4th scale; (k) Reconstruction of 20th angle from 4th scale; (l) Leveling by curvelet transform.

在確定了調平誤差聚集的尺度和角度后，我們可以利用閾值法消除這些誤差.本文中我們將硬閾值法、軟閾值法兩種方法應用于航空電磁數據調平.

硬閾值定義為

(18)

其中，x為曲波系數，T為閾值.硬閾值的主要思想是通過對低于閾值的系數進行置零來消除干擾.硬閾值法雖然易于實現，但其缺點是可能受局部畸變的干擾.為了使處理后的曲波系數更具連續性，我們定義如下形式的軟閾值：

(19)

由(19)式可知，軟閾值對所有系數都進行了修改.與硬閾值法相比，軟閾值法可以得到更加平滑的結果，但可能給數據造成較大的失真.

本文將結合兩種閾值方法的優點.為此，我們采取了2個步驟來消除調平誤差：

(1)我們通過硬閾值法x=hard(x,Ta),x∈Xja,la消除典型的條帶狀誤差；

(2)我們通過軟閾值法x=soft(x,Tb)，x∈X去除上一步驟的殘余誤差及其他短波長誤差.其中，集合X由所有的曲波系數組成，Ta、Tb是兩個步驟中相應的閾值.

上述處理程序可以給出更進一步的說明.事實上，第一步中硬閾值法及相應閾值Ta主要用于調平誤差集中的尺度和角度方向ja,la.由于調平誤差的幅值未知，我們在這一步驟需要設置較大的閾值Ta，并采用硬閾值法以避免信號失真.當ja,la中調平誤差幅值相對有用信號很大時，我們可以將所有ja,la內的曲波系數置零.然而，調平誤差是非常復雜的，可能存在方向性不明顯的局部擾動，或者存在第一步調平后的非線性殘余.在這種情況下，我們可以對所有的尺度和角度應用軟閾值法并設置閾值Tb，進一步去除前一步的調平殘差和短波長誤差.

在實際的數據處理中，首先通過曲波變換將數據轉換至曲波域，然后根據調平誤差的特征，選擇合適的方向閾值或全局閾值進行處理.由于調平誤差的復雜性，我們通過特定尺度和特定角度數據的重構判斷所選閾值的合理性.若重構結果中方向性明顯的誤差消失，說明所選閾值合理可行.最后，我們將數據變換到空間域，如此便可以消除數據中的調平誤差.

1.3 航空電磁數據調平步驟

航空電磁信號對飛行高度很敏感，因此在計算視電阻率或進行反演時，飛行高度往往作為變量.在對航空電磁數據進行調平時，我們需要首先校正飛行高度的影響.為此，我們采用Huang(2008)和Beiki等(2010)提出的方法，即利用收發距與飛行高度之比的立方對測量信號的實虛分量進行歸一化處理：

(20)

(21)

其中，I和Q分別為電磁信號實虛分量，m和n為受飛行高度影響較小的歸一化響應.Huang(2008)指出，當地下介質相對良導時，(20)和(21)式很準確；然而，對于高阻介質，為獲得準確結果，式(21)中的(h/s)的階數應該稍微降低.Beiki等(2010)利用這一近似將航空電磁信號校正到某一參考高度進行調平.本文中我們首先應用(20)和(21)式將數據進行歸一化，進而利用上述基于曲波變換的調平方法對歸一化的響應m和n進行調平，然后我們將數據轉換回到實虛分量.最后,通過已調平的數據計算視電阻率.如果視電阻率中沒有明顯的條帶狀誤差，則調平完成.由于曲波變換的多尺度和多方向特點，基于曲波變換的調平方法只對數據的某些細節尺度和角度進行局部改變，故對數據的原始相位關系影響很小，在去除調平誤差時不容易引起數據畸變.

2 合成數據調平

下面我們首先測試本方法在合成數據調平上的有效性.考慮到調平誤差可存在于各種航空數據中，本節僅給出不具有物理意義的合成模型.對于所有模型，我們在曲波域設置了5個尺度，角度數分別為第2尺度16個、第3、4尺度32個、第5尺度64個.圖6a—c給出了原始數據、加入調平誤差的數據以及調平誤差，圖6d給出了采用我們的曲波變換方法進行調平的結果，為對比起見，圖6e給出了Huang(2008)提出的基于測線相關性調平方法的結果.對比圖6d和圖6e，我們發現基于測線相關性的調平方法更容易造成數據失真，這主要是由于該方法難以利用多項式擬合復雜的非線性調平誤差.相比之下，曲波變換由于出色的二維信號表達能力，很好地還原了原始數據.從圖6f和圖6g展示的原始數據與調平后數據差中，我們也可以清楚地看到，本文方法的殘差比其他方法的殘差小.在點位(4000 m, 3000 m)附近，有一個延伸方向與測線相同的異常體.從理論上講，這類異常體對于調平而言是較為困難的.然而，由于曲波變換的多尺度和多方向特性，我們成功地將調平誤差從有效信號中分離出來，很好地恢復了該異常特征.

圖6 (a) 原始數據； (b) 加入調平誤差的數據； (c) 調平誤差； (d) 曲波變換調平結果； (e) 基于測線相關性的調平結果(Huang，2008)； (f) 圖(a)與(d)的差值； (g) 圖(a)與(e)的差值Fig.6 (a) Original data； (b) Data with leveling errors； (c) Leveling errors； (d) Leveling of curvelet transform; (e) leveling of line-to-line correlation by Huang (2008)； (f) Difference between (a) and (d)； (g) Difference between (a) and (e)

航空電磁中經常會出現測線方向與大地坐標系不一致的情況.此時，傳統的調平方法(例如基于二維濾波的方法)需要對坐標系進行旋轉.由于曲線變換的方向性，我們可以省去這些旋轉，僅通過選擇不同的角度參數很容易地解決問題.這意味著基于曲波變換的調平方法具有更大的靈活性.圖7給出一個測線與大地坐標軸斜交的模型和對應的曲波變換調平結果.由圖可以看出，調平誤差被很好地消除.

圖7 (a) 原始數據； (b) 加入調平誤差的數據； (c) 曲波變換調平結果Fig.7 (a) Original data； (b) Data with leveling errors； (c) Leveling of curvelet transform

航空電磁法中另一個常見的問題是測區形狀不規則.在這種情況下，數據調平不僅費時，而且調平結果也會受到影響.例如，基于測線相關性的調平方法必須對擬合多項式進行外插，而基于二維濾波的調平方法在不規則的邊緣處會產生一定的邊緣效應.然而，基于曲波變換的調平方法由于具有多尺度性，可以有效地處理不規則測區獲取的數據.事實上，在調平過程中，我們用一個常數來填充缺失的區域，然后進行調平.由于缺失區域具有一定的尺度，該區域的信息大多出現在曲波域的粗尺度上，對調平結果影響較小，因此在完成數據調平后我們僅需將該區域的數據剔除即可.圖8給出形狀不規則的合成模型的調平結果.由圖可以看出，盡管測區由于存在三角形缺失變得不規則，調平誤差仍然被很好地消除，數據得到有效地恢復.

圖8 (a) 原始數據； (b) 加入調平誤差的數據； (c) 曲波變換調平結果Fig.8 (a) Original data； (b) Data with leveling errors； (c) Leveling of curvelet transform

3 實測數據調平

為了進一步驗證本文方法的有效性，我們將其應用于愛爾蘭Waterford地區利用頻率域航空系統采集的電磁數據.該項目由愛爾蘭地質調查局(GSI)牽頭，由加拿大Sander地球物理公司與Unicorn礦業有限公司合作組織實施.圖9給出Waterford測區的位置及測線布設情況.其中，測線方向為N15°W，線距200 m，切割線方向為E15°N，線距2000 m.我們選擇圖9中紅色線框標注的面積約為11.5 km×13 km的區域進行測試.我們在曲波域中設置5個尺度，其中第2尺度包含32個角度,第3、4尺度包含64個角度，第5尺度包含128個角度.我們以頻率912 Hz為例對實測數據實虛部進行調平.

圖9 Waterford項目測區位置及測線布設(修改自Bates和Cho, 2016)Fig.9 Survey area and survey lines for the Waterford project(Modified from Bates and Cho, 2016)

圖10a和圖10d為實、虛部原始數據.我們首先利用(20)、(21)式將數據轉換為歸一化響應m和n，然后利用本文方法對m和n進行調平，再將數據轉換到調平后的實虛分量.圖10b和圖10e為調平后的數據.為了便于比較，我們在圖10c和圖10f中還展示了GSI采用的基于微分多項式擬合的調平結果(Beiki et al.，2010).由對比結果可見，對于虛部數據，基于曲波變換和基于GSI的多項式擬合方法均獲得很好的結果；對于實部數據，基于曲波變換的調平方法有效去除了調平誤差，未對數據產生明顯的畸變，然而，GSI利用的多項式擬合方法在消除調平誤差的同時，也去除了同方向上的異常(圖10c).

圖10 愛爾蘭Waterford測區實部(左側)與虛部(右側)數據調平結果(a)與(d)為原始數據； (b)與(e)為曲波變換調平結果； (c)與(f)為GSI的調平結果.Fig.10 Leveling results of in-phase (left column) and quadrature component (right column) for the data from Waterford project, Ireland(a) and (d) Original data； (b) and (e) Leveling by curvelet transform； (c) and (f) Leveling by GSI.

由于本文方法對實虛分量的處理是獨立進行的，因此對于調平結果來說，檢驗其視電阻率和相位變化至關重要.圖11給出圖10中調平前后數據計算的視電阻率和相位.由圖可以看出，由于實虛部數據中含有調平誤差，故利用原始數據計算視電阻率及相位時，誤差被引入其中，導致計算結果中存在明顯的條帶狀誤差.經曲波變換調平后，視電阻率數據中的條帶狀誤差明顯減少，相鄰測線間電性趨于相容、合理.作為對比，由GSI給出的調平結果計算得到的視電阻率在中間相對高阻區域仍存在明顯的條帶狀誤差，且造成了一定程度的數據畸變，特別是相位數據畸變更為嚴重.這主要是由于該方法在分別對實虛部進行調平時對原始數據畸變過多造成的.從圖11d、e、i、j可以得出相似的結論，即本文的調平方法對原始數據的視電阻率和相位畸變較小.

圖11 愛爾蘭Waterford測區調平前后數據計算得到的視電阻率(左側)與相位(右側)(a)與(f)由原始數據得到； (b)與(g)由曲波變換調平結果得到； (c)與(h)由GSI提供的調平結果得到； (d)與(i)分別為(a)與(b)、(f)與(g)的差； (e)與(j)分別為(a)與(c)、(f)與(h)的差.Fig.11 Apparent resistivity (left column) and phase (right column) before and after leveling for the data from Waterford project, Ireland(a) and (f) From original data; (b) and (g) From results by curvelet transform； (c) and (h) From results of GSI; (d) and (i) Differences between (a) and (b), (f) and (g)；(e) and (j) Differences between (a) and (c), (f) and (h).

4 結論

本文我們基于曲波變換成功地開發了一種航空電磁數據調平方法.該方法較之于傳統數據調平方法具有如下優勢：

(1)利用曲波變換多尺度和多方向的優點，可以方便地提取與飛行方向相關的調平誤差，同時有效保留了原始信號.

(2)對于實虛部數據的調平僅在曲波域的部分尺度和方向上進行，保留了實虛部數據中的主要成分，因此由其計算的視電阻率在調平誤差被極大削弱的同時，原始相位關系得到很好保證.

(3)對于傾斜測線及不規則測區，本文調平方法無需坐標旋轉和測區分割，同時避免了傳統方法產生的邊緣效應，具有很好的實用性.

理論和實測數據的實驗結果均表明，本文基于曲波變換的調平方法比傳統調平方法取得了更好的效果.由于調平誤差在航空地球物理數據中廣泛存在，而本文提出的方法實質上與數據類型無關，我們可以期待將該方法進一步推廣應用于其他航空地球物理數據處理中.

致謝感謝愛爾蘭地質調查局(GSI)與加拿大Sander地球物理有限公司允許本文使用他們的Waterford數據并提供技術支持.