基于Python和Fortran程序?qū)崿F(xiàn)瞬變電磁法一維反演軟件的設(shè)計(jì)與開(kāi)發(fā)

孫大利 李方舟 費(fèi)靖淇

中國(guó)地震局第一監(jiān)測(cè)中心 天津 300171

引言

近年來(lái)隨著城市化的不斷發(fā)展，城市規(guī)模的不斷擴(kuò)大，地上有限的空間資源難以滿足人們?nèi)粘Ｉ?、工作的需求。與之相對(duì)的地下空間，大部分城市的開(kāi)發(fā)程度還很低，具有廣闊的發(fā)展前景[1]。地球物理勘探可以有效探測(cè)到地下異常信號(hào)，通過(guò)對(duì)異常信號(hào)的反演得到地下異常結(jié)構(gòu)分布，用以輔助地質(zhì)勘查解釋，為城市地下空間的科學(xué)合理的利用提供保障。

在城市地下空間開(kāi)發(fā)過(guò)程中，地下水的分布對(duì)工程的進(jìn)展有著極大的影響，瞬變電磁法是一種地球物理勘探方法，應(yīng)用于礦產(chǎn)勘探和地下水調(diào)查等，對(duì)于地下水分布具有良好的反映，但是瞬變電磁法勘探原理較為復(fù)雜，數(shù)據(jù)處理更是專業(yè)性強(qiáng)、難度大，需要專業(yè)的技術(shù)人員進(jìn)行分析處理，難以普及。針對(duì)這一問(wèn)題，亟須設(shè)計(jì)和開(kāi)發(fā)一款瞬變電磁法反演軟件，實(shí)現(xiàn)數(shù)據(jù)處理流程的自動(dòng)化，降低應(yīng)用門檻。目前，主流的編程語(yǔ)言Python具有豐富的標(biāo)準(zhǔn)庫(kù)，擴(kuò)展性和可嵌入性強(qiáng)，可以提供簡(jiǎn)潔的交互界面，用以滿足普通用戶的使用，但是對(duì)于復(fù)雜的數(shù)值計(jì)算執(zhí)行效率低，占用內(nèi)存較大，耗時(shí)長(zhǎng)。而Fortran語(yǔ)言語(yǔ)法格式接近數(shù)學(xué)公式的自然描述，執(zhí)行效率高，在數(shù)值計(jì)算領(lǐng)域運(yùn)用廣泛，可以滿足瞬變電磁反演算法的復(fù)雜計(jì)算。

綜上所述，本文基于Python和Fortran語(yǔ)言聯(lián)合設(shè)計(jì)開(kāi)發(fā)，核心算法采用瞬變電磁法中心回線系統(tǒng)Occam一維反演算法，使用Fortran語(yǔ)言編寫，外部交互界面使用Python語(yǔ)言編寫，最后在將兩部分銜接形成完整的反演軟件，實(shí)現(xiàn)自動(dòng)化處理功能。

1 算法介紹

瞬變電磁法反演算法的核心是構(gòu)建目標(biāo)函數(shù)，選取最優(yōu)算法，通過(guò)調(diào)整搜索策略，使得實(shí)際探測(cè)數(shù)據(jù)與理論數(shù)據(jù)具有較好的擬合度，來(lái)到達(dá)探測(cè)地下電性結(jié)構(gòu)分布特征的目的。Occam反演算法是1987年StevenC.Cconstable首次提出的，由于反演過(guò)程是非唯一不確定的，具有多解性，在求解策略上，Occam反演算法通過(guò)引入粗糙度，選擇模型盡可能簡(jiǎn)單平滑的解，從而減少對(duì)實(shí)測(cè)數(shù)據(jù)的過(guò)度解釋，并且消除簡(jiǎn)單分成模型中的任意不連續(xù)性[2]，在保證算法穩(wěn)定的同時(shí)，加快收斂速度，提高處理效率。本文模擬的數(shù)據(jù)采集裝置是中心回線系統(tǒng)，通過(guò)移動(dòng)線圈進(jìn)行工作，適合城市地下空間探測(cè)。

2 軟件設(shè)計(jì)開(kāi)發(fā)

程序主要分為3大部分，包括主程序、Occam反演算法程序和操作程序。

主程序：采用Pyside2內(nèi)置可執(zhí)行程序QtDesigner設(shè)計(jì)界面組件和布局，生成ui文件，然后在主程序中調(diào)用PySide2模塊中的QtWidgets、QtUiTools、QtCore用于導(dǎo)入ui文件，生成界面窗口、加載界面組件，并初始化事件觸發(fā)信號(hào)和按鈕。

Occam反演算法程序：采用Fortran語(yǔ)言編寫，使用visual studio編譯器生成可執(zhí)行文件保存到軟件所在目錄下，等待調(diào)用。

操作程序：包括輸入部分、運(yùn)行Fortran程序和輸出部分，輸入部分調(diào)用shutil模塊和QUiLoader函數(shù)，用于導(dǎo)入實(shí)測(cè)數(shù)據(jù)文件，以及將界面內(nèi)用戶輸入信息傳入相應(yīng)變量中。運(yùn)行Fortran程序部分調(diào)用subprocess和threading模塊，創(chuàng)造子線程調(diào)用外部可執(zhí)行程序。最后輸出部分調(diào)用numpy、matplotlib模塊，將計(jì)算后輸出數(shù)據(jù)數(shù)組化并繪制折線圖[3]。

2.1 基于Fortran語(yǔ)言實(shí)現(xiàn)Occam反演算法

電磁法勘探一般實(shí)測(cè)數(shù)據(jù)和地下模型之間的關(guān)系可以表示為：

反演問(wèn)題是已知實(shí)測(cè)數(shù)據(jù)d反過(guò)來(lái)求解真實(shí)模型空間m，根據(jù)Occam反演算法我們需要在求解式中加入正則化項(xiàng)，從而在眾多滿足數(shù)據(jù)擬合條件的解中選取最簡(jiǎn)單光滑的解模型?；谏鲜鲈?，可以得到目標(biāo)函數(shù)，同時(shí)定義：

式中，d表示實(shí)測(cè)數(shù)據(jù)向量，m表示模型參數(shù)空間向量，X為期望擬合差，W為對(duì)角加權(quán)矩陣為拉格朗日乘子，為粗糙度矩陣。

求取目標(biāo)函數(shù)的極小值，即可直接獲得此次迭代的新模型：

2.2 基于Pyside2模塊開(kāi)發(fā)交互界面

軟件交互界面主要框架采用由諾基亞公司開(kāi)發(fā)的Python圖形化界面開(kāi)源庫(kù)Pyside，所有操作均通過(guò)用戶界面上創(chuàng)建的相應(yīng)輸入?yún)^(qū)和按鈕完成，主要包括參數(shù)設(shè)計(jì)區(qū)、運(yùn)行日志區(qū)和結(jié)果輸出區(qū)。

其中參數(shù)設(shè)計(jì)設(shè)置“選擇文件”按鈕，用于快速導(dǎo)入實(shí)際采集到的數(shù)據(jù)，設(shè)置“修改參數(shù)”按鈕，驅(qū)動(dòng)內(nèi)部程序根據(jù)用戶在文本框內(nèi)輸入的內(nèi)容，修改模型初始化參數(shù)；運(yùn)行日志區(qū)會(huì)根據(jù)用戶點(diǎn)擊操作實(shí)時(shí)反饋內(nèi)部程序運(yùn)行信息，提示進(jìn)程并反饋錯(cuò)誤操作；結(jié)果輸出區(qū)會(huì)根據(jù)算法程序得到的反演記過(guò)繪制折線圖，直觀反映地下電性分布結(jié)構(gòu)，同時(shí)所在路徑下會(huì)生成結(jié)果文件，輸出并保存計(jì)算結(jié)果。圖1展示運(yùn)行初始化和模型參數(shù)加載成功界面。

圖1 界面展示

圖2 反演結(jié)果示例

2.3 實(shí)例演示

使用瞬變電磁中心回線裝置探測(cè)某區(qū)域，裝置直徑10m，電流1A，得到單點(diǎn)的實(shí)測(cè)數(shù)據(jù)。使用一維Occam反演軟件進(jìn)行處理，如圖4所示，在界面上生成折線圖，并將反演數(shù)據(jù)保存。可以看到經(jīng)過(guò)軟件處理后的結(jié)果呈現(xiàn)出多層狀的分布，代表著地下不同電阻率地層的分布情況，電阻率低的分層是含有地下水的可能性最高，可以作為地下空間開(kāi)發(fā)時(shí)的輔助意見(jiàn)。

3 結(jié)束語(yǔ)

本文設(shè)計(jì)的軟件，通過(guò)簡(jiǎn)潔的輸入模塊和直觀的輸出模塊，將復(fù)雜的參數(shù)調(diào)整、核心計(jì)算和數(shù)據(jù)可視化等操作封裝在內(nèi)部模塊中，用戶只需要進(jìn)行簡(jiǎn)單的操作，就可以得到瞬變電磁法探測(cè)數(shù)據(jù)的處理結(jié)果，省去了中間專業(yè)化的處理，極大地簡(jiǎn)化了實(shí)際操作流程，可以更快更高效的將該方法運(yùn)用到實(shí)際地下空間開(kāi)發(fā)中。通過(guò)對(duì)實(shí)際案例的應(yīng)用和測(cè)試，取得了很好的實(shí)踐效果。