MATLAB在深層構造預測中的應用
——以房山張坊水源地為例

孫少游，馬靜晨，劉全國，楊 楠，潘 懋

（1. 北京市地質工程勘察院，北京 100048； 2. 北京大學，北京 100087；3. University of British Columbia, Vancouver, Canada V6P4Y 7）

MATLAB在深層構造預測中的應用
——以房山張坊水源地為例

孫少游1，馬靜晨1，劉全國1，楊 楠2,3，潘 懋2

（1. 北京市地質工程勘察院，北京 100048； 2. 北京大學，北京 100087；3. University of British Columbia, Vancouver, Canada V6P4Y 7）

利用MATLAB強大的分析能力，基于北京房山巖溶區淺層節理裂隙數據的進行性半自動化的統計與分析，繪制淺層節理裂隙玫瑰花圖，并提取深層構造的走向線，最后與深層斷層、向斜與背斜等構造進行擬合。擬合結果顯示，北京房山巖溶地區淺層節理裂隙的統計走向與深層構造走向基本一致，北東向與北西向兩條主構造方向的擬合誤差分別為6.4%和2.1%。這一誤差遠小于設置的15%的誤差容限，表明此代碼在利用預測深層構造的應用中具有一定的正確性與潛力，能夠為野外數據的分析統計提供簡潔有效的幫助。但方法的有效性，仍需要通過更多的實例進行驗證。

MATLAB; 巖溶水; 深層構造預測; 房山

0 前言

北京地區水資源匱乏，其地下水資源在城市供水中占據十分重要的位置。巖溶水作為北京地下水的重要組成部分因賦存條件復雜、勘探難度大、開發利用成本高，勘察評價研究起步較晚（Pyne，1995）。根據《北京山前淺埋區巖溶水資源勘察評價》項目初步成果，北京市巖溶水天然補給資源量可達5.83×108m3/a（韓再生， 2003；Herman，2002）。北京地區碳酸鹽巖分布面積約4900km2，其中山區碳酸鹽巖面積約2900km2，占山區總面積的28%；平原隱伏區約2000km2，以奧陶系灰巖、寒武系灰巖和薊縣系白云巖為主（北京市水文地質工程地質大隊，2001(a)；北京市水文地質工程地質大隊，2001(b)；British Geology Survey, 2015；U.S.Department of the Interior et al，2008）。

北京市房山區張坊水源地是北京市建成的第一個大型應急供水工程，地處房山西南山前地區的霞云嶺—龍門臺復向斜水文地質單元，受下游斷裂構造控制，形成了相對獨立的巖溶水系統。區內巖溶裂隙發育、含水層厚度大，地下水易采易補。大石河從水源地穿過，在豐水期河水直接補給地下水，使該地區具備地表水和地下水聯合調蓄供水條件。北京市房山地區多年勘查已經積累了一定的基礎資料（Bear et al，1987），但由于巖溶區構造裂隙網絡空間分布的復雜性、巖溶水分布的不均一性以及巖溶區地下水與地表水之間相互轉化的復雜性，目前對北京房山地區仍缺乏綜合系統的研究，且尚未制定巖溶水資源合理開發保護方案、法規或條例。

Matlab是矩陣（Matrix）和實驗室（Laboratory）的組合，不僅自帶有強大的函數庫，實現各種運算，而且還允許用戶進行自編程，實現特定的功能。Matlab以其強大的數值運算與處理能力，在水資源規劃、水質評價、土壤重金屬污染空間分布（李鵬等，2014）等方面得到了廣泛的應用。本次利用MATLAB代碼實現的功能：包括提取三組深層節理數據；提取節理裂隙玫瑰花圖指示的淺層節理裂隙數據；擬合淺層與深層節理裂隙數據。本文研究的主要原因在于：程序化的方法能夠加快數據分析，尤其是大量數據存在的情況下每一個微小的錯誤都會導致結果的偏差，而重復分析十分耗費時間；其次，程序化的方法使分析過程更簡潔，分析結果的視覺化效果更直觀更完善；最后，本文的代碼可用于其他區域的檢測，最終用于深層構造的預測。這一代碼的實現將會加快分析和統計節理裂隙的進度并減小繪圖的工作量，同時能夠完善視覺效果，最終為科學有效地管理和保護巖溶水資源提供有效的幫助。

1 研究區域及背景

房山區位于北京市西南，其地理坐標為東經115°25′—116°15′，北緯39°30′—39°55′。西、南部與河北省接壤，北接門頭溝區，東接大興區，東北部與豐臺區相接。區內主要道路有房易路、周張路，京石高速公路從研究區的東部通過，成為房山區與北京城聯系的重要通道（Coppana et al， 1985；U.S. Geological Survey，2016）。

房山區地處北京西部山區和永定河沖積平原交匯處（Bouwer，2002），地貌主要為山地和平原兩種類型，山地面積廣大，位于房山區的西北部，約占總面積的三分之二，如圖1所示。房山區屬于大清河水系，境內地勢西北高，東南低。以鳳凰山—上寺嶺—霞云嶺分水嶺為界，北部為大石河，南部為拒馬河。北部大石河切割較淺，由霞云嶺向北經薊縣系地層進入寒武系、奧陶系灰巖地層之中，沿岸現巖溶溶蝕作用形成的洼地，并受構造作用影響形成佛子莊泉、河北泉等多處泉。南部低山丘陵區廣泛分布薊縣系霧迷山組碳酸鹽巖，受拒馬河切割影響，形成兩岸陡峭的峰叢、峰林等巖溶地貌，且巖溶洞穴較發育，形成石花洞、云水洞等重點景觀巖溶洞穴。山區巖溶裂隙發育，地下水入滲補給條件好。

圖1 北京市房山區地貌圖Fig.1 Terrain of Fangshan District

2 方法與材料

2012年4月—2013年4月，野外調查小組通過對北京市房山區地表各時期地層野外出露的節理進行詳細的觀察與測量，共采集來自295個觀測點的2095條節理裂隙數據，其中忽略延伸小于1m的節理，共計有效節理2036條，裂隙數據采樣點分布如圖2所示。

圖2 房山區節理裂隙數據分布點Fig.2 Fangshan Distribution of data collected points

3 結果與分析

將2036條有效的節理裂隙數據根據所屬地層劃分為32個組，并根據其走向數據在MATLAB中繪制玫瑰花圖，這些玫瑰花圖指示的總體方向如圖3所示。但基于所有的走向分布，我們可以利用MATLAB代碼提取出主要的走向，總體來說，北京房山地區以低應變脆性破裂變形為主，形成X形共軛節理。節理面平直、傾角較陡、延伸中等。利用獲得的節理數據繪制節理走向玫瑰花圖可見，房山區節理發育指示出的優勢方位如圖4所示。

在獲取了淺層節理裂隙所指示的主要走向之后，利用MATLAB去擬合淺層與深層構造的基本走向。首先導入基本構造圖（圖5）并提取基本構造圖中所顯示的深層構造走向。隨后，MATLAB代碼將進行淺層節理裂隙與深層構造走向的擬合，擬合的基本準則為：在15%的誤差范圍之內尋找與淺層節理裂隙提取的主要走向吻合的深層構造。結果顯示，兩條淺層節理裂隙的主方向與深層構造的擬合誤差分別為6.4%和2.1%。

圖3 玫瑰花圖指示的走向分布Fig.3 All strike directions in rose diagram

圖4 兩條主要的走向Fig.4 Two major strike directions

圖5 房山區節理裂隙與區域斷裂構造關系示意圖Fig.5 The base map for deep structure

較小范圍的誤差，表明淺層節理裂隙所反映的構造走向與深層構造（如斷層、背斜和向斜）的實際走向基本吻合。因為區域地質構造是多期次不同構造運動的結果，而這些大型的深層構造均為地質歷史時期多起構造應力場作用的結果，因此會疊加反映多期構造運動的應力場方向，但是簡而言之，無論在淺層節理裂隙還是深層構造中，區域主應力與構造走向垂直。

構造應力場，就是在一個空間范圍內構造應力的分布。一個地區的構造應力可以用主應力方向來表示，即剪應力等于零的方向。若要能夠產生構造變形，一般需要多向應力，其中3個主應力方向上的應力值一般應該是不相等的，其中必有一個為最大主壓應力方向（σ1），一個為中間主應力方向，一個為最小主應力方向。當3個主應力值相等時，就處于靜巖壓力狀態，此時則不會發生構造變形。

4 分析與后續研究

MATLAB具有強大的數值運算與處理能力，本研究雖然整個分析是在MATLAB中開發的，但整個過程并非真正的自動化，在某些步驟仍需要人機交互進行。同時，獲得主要構造走向，15%誤差容限和獲取玫瑰圖過程的平均方法都對結果具有重要影響。基于上述信息，來自該文章的結果需要通過更多應用程序和實踐進行驗證：

（1）誤差的疊加可能大于我們的預期，因為誤差產生于過程的多個節點；

（2）方法的有效性需要通過更多的應用和實踐來驗證；

（3）MATLAB中的編碼方法可能不是最簡潔的；

（4）鑒于作者在MATLAB上的知識有限，MATLAB編碼可能不如我們預期的有效；

（5）為MATLAB在深層構造預測中的應用提供了范例，也為后續的研究進行了有益的嘗試。

Bear J, Berkowitz B, 1987. Groundwater flow and pollution in fractured rock aquifers[J]. Developments In hydraulic Engineering, (4)： 175-238.

Bouwer H, 2002. Artificial recharge of groundwater：Hydrogeology and engineering[J]. Hydrogeology Journal，10(1)：121-142.

British Geology Survey. Elucidating the hydrogeological issues associated with aquifer storage and recovery in the UK[EB/OL].(2015-11-23)[2016-10-22]. http：//nora.nerc.ac.uk/512266/.

Campana ME，Mahin DA, 1985. Model-derived estimates of groundwater mean ages，recharge rates，effective porosities and storage in a limestone aquifer[J]．Journal of Hydrology，76(3-4)：247-264．

Herman Bouwer, 2002. Artificial recharge of groundwater：hydrogeology and engineering[J]. Hydrogeology Journal, 10(1)： 121-142.

Pyne RDG,1995. Groundwater recharge and wells： A guide to Aquifer Storage Recovery[M]. Florida：Lewis Publishing.

U.S.Department of the Interior,U.S.Geological Survey. Ground-water recharge in the arid and semiarid southwestern United States∥USGS Professional Paper 1703： 375[EB/OL]. (2008-8-8)[2016-10-22]. https：//pubs.usgs.gov/pp/pp1703/.

U.S. Geological Survey.U.S. Geological Survey artificial recharge workshop proceedings,[USGS Open-File Report 02-89] [EB/OL] (2016-11-28)[2016-12-5]. https：//water.usgs.gov/ogw/pubs/ofr0289/.

北京市水文地質工程地質大隊，2001. 北京地區巖溶水資源研究與開發利用工程可行性研究報告[R].

北京市水文地質工程地質大隊，2001. 北京地區巖溶水資源研究與開發利用工程工作總體設計[R].

韓再生，2003. 為可持續利用而管理含水層補給--第四屆國際地下水人工補給會議綜述[J]. 地質通報，22(2)：142-143.

李鵬，陳永當，葉宏武，等． 基于Matlab 的城市表層土壤重金屬污染空間分布的研究[J]．環境保護科學，2014，40(2)：51-54.

Application of MATLAB in Deep Structure Prediction——Application in Fangshan, Beijing

SUN Shaoyou1, YANG Nan2,3, PAN Mao2

(1. Beijing Institute of Geological & Prospecting Engineering, Beijing100048; 2. Peking University, Beijing 100871; 3. University of British Columbia, Vancouver, Canada, V6P 4Y7)

Taking advantage of the robust analytic ability of MATLAB, the analysis of shallow fractured joints in Fangshan karst area of Beijing were analyzed by a semi-automatic coding to draw the shallow joints, and extracted the deep structure trending lines. The results showed that the statistical trending of shallow joints in the Fangshan karst area of Beijing is basically the same as that of the deep tectonics, and the fi tting errors of the two main tectonic directions are 6.4% and 2.1%, respectively. This error is far less than the 15% error margin of the set, indicating that the code has some correctness and potential in the application of predicting deep structure, which can provide a simple and effective help for the analysis and statistics of the fi eld data, but the method’s effectiveness still need to be verif i ed through more examples.

MATLAB; Karst Water; Deep Structure Prediction; Fangshan District

P542.3

A

1007-1903（2017）01-0117-04

10.3969/j.issn.1007-1903.2017.01.021

孫少游（1983- ），男，工程師，主要研究方向：地下水資源。E-mail：418171619@qq.com