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