CN105911603A - 基于天然电场的四维物探方法 - Google Patents

Abstract

本发明公开一种基于天然电场的四维物探方法，它是在现场数据采集前对探测工程进行完备的设计，通过探测仪器在地面同步采集每条测线上多根探针接收的天然电磁波穿透地层结构后在地表形成的分布电场信号并实时进行频谱细化、特征信息识别后保存，上位机读取探测数据经综合统计后绘制测线四维物探剖面图并对存在特征信息的测线进行定性、定量分析，然后综合多条测线的探测分析结果得出探测工程结论。其有益效果是：在有效运用频域三维物探方法的基础上引入随时间变化的特征参数而形成的四维物探方法，与现有物探方法相比：探测仪器所需存储空间少、同步探测测点数据的关联性好且抗干扰能力强，能对探测目标地质体进行定性、定量的综合分析。

Description

基于天然电场的四维物探方法

技术领域

本发明属于地球物理探测技术，具体地说是一种基于天然电场的四维物探方法。

背景技术

天然电场法是大地电磁测深法（Magnetotelluric Sounding，简称MT）中的一种，MT法是苏联学者Tikhonov(1950)和法国学者Cagniard(1953)50年代初提出来的利用天然交变电磁场研究地球电性结构的一种地球物理勘探方法；由于它不用人工供电、成本低、工作方便、不受高阻层的屏蔽、对低阻层分辨率高，而且勘探深度随电磁场的频率而异，浅可以几十米，深可达数百公里，因此，近年来在许多领域都得到了成功的应用，引起了地球物理学家的广泛兴趣和极大的重视。天然电场法利用雷电场和游散电流场作为场源，以岩矿石电阻率差异为基础，在地面沿特定的剖面测量大地电场的水平分量，并通过对其变化规律的研究来解决水文、工程地质和找矿勘探等地质任务；由于不需要人工供电，可省去笨重的供电设备，因此具有简单、轻便、经济、效率高等特点。天然电场探测方法类似于音频大地电磁法，所不同的是野外工作中只测量单方向的大地电场的水平分量，在场源上除了利用天然电磁场变化在大地中感应的大地电流之外，还利用了工业游散电流场。随时间变化的天然电磁场在大地中感应出电流，这种地电流频谱宽、分布广、能量丰富，是一种能解决地质任务的天然场源。天然电场产生的原因很复杂，通常认为频率大于1Hz的天然电场主要是由赤道地区雷电放电产生的，局部地区是由工业游散电流及谐波引起的。大地电磁场的振幅、频率、方向均随时间而随机变化，但这种随机有一定的规律性；如具有同源性，这表现在某一瞬间，在几百平方公里范围内，场的振幅、频率和强度均能同步变化。天然电场的场源是由雷电和游散电流场组成的，由远处传输来的雷电电磁场，可近似看成一个平面电磁波，它的传输特性满足麦克斯威方程组；至于游散电流场，当它刚入地时电流线呈放射状，而远离接地点时，可近似认为电流线互相平行，可视为似稳的均匀电场，它在物探方面的应用原理与直流电法中的中间梯度相类似。

现有的天然电场探测中，一般是“探针+电缆+仪器”探测方式，探测仪器通过探测电缆与探针相连以接收探针获得的地表每个测点的电磁场信息，每个测点信息的获取至少需要2根探针，通常称探针与探针之间的距离为“极距”、测点与测点之间的距离为“点距”，极距大小取决于勘查纵深范围要求，而点距大小决定了勘查的平面精度。现有方法的实际探测工程中，为得到满足要求的极距和点距，需频繁移动探针或布置大量的探针和探测电缆（或探针与仪器的连接线），存在以下不足：（1）效率低 对于仅配2根探针的仪器，采集一个测点的数据只需几秒到几分钟的时间（因探测深度不同而不同），而移动探针并布置好一个测点往往需要几十分钟，对于配置多探针的仪器，每探针在满足极距和点距要求的同时还需布置一根与仪器的连接线，布线繁琐且工程量大；（2）探测数据的时空差异大 现有探测方式需逐个测点进行数据采集，测点之间采集数据的起始时刻相隔较长，造成各测点数据取样时的时空条件不一，增加了后续数据分析与解释的困难；（3）数据解析能力差 现有的探测方式中，大都采用基于时域数据分析的二维或三维方法，难以解决较为复杂地质工程问题。现有天然电场探测方法中，对地下地质体的频域特征研究较少，能否发明一种根据地下地质体频域特征来分析地层结构的物探方法，目前国内外科技界一直没有解决，尚无相关的研究成果和产品报道。

发明内容

针对目前物探方法存在的不足，本发明为解决该问题提供了一种基于天然电场的四维物探方法。

本发明的术方案是：在频域三维物探基础上引入随时间变化的特征参数，通过探测仪器在地面同步采集每条测线上多根探针接收的天然电磁波穿透地层结构后在地表形成的分布电场信号并实时进行频谱细化、特征信息识别后保存，由上位机读取探测数据经综合统计后绘制测线四维物探剖面图并对存在特征信息的测线进行定性、定量分析，然后综合多条测线的探测分析结果得出探测工程结论；分以下四个步骤进行：（1）按照探测工程的具体要求，设计能满足探测要求的K条测线以及每条测线的测深个数J、重复探测次数I和相应的探测参数；（2）通过探测仪器同步采集每条测线上多个测点的天然电场信息，每个测深重复采集I次并实时进行频谱细化、特征信息识别后保存，直到采集完第J个测深后再进行下一条测线探测；（3）上位机读取探测装置保存的测线探测数据，综合统计测线上每个测深的每个测点的I次重复探测数据及其动态信息、裂隙信息的识别结果，绘制测点信号谱图和测线四维物探剖面图并进行测线低值异常判断，若测线存在动态信息、裂隙信息、低值异常则相应地进行地下径流、地下空区、断裂构造的定性、定量分析，否则进行下一条测线分析，直到完成第K条测线分析；（4）综合K条测线探测分析结果得出探测工程结论，对探测目标地质体进行综合评估。

在本发明中，测深h与探测仪器工作频率关系为f _c＝ (1591.58/h)²；所述的探测参数包括采样频率f _s=k₁+k₂*(f _c-k₃)、分析半带宽B=k₄+k₅* (f _c-k₃)、细化倍数D=f _s/2B、滤波器的半阶数M=4D、采样数据个数N=2M+511D；当1Hz <f _c≤10Hz时，取系数k₁=12.8、k₂=5.6889、k₃=1.0、k₄=0.6、k₅=0.26667；当10Hz <f _c≤100Hz时，取系数k₁=64.0、k₂=3.5556、k₃=10.0、k₄=3.0、k₅=0.03333；当100Hz <f _c≤500Hz时，取系数k₁=384.0、k₂=2.24、k₃=100.0、k₄=6.0、k₅=0.0175；当500Hz <f _c≤1200Hz时，取系数k₁=1280.0、k₂=1.77、k₃=500.0、k₄=13.0、k₅=0.01143；当1200Hz <f _c<6000Hz时，取系数k₁=3648.0、k₂=1.95、k₃=1200.0、k₄=16.0、k₅=0.002917。

在本发明中，探测仪器由16根探针、2段信号电缆、1个探测装置组成；16根插入地表的探针用于接收天然瞬变电磁波穿透地层在地面形成的分布电场信号，分别经2段信号电缆（每段电缆软件8根探针、探针间距大小可探测任务要求选择）连接到探测装置；探测装置采用以DSP+CPLD为核心的结构设计，具有信号调理、16通道信号同步采集、ZFFT频谱细化、特征信息识别和数据保存等功能。

在本发明中，特征信息识别方法是：针对探测仪器接收天然电场信号并细化得到分析带宽f _c-B～f _c+B的分辨率为∆f的256条谱线，在频谱中搜索出最大幅值F_MAX对应的谱线L_MAX及另外12个次大幅值谱线L_MAX1、L_MAX2、…、L_MAX12并同时计算频谱的平均幅值F_J、最大幅值谱线L_MAX及其前后各3条谱线一起的平均幅值F_J0；接着缓存256条谱线到另一指定存储区，在其中按从前往后两两比较幅值低者置0、从后往前两两比较幅值低者置0的原则找出频谱中的所有极值点，再将幅值小于F_MAX/30的极值点置0，结合原频谱计算次大极值点谱线J_MAX1、J_MAX2及其前后各3条谱线一起的平均幅值F_J1、F_J2；当极点数＜3进行裂隙信息识别，当3≤极点数＜30进行动态信息识别，极点数≥30无特征信息退出识别，此时的主频频率为仪器工作频率f _c，幅值为谱线L_MAX、L_MAX1、L_MAX2的平均幅值。

所述的裂隙信息识别为：当存在频谱能量衰落且极点谱线不分散时，若最大幅值谱线L_MAX频率非工频或其高次谐波、也不是测点区域已知的固有频率成分，则该频谱信息为裂隙信息，主频信息为最大幅值谱线L_MAX的幅值、频率，反之则无特征信息退出识别；频谱能量衰落判据是：若谱线L_MAX的幅值F_MAX≥20F_J并谱线L_MAX12的幅值F_MAX12≥F_J且12个次大幅值谱线L_MAX1、…、L_MAX12分布在最大幅值谱线L_MAX前后各6条谱线范围内，则频谱能量存在衰落现象；极点谱线分散判据是：若在最大幅值谱线L_MAX前后各6条谱线范围内存在2个以上极值点，则极点谱线分散。

所述的动态信息识别为：在满足F_J0≥10F_J则频谱能量集中，判断最大幅值谱线L_MAX的幅值F_MAX0≥3F_J0则L_MAX离散、次大极值点谱线J_MAX1的幅值F_MAX1≥3F_J1则J_MAX1离散、次大极值点谱线J_MAX2的幅值F_MAX2≥3F_J2则J_MAX2离散，通过边频识别规则一、规则二、规则三判断是否存在关于L_MAX、J_MAX1及J_MAX2的对称极值谱线来确定地下径流动态信息的存在，否则无特征信息退出识别；边频识别规则一是：L_MAX离散且极值点分布于其左、右两边，若存在关于L_MAX对称的边频极值谱线，则主频为L_MAX，动态信息为一阶边频的信息，若不存在关于L_MAX对称的边频极值谱线、但存在关于J_MAX1或J_MAX2对称的边频极值谱线且相应的J_MAX1或J_MAX2离散，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，否则无特征信息退出识别；边频识别规则二是：L_MAX离散但主要极值点分布在其一边，若J_MAX1或J_MAX2离散且存在关于J_MAX1或J_MAX2对称的边频极值谱线，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，否则无特征信息退出识别；边频识别规则三是：L_MAX不离散，若J_MAX1或J_MAX2离散且存在关于J_MAX1或J_MAX2对称的边频极值谱线，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，若J_MAX1或J_MAX2均不离散但存在包括L_MAX、J_MAX1、J_MAX2中至少2个在内的4个以上的点距误差在±2条谱线内的极值谱线，则主频频率为L_MAX谱线频率、幅值为L_MAX与J_MAX1的均值，动态信息的频率为一阶边频频率、幅值为J_MAX1与J_MAX2的均值，否则无特征信息退出识别。

在本发明中，低值异常是指四维物探曲线中部分曲线在一定连续测深范围内或全部曲线在整个测深范围内出现部分测点的探测电性值比测线上其它测点的探测电性值低2个及2个以上数量级。

在本发明中，定性分析是根据多测点主频f _c幅值F_MAX的变化规律来评估地下空区的塌陷、充填情况；通过多测线的四维物探曲线低值异常来确定断裂构造走向并结合动态信息、裂隙信息参数来综合分析断裂构造的存在状态。

在本发明中，定量分析是根据主频频率f _c、相位θ _MAX (f)、幅值F_MAX，动态信息的边频频率f _b、边频综合幅值e (t)和裂隙信息的幅值小于频谱平均幅值F_J的第一条衰落谱线的相位θ _MIN (f)、幅值F_MIN以及传输路径公共损耗K来确定地下径流、地下空区、断裂构造的埋藏深度h=1591.58(1/f _c)^0.5，

地下径流或导水断层的水量，

地下空区的横向尺寸B_x= 4.775*10⁷*[θ _MAX (f)-θ _MIN (f)] / f _c，

地下空区的纵向尺寸B_x=10^{[lg (FMAX/FMIN)+K/20-1.172-lg fc]}。

在本发明中，目标地质体综合评估是对四维物探剖面图进行分解，综合运用低值异常、动态信息、裂隙信息的分析方法对其做出的定性、定量分析。

本发明的有益效果是：在有效运用频域三维物探方法的基础上引入随时间变化的特征参数而形成一种基于天然电场的四维物探方法，与现有物探方法相比：探测仪器所需存储空间少、同步探测测点数据的关联性好且抗干扰能力强，能对探测目标地质体进行定性、定量的综合分析；能满足不同探测目的的物探工程需要。

附图说明

图1是本发明的四维物探方法流程图；

图2是本发明的实施例工程设计测线分布图；

图3是本发明的探测仪器现场测线布置示意图；

图4是本发明的特征信息识别流程图；

图5是本发明的动态信息特征图；

图6是本发明的裂隙信息特征图；

图7是本发明的实施例测线信号谱图；

图8是本发明的实施例测线数据综合统计表；

图9是本发明实施例测线一的四维物探剖面图；

图10是本发明实施例测线二的四维物探剖面图；

图11是本发明实施例测线三的四维物探剖面图；

图12是本发明实施例测线四的四维物探剖面图；

图13是图12的横向分解图；

图14是图12的纵向分解图；

图15是图12的单层次分解图。

具体实施方式

下面结合附图和实施例对本发明的方案作进一步详细说明。

参见附图1，一种基于天然电场的四维物探方法，它是在频域三维物探基础上引入随时间变化的特征参数，通过探测仪器在地面同步采集每条测线上多根探针接收的天然电磁波穿透地层结构后在地表形成的分布电场信号并实时进行频谱细化、特征信息识别后保存，由上位机读取探测数据经综合统计后绘制测线四维物探剖面图并对存在特征信息的测线进行定性、定量分析，然后综合多条测线的探测分析结果得出探测工程结论；分以下四个步骤进行：

步骤一 按照探测工程的具体要求，设计能满足探测要求的K条测线以及每条测线的测深个数J、重复探测次数I和相应的探测参数。

大地电磁法（MT）是深部地球物理探测必不可少的手段，已在国内外的矿产普查和勘探、地热田调查、地壳和上地幔电性结构研究、海洋地球物理、环境地球物理和地质工程中广为应用，当前的MT法大都是基于时域数据分析的二维或三维方法，难以解决较为复杂地质问题。本发明方法是在自有的频域三维物探方法基础上引入了随时间变化的特征参数而形成的基于天然电场的四维物探方法，分探测工程设计、现场数据采集、测线数据分析、综合结论等四步进行，具体过程如下：参见附图2，为考察某矿区以及其周边含水层和断裂构造的分布情况和具体位置，根据该区域的具体地形和踏勘情况，设计探测工程的测线布置如图2所示。计划在地热异常点周边安排A、B、C、等3组测线，测线水平间距约60m左右（具体根据实际地形决定），各测线的功能及作用如下：A组测线，西北→东南（每条测线平均长度约900m），主要探测该地区西南→东北方向的含水和岩层断裂情况，为判断地下径流、地下矿区走向提供依据；受地面山体起伏的影响，各测线的地面水平间距取40～80m，估计需40条测线左右。B组测线，西南→东北（每条测线平均长度约2500m），主要探测该地区东部的含水和岩层断裂分布情况，查询该区域多个泉水出露点的来源；受地面山体起伏的影响，各测线的地面水平间距取40～80m，估计需15条测线左右。C组测线，西南→东北（每条测线平均长度约800m），主要探测该地区在西北方向的含水和岩层断裂情况，为判断地下径流、地下矿区走向提供依据；受地面山体起伏的影响，各测线的地面水平间距取40～80m，估计需25条测线左右。因此，本实施例中，取测线条数K=80；根据工程的具体探测要求和踏勘的实际情况分析，在本实施例探测的垂深为地下48m到地下280m，每间隔8m做一次探测，则相应的探测深度分别为48m、56m、64m、72m、80m、88m、96m、104m、112m、120m、128m、136m、144m、152m、160m、168m、176m、184m、192m、200m、208m、216m、224m、232m、240m、248m、256m、264m、272m、280m，即测深个数J＝30；令探测区域的电阻率ρ＝10，根据电法勘探的深度与频率对应关系可得：测深h与探测仪器工作频率f _c的关系为f _c＝ (1591.58/h)²，则30个测深h对应探测仪器的工作频率f _c分别为1099.44Hz、807.75Hz、618.43Hz、488.64Hz、395.80Hz、327.11Hz、274.86Hz、234.20Hz、201.94Hz、175.91Hz、154.61Hz、136.95Hz、122.16Hz、109.64Hz、98.95Hz、89.75Hz、81.78Hz、74.82Hz、68.71Hz、63.33Hz、58.55Hz、54.29Hz、50.48Hz、47.06Hz、43.98Hz、41.19Hz、38.65Hz、36.35Hz、34.24Hz、32.31Hz，探测仪器的工作频率范围为1Hz ~6000Hz。现有三维物探一般是在X方向布置测线，探测Y方向的电场或磁场分量以进一步分析Z方向的地下地质体情况，且大都采样时域数据分析方法，数据量大、针对性不强、抗干扰性能差；本发明方法是在自有的频域三维物探方法基础上引入了随时间变化的特征参数而形成的基于天然电场的四维物探方法，对测点数据需要进行多次重复采集，综合考虑探测工作量和现场数据采集的可靠性，重复探测次数I一般在5~9次之间选取。本发明方法为减少数据的保存量，每次采集完成后实时进行频谱细化与分析，仅保存细化前的512点时域数据、细化后256条有效频域谱线值和特征信息；本发明的频谱细化采用ZFFT方法，因此，在数据采样前需确定以下探测参数：采样频率f _s=k₁+k₂*(f _c-k₃)，分析半带宽B=k₄+k₅* (f _c-k₃)，细化倍数D=f _s/2B，滤波器的半阶数M=4D，采样数据个数N=2M+511D；当1Hz <f _c≤10Hz时，取系数k₁=12.8、k₂=5.6889、k₃=1.0、k₄=0.6、k₅=0.26667；当10Hz <f _c≤100Hz时，取系数k₁=64.0、k₂=3.5556、k₃=10.0、k₄=3.0、k₅=0.03333；当100Hz <f _c≤500Hz时，取系数k₁=384.0、k₂=2.24、k₃=100.0、k₄=6.0、k₅=0.0175；当500Hz <f _c≤1200Hz时，取系数k₁=1280.0、k₂=1.77、k₃=500.0、k₄=13.0、k₅=0.01143；当1200Hz <f _c<6000Hz时，取系数k₁=3648.0、k₂=1.95、k₃=1200.0、k₄=16.0、k₅=0.002917。

步骤二 通过探测仪器同步采集每条测线上多个测点的天然电场信息，每个测深重复采集I次并实时进行频谱细化、特征信息识别后保存，直到采集完第J个测深后再进行下一条测线探测。

参见附图3，探测仪器由16根探针、2段信号电缆、1个探测装置组成；16根插入地表的探针用于接收天然瞬变电磁波穿透地层在地面形成的分布电场信号，分别经2段信号电缆（每段电缆软件8根探针、探针间距大小可探测任务要求选择）连接到探测装置；探测装置采用以DSP+CPLD为核心的结构设计，具有信号调理、16通道信号同步采集、ZFFT频谱细化和数据保存等功能；上位机读取探测装置保存的探测数据经综合统计、特征信息识别与分析后绘制测线四维物探剖面图，然后综合多条测线的探测分析结果得出结论并对探测目标地质体进行定性、定量评估。结合附图1，探测仪器每次进行1个测深的测线上16个测点的同步数据采集，探针的间距即是极距也是点距，对于普查时一般选用极距较大的信号电缆以加大测线探测的控制范围，祥查、精查选择极距较小的信号电缆以提高探测精度；数据采集完成，探测装置的DSP根据步骤一设置的探测参数对采集的时域信号实时进行ZFFT频谱细化、动态信息与裂隙信息识别并保存，每个测深重复进行I次探测后，探测装置自动根据步骤一设置进行下一个测深的探测，直到完成第J个测深的探测后结束该测线的探测任务；接着按相同的方法进行下一条测线探测。

本发明的频谱细化方法采用基于解析带通滤波的ZFFT方法，解析带通滤波器选用FIR非递归复解析带通滤波器，其细化过程按如下5步进行：

（1）确定中心频率及细化倍数 每次探测时欲细化频段的中心频率即为探测仪器的工作频率f _c，根据步骤一中的方法计算采样频率f _s、分析半带宽B、细化倍数D、滤波器的半阶数M=4D、采样数据个数N。

（2）构造复解析带通滤波器 设计一个复解析带通滤波器，滤波器的宽度为f _s/2D。

（3）选抽滤波 用设计好的复解析带通滤波器对N个采样数据x (n)作选抽滤波，并将滤波与选抽结合起来以提高计算效率，且实信号经过复滤波器以后就成为频率在通带以内的复解析信号；由于分析带宽内所要求的谱线数为256条，故选抽点数和FFT运算点数均为512点。

（4）复调制移频 对选抽后的g(m)进行复调制移频，将细化的起始频率f _c-B移到零频点后得到g ^/ (m)。

（5）FFT和谱分析 作512点FFT和谱分析，取正频率部分，不需要进行频率调整就可以得到具有256条独立谱线的细化频谱。

结合附图4、附图5、附图6，本发明特征信息识别方法是：针对探测仪器接收天然电场信号并细化得到分析带宽f _c-B～f _c+B的分辨率为∆f的256条谱线，在频谱中搜索出最大幅值F_MAX对应的谱线L_MAX及另外12个次大幅值谱线L_MAX1、L_MAX2、…、L_MAX12并同时计算频谱的平均幅值F_J、最大幅值谱线L_MAX及其前后各3条谱线一起的平均幅值F_J0；接着缓存256条谱线到另一指定存储区，在其中按从前往后两两比较幅值低者置0、从后往前两两比较幅值低者置0的原则找出频谱中的所有极值点，再将幅值小于F_MAX/30的极值点置0，结合原频谱计算次大极值点谱线J_MAX1、J_MAX2及其前后各3条谱线一起的平均幅值F_J1、F_J2；当极点数＜3进行裂隙信息识别，当3≤极点数＜30进行动态信息识别，极点数≥30无特征信息退出识别，此时的主频频率为仪器工作频率f _c，幅值为谱线L_MAX、L_MAX1、L_MAX2的平均幅值。

裂隙信息识别：当存在频谱能量衰落且极点谱线不分散时，若最大幅值谱线L_MAX频率非工频或其高次谐波、也不是测点区域已知的固有频率成分，则该频谱信息为裂隙信息，主频信息为最大幅值谱线L_MAX的幅值、频率，反之则无特征信息退出识别；频谱能量衰落判据是：若谱线L_MAX的幅值F_MAX≥20F_J并谱线L_MAX12的幅值F_MAX12≥F_J且12个次大幅值谱线L_MAX1、…、L_MAX12分布在最大幅值谱线L_MAX前后各6条谱线范围内，则频谱能量存在衰落现象；极点谱线分散判据是：若在最大幅值谱线L_MAX前后各6条谱线范围内存在2个以上极值点，则极点谱线分散。参见附图6，最大谱线L_MAX的幅值F_MAX≥20F_J且不属干扰信号，次大谱线L_MAX12的幅值F_MAX12≥F_J，12个次大幅值谱线L_MAX1、…、L_MAX12分布在最大幅值谱线L_MAX前后各6条谱线范围内，故附图6的频谱变化属于衰落现象，是地下空区形成的裂隙信息。

动态信息识别：在满足F_J0≥10F_J则频谱能量集中，判断最大幅值谱线L_MAX的幅值F_MAX0≥3F_J0则L_MAX离散、次大极值点谱线J_MAX1的幅值F_MAX1≥3F_J1则J_MAX1离散、次大极值点谱线J_MAX2的幅值F_MAX2≥3F_J2则J_MAX2离散，通过边频识别规则一、规则二、规则三判断是否存在关于L_MAX、J_MAX1及J_MAX2的对称极值谱线来确定地下径流动态信息的存在，否则无特征信息退出识别；边频识别规则一是：L_MAX离散且极值点分布于其左、右两边，若存在关于L_MAX对称的边频极值谱线，则主频为L_MAX，动态信息为一阶边频的信息，若不存在关于L_MAX对称的边频极值谱线、但存在关于J_MAX1或J_MAX2对称的边频极值谱线且相应的J_MAX1或J_MAX2离散，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，否则无特征信息退出识别；边频识别规则二是：L_MAX离散但主要极值点分布在其一边，若J_MAX1或J_MAX2离散且存在关于J_MAX1或J_MAX2对称的边频极值谱线，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，否则无特征信息退出识别；边频识别规则三是：L_MAX不离散，若J_MAX1或J_MAX2离散且存在关于J_MAX1或J_MAX2对称的边频极值谱线，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，若J_MAX1或J_MAX2均不离散但存在包括L_MAX、J_MAX1、J_MAX2中至少2个在内的4个以上的点距误差在±2条谱线内的极值谱线，则主频频率为L_MAX谱线频率、幅值为L_MAX与J_MAX1的均值，动态信息的频率为一阶边频频率、幅值为J_MAX1与J_MAX2的均值，否则无特征信息退出识别。在附图5中，明显可见最大幅值谱线L_MAX0的左边存在二阶边频极值点、右边存在三阶边频极值点，且它们都关于L_MAX0谱线对称；因此，主频信息为L_MAX0谱线所对应的频率f _c（＝f ₁＋f ₁与f _c的点距×∆f）、幅值（可从频谱中直接读取），动态信息为一阶边频的频率f _b（＝0.5×两个一阶边频的点距×∆f）、幅值（＝0.5×左一阶边频幅值＋0.5×右一阶边频幅值）。

步骤三，上位机读取探测装置保存的测线探测数据，综合统计测线上每个测深的每个测点的I次重复探测数据及其动态信息、裂隙信息的识别结果，绘制测点信号谱图和测线四维物探剖面图并进行测线低值异常判断，若测线存在动态信息、裂隙信息、低值异常则相应地进行地下径流、地下空区、断裂构造的定性、定量分析，否则进行下一条测线分析，直到完成第K条测线分析。

完成步骤二的测线探测任务后，上位机读取探测装置保存的探测数据，针对每条测线进行如下的统计与分析：

（1）测点信息观察验证 每个测点每次采集的存储数据包括：512 g ^/ (m)点值、256条谱线值和特征信息识别结果，通过每个测点的256条谱线值和特征信息识别结果可绘制测线上16测点每次采集的细化信号谱图如附图7所示；为验证探测仪器识别动态信息、裂隙信息的准确性或观察其细微变化特征，可双击对应窗口得到附图5或附图6所示的动态信息特征图或裂隙信息特征图，该图的上半部分为512点时域信号g ^/ (m)、下半部分为256条谱线，图中可清楚的观察到特征信息的细微变化。

（2）测线数据综合统计 如附图7所示每条测线包括16个测点，但其中的A8B8是整个测线范围的综合信息不代表一个测点，因此每条测线只有15个有效测点；附图8所示为本发明实施例的测线数据综合统计表，为消除探测过程中随机干扰、系统白噪声等的影响，对每一测点I次重复探测数据的识别结果中（由步骤二可知在识别结果中不管是否存在特征信息都会具有主频的幅值和频率）的主频的幅值和频率去掉一个最大值、去掉一个最小值、剩余的数求平均的方法求出主频的综合幅值和频率，附图8表中的第1列为各探测深度的主频综合频率，表中的上半部分为1条测线的15个测点在各个深度的主频幅值综合值，在综合数据统计表中的数据上、下、左、右间作比较可清楚地了解到测点范围内地下纵向及横向地质体变化的情况；若I次重复探测数据结果中有I ^/ 次识别出有特征信息且I ^/ /I≥0.6，则判定该测点在探测深度存在特征信息，附图8表中的下半部分为对应15个测点在各个深度的特征信息识别情况，特征信息中的动态信息用“Δ”表示、裂隙信息用“О”表示。

（3）绘制测线四维物探剖面图 附图8表中综合统计测线上每个测深的每个测点的I次重复探测数据及其动态信息、裂隙信息的识别结果，将表中的信息用图形方式表示出来则得到了本发明的四维物探剖面图；附图9所示为本发明的一个实施例的测线四维物探剖面图，横坐标为探针A8……B8的间距（单位：m），纵坐标为天然电磁波穿透不同地层后反映到地表的电性值，图中每一条曲线表示一个深度的地层电性值变化情况（用不同颜色、不同线型的线段加以区分），再将某一深度的地下特征信息（动态信息“Δ”、裂隙信息“О”等）标注到相应的地层电性值曲线上；图中包含了四维信息，故称四维物探剖面图；附图9测线一的探测频段频率为35Hz、40Hz、45Hz、55Hz、65Hz、75Hz、85Hz、95Hz、125Hz、165Hz、205Hz、255Hz、335Hz、385Hz、505Hz、625Hz等16个频段，可由公式深度h＝1591.58(1/f _z )^0.5换算为对应的探测深度；从附图9所示的测线四维物探剖面图可知，测线范围内无明显地下水动态信息、裂隙信息和电性低值异常区域，该测线在探测深度范围内应属完好地层结构。

（4）地下径流分析 对存在动态信息的测点进行地下径流分析。地下径流又称地下水径流、地下径流水，是重力作用下地下水在自然界陆地水循环过程中的流动。岩溶、裂隙中脉状流的地下水，按水文学规律流动在地底下会切割地磁场磁力线产生感应电动势f _b，但该感应电动势较微弱，无法穿透地层传送到地面；但它会对经过该处的天然电磁波f _c产生持续不断的干扰，最终附加在该电磁波信号上传输到地表，把这种附加有地下水流动信息的信号定义为地下径流动态信息。探测仪器通过在地表接收该天然电磁波的分布电场信号，并据此分析地下径流的赋存及状态。附图5是本发明实施例的动态信息，图中上半部分为天然瞬变电磁波的时域波形，下半部分为天然瞬变电磁波中指定分析带宽f ₁-B～f ₂+B的细化频谱图。频谱图中间的f _c谱线称为主频谱线，分布于主频谱线两边并关于主频对称的谱线f _b称为边频谱线。根据电磁场理论，天然瞬变电磁波在地层中的传导电流远远大于位移电流，则其穿透深度为

h＝1/(πf _c μ/ρ)^0.5(单位：m)

式中，ρ为探测点的地表电阻率，一般地表大都是松散的土层，取其平均电阻率ρ≈10（Ω.m）；假设天然瞬变电磁波穿透的地层为无磁性介质，则磁导率取 μ＝4π×10^-7 H/m；据此得到地下水径流的埋藏深度与主频频率的关系如下

h＝1/(πf _c μ/ρ)^0.5＝1591.58(1/f _c)^0.5(单位：m)

流动的地下水会切割地磁场的磁力线而产生感应电动势e，根据电磁流量原理，感应电动势e（单位：mV）与流量的方程式为

e＝（4BK_s/L）Q

式中，B为探测地点的大地电磁场感应强度（单位：mV/m ² ）、K_s为常数，对于一个已知的探测区域来说它们都是常数，但需要在探测区域附近找已知的水量的参考测试点进行整定；L为地下水流通道截面的周长（单位：m）。实际探测工作表明，地下水在岩溶裂隙中流动的流速、流量都是随时间变化的，其产生的感应电动势e是时间t的函数，计为e (t)＝0.5〔e ₁(t)＋e ₂(t)〕。e (t)与地下水储量Sw之间的关系为

Sw=∫₀ ^T[L.e (t)/4BK]dt=( L/4BK)∫₀ ^T e (t) dt(单位：m ³/h)

式中，T为地下水径流的变化周期，探测仪器能探测出地下水径流的变化周期T（边频频率f _b的倒数）和周期T内感应电动势e (t)的综合值，由上式可简单地估算出地下水储量Sw。但由于地下水径流的变化周期T受地质构造、补给情况、开采量等因素的影响而各不相同，且该周期也是随时间变化的，这给准确评估地下水储量或开采量带来一定困难。通常采用离散化方法来计算上式，具体做法是：由探测仪器对该固定深度连续重复探测1小时，得到M组探测数据，而这其中只有N组探测数据具有地下水径流特征信息，即边频综合幅值e (t ₀)、e (t ₁)、e (t ₂)……e (t _N-1)，边频频率f _b0、f _b1、f _b2……f _bN-1；因此，上式的离散化计算公式为：

(单位：m ³/h)

因此，本发明方法通过步骤一得出有无地下径流动态信息的结论并相应给出主频信息、动态信息的幅值、频率，即可由识别的主频f _c的频率与幅值、边频f _b的频率与幅值来定量评估地下径流。

（5）地下空区分析 对存在裂隙信息的测点进行地下空区分析。地下空区是地下溶洞、断裂构造形成的地下空间、矿产资源开采后形成的老空区以及各种地下掩埋空洞的总称；所谓老空区是指矿产资源开采后形成的采空区、老窑和已经报废井的总称。根据天然瞬变电磁波的频率和不同地质体的电磁波穿透情况，由于地下空区与完好地层间存在密度突变，天然瞬变电磁波在这里会引起信号的反射、折射、绕射、散射和吸收等现象，导致信号产生衰落，从而降低了信号的传输质量；这个由地下空区所形成的天然瞬变电磁波传输方式相当于一个窄带随机过程：只允许靠近中心频率ω _c(＝2πf _c)附近很窄范围(Δω<<ω _c)的频率成分通过窄带系统，随机信号通过窄带系统后，输出信号即为窄带随机信号，而且是存在衰落现象的窄带随机信号，这就是步骤二中裂隙信息识别的理论依据。天然电磁波穿过地下空区的窄带随机过程的垂直分解可知X(t)和Y(t)在几何上是正交的，变换到频域后其实部X(f)和虚部Y(f)也是正交的；附图6下半部分所示是探测信号的细化幅值谱，图中的每一条的幅值为A(f)=[ X(f)²+Y(f)²]^0.5，相位为：θ(f)=tan^－1 [Y(f)/X(f)]。如附图6所示，设最大幅值谱线L_MAX（即主频f _c）的相位为θ _MAX (f)，其衰落到幅值小于频谱平均幅值F_J的第一条谱线的相位为θ _MIN (f)，则衰落相位差为Δθ(f)=θ _MAX (f)-θ _MIN (f)；由此可计算出主频f _c信号衰落时间差为

Δt(f)= Δθ(f)/2πf _c=[θ _MAX (f)-θ _MIN (f)] /2πf _c(单位：s)

天然电磁波的波速与光速C=3*10⁸m/s相当，故可通过下式来评估地下空区的横向尺寸B_x大小

B_x=C*Δt(f) =4.775*10⁷*[θ _MAX (f)-θ _MIN (f)] / f _c(单位：m)

本发明方法的探测仪器是在地面天然电磁波信号，则天然电磁波的各频率成分穿过地下空区传输到地面的路径是基本一致的，令其传输路径损耗为K(dB) 。根据电波在自由空间（所谓自由空间传播系指天线周围为无限大真空时的电波传播，电波在自由空间传播时，其能量既不会被障碍物所吸收，也不会产生反射或散射）传播损耗（亦称衰减）计算公式

Los=32.44+20*lg d+20*lg f

式中，Los是传播损耗，单位为dB；d是传播距离，单位是Km；f是工作频率，单位是MHz。结合本发明方法的应用实际，上式可改写为

Los+K=23.44+20*lg d+20*lg f _c

式中，d是传播距离，单位是m；f _c是主频频率，单位是Hz；Los和K的单位均为dB，K是传输路径公共损耗，对于同一探测区域的多个测点在相同深度探测时其值是基本一致的，Los是天然电磁波穿过地下空区的衰落损耗，在本发明方法中，可由最大幅值谱线L_MAX（即主频f _c）的F_MAX和主频f _c衰落到幅值小于频谱平均幅值F_J的第一条谱线的F_MIN求得

Los=20*lg (F_MAX/F_MIN)

因此，可由该衰落损耗来评估地下空区的纵向尺寸B_y大小(单位：m)

B_y=10^{[20*lg (FMAX/FMIN)+K-23.44-20*lg fc]/20}=10^{[lg (FMAX/FMIN)+K/20-1.172- lg fc]}

以附图6实施例的探测结果为例，主频频率f _c=169.85Hz，主频幅值实部X_MAX (f) =0.69856mV、虚部Y_MAX (f) =0. 0.66429mV、模A_MAX (f) =F_MAX=0.963985mV，衰落到幅值小于频谱平均幅值F_J的第一条谱线的幅值实部X_MIN (f) =0.003698mV、虚部Y_MIN (f) =0.003516mV、模A_MIN (f) =F_MIN=0.0051027mV，令传输路径公共损耗K=30dB，根据前述计算公式可得该测点地下空区的定量分析结果如下：

地下空区的埋藏深度h=1591.58(1/f _c)^0.5=122.12m；

地下空区的横向尺寸B_x= 4.775*10⁷*[θ _MAX (f)-θ _MIN (f)] / f _c=23.287m；

地下空区的纵向尺寸B_y=10^{[lg (FMAX/FMIN)+K/20-1.172-lg fc]} =2.367m。

通过本发明方法的大量试验和现场数据分析，发现地下空区具有以下特点：①当地下空区未塌陷且不含水时，将引起特征信息主频幅值的高值异常，而保存完好的地层主频幅值相对较低；②当地下空区塌陷且充填了泥砂或者其它较软杂质时，会引起相对低的特征信息主频幅值，而保存完好的地层主频幅值相对稍高；③当地下空区塌陷（或未塌陷）且含水时，将引起主频幅值低值异常，而保存完好的地层主频幅值较高。因此，根据多测点特征信息主频幅值F_MAX的变化规律可有效评估地下空区的存在状态；另外对比特征信息主频幅值F_MAX的综合值与周围测点的变化情况，结合附近已知的钻孔资料，还可进一步评估测点在该深度的地下岩性。

（6）断裂构造分析 对四维物探曲线中整体出现电性低值异常的测线进行断裂构造分析。断裂构造又称断裂，是指组成地壳的岩体在地应力作用下发生变形，当应力超过岩石的强度，岩体的完整性和连续性受到破坏而产生的大小不一的断裂；断裂构造是岩石破裂的总称，包括劈理、节理、断层、深大断裂和超壳断裂等；断裂可以作为石油天然气二次运移的良好通道，油气沿断裂通道运移比在岩石孔隙中运移更加容易；研究断裂构造对找矿勘探、水文地质与工程地质以及了解区域构造特点均有实际意义。附图10测线二的探测频段频率为30Hz、40Hz、55Hz、75Hz、95Hz、125Hz、165Hz、225Hz、305Hz、355Hz、475Hz、685Hz、885Hz、1085Hz等14个频段，在测线二的探针A7A6之间四维物探曲线出现电性低值异常区域，且在75Hz（对应探测深度约183米）频段存在地下水动态信息、305Hz（对应探测深度约91米）频段存在裂隙信息。所谓的电性低值异常是指四维物探曲线中部分曲线在一定连续测深范围内或全部曲线在整个测深范围内出现部分测点的探测电性值比测线上其它测点的探测电性值低2个及2个以上数量级。附图11测线三的探测频段频率为30Hz、40Hz、55Hz、75Hz、95Hz、125Hz、165Hz、225Hz、305Hz、385Hz、505Hz、625Hz等12个频段，在测线三的探针B1B2之间四维物探曲线出现电性低值异常区域，且在95Hz（对应探测深度约163米）频段存在地下水动态信息、305Hz（对应探测深度约91米）频段存在裂隙信息。当地层中出现断裂时，在断裂带中充填的物质与两边岩层相比较松软、或者断裂带内含水时，会引起相对低的主频幅值反映，而两边岩层主频幅值相对较高；断裂带电阻率与岩石孔隙度的定量关系可以用阿尔奇公式来表述

R ₀ = aR _W Φm

式中，R ₀为岩石的宏观电阻率，R _W为孔隙中水的电阻率，Φ为岩石孔隙度，a和m为待定系数。砂岩的电阻率一般在几欧姆米至几千欧姆米变化，分选差、颗粒粗及胶结程度高的致密砂岩，其电阻率高； 反之，分选好、颗粒细及胶结程度低的疏松砂岩往往具有低电阻率。砾岩由于颗粒粗、分选性差，故常具有比砂岩较高的电阻率。一般土层结构疏松，孔隙度大，且与地表水密切相关，因而它们的电阻率均较低，一般为几十Ω·m。在本发明仪器的探测中：被测地质体的电阻率与四维物探曲线的电性值成正比例关系。所以当在某一测点位置出现纵向上所有曲线出现连续性低值异常时，则表明此处有断层存在。综合分析附图10和附图11，由于测线二和测线三的布置相距仅300米，可初步判定在测线二的A7A6测点到测线三的B1B2测点之间存在一地下断裂构造，且断裂带宽度小于10m。具体断裂带宽度和结构需通过测线平移或更换极距较小的信号电缆进行精确定位探测、并多次重复探测才能得到准确结果。因此，识别地下断裂构造走向需要进行2条以上的测线探测；要更精确的掌握断裂构造宽度、含水情况等，则需进行精确定位探测、多次重复探测。单凭1条测线的数据或地层中某一段深度的数据进行判断是不完整的、不全面的。由于测线二的A7A6测点在183米处存在地下水动态信息、垂深91米存在裂隙信息，测线三的B1B2测点在163米处存在地下水动态信息、垂深91米存在裂隙信息；可知该断裂构造为一导水断层，导水流向是从测线三的B1B2测点到测线二的A7A6测点，且断裂带的顶部结构疏松，存在孔隙。

步骤四，综合K条测线探测分析结果得出探测工程结论，对探测目标地质体进行综合评估。

在完成步骤二的测线探测、步骤三的测线数据综合与分析后，通过综合K条测线探测分析结果，即可得出探测工程结论，对存在特征信息或探测目标地质体做进一步详细分析与评估。

参见图12，为本发明方法实施例的测线四的四维物探剖面图，由步骤二中的分析可知四维物探剖面图的表示方法对分析地下断裂构造是非常有用，但在该图中它将所有的地层电性曲线画在一个图中，不能清楚地表示出地层结构的细微变化。为此，本发明方法还可根据综合分析的需要、针对不同的探测目的对四维物探剖面图进行横向、纵向分解，图13所示为图12的横向分解图，图14所示为图12的纵向分解图。综合分析图12～图14可知：在垂深68m层次测线的探针B1到探针B5之间曲线存在低值异常，且B2B3测点有地下水动态信息，B4B5测点有地下裂隙信息。为更清楚的了解该层次的结构变化，本发明的软件还可对图12中垂深68m的单层次分解，如图15所示。从图15可清楚看到：测线的探针B1到探针B5之间垂深68m处存在老空区，且该老空区已部分塌陷，底部存在流动水，顶部还有空隙；因此，还需结合步骤三中的关于低值异常、动态信息、裂隙信息的分析方法对其做进一步的定性、定量综合分析。

综上所述，本发明方法是在现场数据采集前对探测工程进行完备的设计，通过探测仪器在地面同步采集每条测线上多根探针接收的天然电磁波穿透地层结构后在地表形成的分布电场信号并实时进行频谱细化、特征信息识别后保存，上位机读取探测数据经综合统计后绘制测线四维物探剖面图并对存在特征信息的测线进行定性、定量分析，然后综合多条测线的探测分析结果得出探测工程结论。本发明的有益效果是：在有效运用频域三维物探方法的基础上引入随时间变化的特征参数而形成一种基于天然电场的四维物探方法，与现有物探方法相比：探测仪器所需存储空间少、同步探测测点数据的关联性好且抗干扰能力强，能对探测目标地质体进行定性、定量的综合分析；能满足不同探测目的的物探工程需要。

Claims (10)

1.一种基于天然电场的四维物探方法，其特征在于：在频域三维物探基础上引入随时间变化的特征参数，通过探测仪器在地面同步采集每条测线上多根探针接收的天然电磁波穿透地层结构后在地表形成的分布电场信号并实时进行频谱细化、特征信息识别后保存，由上位机读取探测数据经综合统计后绘制测线四维物探剖面图并对存在特征信息的测线进行定性、定量分析，然后综合多条测线的探测分析结果得出探测工程结论；分以下四个步骤进行：

步骤一，按照探测工程的具体要求，设计能满足探测要求的K条测线以及每条测线的测深个数J、重复探测次数I和相应的探测参数；

步骤二，通过探测仪器同步采集每条测线上多个测点的天然电场信息，每个测深重复采集I次并实时进行频谱细化、特征信息识别后保存，直到采集完第J个测深后再进行下一条测线探测；

步骤三，上位机读取探测装置保存的测线探测数据，综合统计测线上每个测深的每个测点的I次重复探测数据及其动态信息、裂隙信息的识别结果，绘制测点信号谱图和测线四维物探剖面图并进行测线低值异常判断，若测线存在动态信息、裂隙信息、低值异常则相应地进行地下径流、地下空区、断裂构造的定性、定量分析，否则进行下一条测线分析，直到完成第K条测线分析；

步骤四，综合K条测线探测分析结果得出探测工程结论，对探测目标地质体进行综合评估。

2.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的测深h与探测仪器工作频率关系为f _c＝ (1591.58/h)²；所述的探测参数包括采样频率f _s=k₁+k₂*(f _c-k₃)、分析半带宽B=k₄+k₅* (f _c-k₃)、细化倍数D=f _s/2B、滤波器的半阶数M=4D、采样数据个数N=2M+511D；当1Hz <f _c≤10Hz时，取系数k₁=12.8、k₂=5.6889、k₃=1.0、k₄=0.6、k₅=0.26667；当10Hz <f _c≤100Hz时，取系数k₁=64.0、k₂=3.5556、k₃=10.0、k₄=3.0、k₅=0.03333；当100Hz <f _c≤500Hz时，取系数k₁=384.0、k₂=2.24、k₃=100.0、k₄=6.0、k₅=0.0175；当500Hz <f _c≤1200Hz时，取系数k₁=1280.0、k₂=1.77、k₃=500.0、k₄=13.0、k₅=0.01143；当1200Hz <f _c<6000Hz时，取系数k₁=3648.0、k₂=1.95、k₃=1200.0、k₄=16.0、k₅=0.002917。

3.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的探测仪器由16根探针、2段信号电缆、1个探测装置组成；16根插入地表的探针用于接收天然瞬变电磁波穿透地层在地面形成的分布电场信号，分别经2段信号电缆（每段电缆软件8根探针、探针间距大小可探测任务要求选择）连接到探测装置；探测装置采用以DSP+CPLD为核心的结构设计，具有信号调理、16通道信号同步采集、ZFFT频谱细化、特征信息识别和数据保存等功能。

4.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的特征信息识别方法是：针对探测仪器接收天然电场信号并细化得到分析带宽f _c-B～f _c+B的分辨率为∆f的256条谱线，在频谱中搜索出最大幅值F_MAX对应的谱线L_MAX及另外12个次大幅值谱线L_MAX1、L_MAX2、…、L_MAX12并同时计算频谱的平均幅值F_J、最大幅值谱线L_MAX及其前后各3条谱线一起的平均幅值F_J0；接着缓存256条谱线到另一指定存储区，在其中按从前往后两两比较幅值低者置0、从后往前两两比较幅值低者置0的原则找出频谱中的所有极值点，再将幅值小于F_MAX/30的极值点置0，结合原频谱计算次大极值点谱线J_MAX1、J_MAX2及其前后各3条谱线一起的平均幅值F_J1、F_J2；当极点数＜3进行裂隙信息识别，当3≤极点数＜30进行动态信息识别，极点数≥30无特征信息退出识别，此时的主频频率为仪器工作频率f _c，幅值为谱线L_MAX、L_MAX1、L_MAX2的平均幅值。

5.根据权利要求4所述的基于天然电场的四维物探方法，其特征在于：所述的裂隙信息识别为：当存在频谱能量衰落且极点谱线不分散时，若最大幅值谱线L_MAX频率非工频或其高次谐波、也不是测点区域已知的固有频率成分，则该频谱信息为裂隙信息，主频信息为最大幅值谱线L_MAX的幅值、频率，反之则无特征信息退出识别；频谱能量衰落判据是：若谱线L_MAX的幅值F_MAX≥20F_J并谱线L_MAX12的幅值F_MAX12≥F_J且12个次大幅值谱线L_MAX1、…、L_MAX12分布在最大幅值谱线L_MAX前后各6条谱线范围内，则频谱能量存在衰落现象；极点谱线分散判据是：若在最大幅值谱线L_MAX前后各6条谱线范围内存在2个以上极值点，则极点谱线分散。

6.根据权利要求4所述的基于天然电场的四维物探方法，其特征在于：所述的动态信息识别为：在满足F_J0≥10F_J则频谱能量集中，判断最大幅值谱线L_MAX的幅值F_MAX0≥3F_J0则L_MAX离散、次大极值点谱线J_MAX1的幅值F_MAX1≥3F_J1则J_MAX1离散、次大极值点谱线J_MAX2的幅值F_MAX2≥3F_J2则J_MAX2离散，通过边频识别规则一、规则二、规则三判断是否存在关于L_MAX、J_MAX1及J_MAX2的对称极值谱线来确定地下径流动态信息的存在，否则无特征信息退出识别；边频识别规则一是：L_MAX离散且极值点分布于其左、右两边，若存在关于L_MAX对称的边频极值谱线，则主频为L_MAX，动态信息为一阶边频的信息，若不存在关于L_MAX对称的边频极值谱线、但存在关于J_MAX1或J_MAX2对称的边频极值谱线且相应的J_MAX1或J_MAX2离散，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，否则无特征信息退出识别；边频识别规则二是：L_MAX离散但主要极值点分布在其一边，若J_MAX1或J_MAX2离散且存在关于J_MAX1或J_MAX2对称的边频极值谱线，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，否则无特征信息退出识别；边频识别规则三是：L_MAX不离散，若J_MAX1或J_MAX2离散且存在关于J_MAX1或J_MAX2对称的边频极值谱线，则主频为J_MAX1或J_MAX2，动态信息为一阶边频的信息，若J_MAX1或J_MAX2均不离散但存在包括L_MAX、J_MAX1、J_MAX2中至少2个在内的4个以上的点距误差在±2条谱线内的极值谱线，则主频频率为L_MAX谱线频率、幅值为L_MAX与J_MAX1的均值，动态信息的频率为一阶边频频率、幅值为J_MAX1与J_MAX2的均值，否则无特征信息退出识别。

7.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的低值异常是指四维物探曲线中部分曲线在一定连续测深范围内或全部曲线在整个测深范围内出现部分测点的探测电性值比测线上其它测点的探测电性值低2个及2个以上数量级。

8.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的定性分析是根据多测点主频f _c幅值F_MAX的变化规律来评估地下空区的塌陷、充填情况；通过多测线的四维物探曲线低值异常来确定断裂构造走向并结合动态信息、裂隙信息参数来综合分析断裂构造的存在状态。

9.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的定量分析是根据主频频率f _c、相位θ _MAX (f)、幅值F_MAX，动态信息的边频频率f _b、边频综合幅值e (t)和裂隙信息的幅值小于频谱平均幅值F_J的第一条衰落谱线的相位θ _MIN (f)、幅值F_MIN以及传输路径公共损耗K来确定地下径流、地下空区、断裂构造的埋藏深度h=1591.58(1/f _c)^0.5，

地下径流或导水断层的水量，

地下空区的横向尺寸B_x= 4.775*10⁷*[θ _MAX (f)-θ _MIN (f)] / f _c，

地下空区的纵向尺寸B_x=10^{[lg (FMAX/FMIN)+K/20-1.172-lg fc]}。

10.根据权利要求1所述的基于天然电场的四维物探方法，其特征在于：所述的目标地质体综合评估是对四维物探剖面图进行分解，综合运用低值异常、动态信息、裂隙信息的分析方法对其做出的定性、定量分析。