CN102116872B - 一种stratagem大地电磁测量系统阻抗张量的稳健估算方法 - Google Patents
一种stratagem大地电磁测量系统阻抗张量的稳健估算方法 Download PDFInfo
- Publication number
- CN102116872B CN102116872B CN200910216980.8A CN200910216980A CN102116872B CN 102116872 B CN102116872 B CN 102116872B CN 200910216980 A CN200910216980 A CN 200910216980A CN 102116872 B CN102116872 B CN 102116872B
- Authority
- CN
- China
- Prior art keywords
- impedance
- value
- calculate
- frequency
- independent observation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明属于一种地球物理数据处理方法,具体涉及一种STRATAGEM大地电磁测量系统(EH4测量系统)阻抗张量的稳健(Robust)估算方法。该方法把统计学中的最大似然估计(M估计)引入到STRATAGEM系统的独立观测中,使多次独立观测的数据均参与最终阻抗运算,但是每次观测参与计算的权值(重)不同。从而避免了在最小二乘估算中某次独立观测影响整体计算结果的可能性,提高了估算的可靠性与稳健性。
Description
技术领域
本发明属于一种地球物理数据处理方法,具体涉及一种STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法。
背景技术
STRATAGEM(EH4)大地电磁测量系统是上世纪末期引入我国的先进大地电磁(音频)测量设备。该设备目前在国内有一百余套,在寻找地下水、资源勘探中有广泛的用途。但是该系统所配的数据处理软件(IMAGEM)仍然采用较原始的大地电磁阻抗估算方法--最小二乘方法。由于个别时段的干扰信号影响最小二乘估算结果,不同时段组合运算所得到的阻抗曲线可能有较大差异,从而导致两个方面结果:
(1)阻抗运算结果飞点数量多;
(2)阻抗运算结果不稳定。
而调整系统阻抗估算流程,并融合稳健统计方法所计算出的阻抗曲线则可以有效地抑制这种情况的发生,为大地电磁资料的解释提供较为可靠的依据,是长时间以来广大用户迫切希望解决的实际应用问题。
发明内容
本发明的目的在于针对现有大地电磁阻抗估算方法的缺陷,提供一种STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,以提高估算的可靠性。
为实现上述目的,本发明的技术方案如下:一种STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,包括如下步骤:
(1)读取原始时间序列数据,获取单个测点的多次独立观测时间序列原始结果;
(2)用FFT变换方法,实现由时间域到频率域的转化,得到每个频点多次独立观测的频谱值;
(3)采用最大似然估计的方法,对同一频点的多个独立观测结果给出不同的权重;
(4)进行原始频点谱的叠加,实现由FFT变换得到的原始频点到阻抗频点的叠加;
(5)由每个阻抗频点叠加频谱值计算不同测量场的互功率谱;
(6)通过融合标定使测量场转化为地球物理单位意义的电磁场数值,从而计算出阻抗值。
进一步,如上所述的STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,步骤(3)中所述的最大似然估计的方法包括如下步骤:
①用最小二乘法计算同一频点的多个独立观测结果的初始阻抗;
②计算初始尺度估计参数;
③计算权函数w(r),并计算迭代更新电场值数据;
④采用最小二乘法计算同一频点的多个独立观测结果的新的阻抗张量与误差尺度;
⑤判断残差是否达到设定值以及是否达到设定的迭代次数;如果否,则返回步骤③,如果是则进行后续的阻抗计算。
更进一步,如上所述的STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,步骤⑤中残差r的表达式采用如下形式:
其中,Ei为实测电场,Hi为实测磁场,为阻抗张量,σ为尺度估计参数。
更进一步,如上所述的STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,步骤③中以权函数w(r)的权,迭代后的电场值表示为:
其中,预测电场值 riO为预测电场值及其与实测电场值之间的残差,分别为阻抗张量与尺度估计参数σ初始值由最小二乘方法对N次独立观测所得的结果,HT i为估算阻抗张量。
本发明的有益效果如下:本发明把统计学中的最大似然估计(M估计)引入到STRATAGEM系统的独立观测中,使多次独立观测的数据均参与最终阻抗运算,但是每次观测参与计算的权值(重)不同。从而避免了在最小二乘估算中某次独立观测影响整体计算结果的可能性,提高了估算的可靠性。
附图说明
图1为本发明的总体方法流程图;
图2为最大似然估计(M估计)的方法流程图;
图3为对存在干扰的STRATAGEM系统的时序数据采用本发明的处理方法得到的阻抗曲线;
图4为STRATAGEM系统配置的IMAGEM处理方法得到的阻抗曲线。
具体实施方式
下面结合附图和实施例对本发明进行详细的描述。
本发明的主要创新点是把统计学中的最大似然估计(M估计)引入到STRATAGEM系统的独立观测中,使多次独立观测的数据均参与最终阻抗运算,其它有关FFT变换(傅里叶变换),原始频点谱的叠加,互功率谱的计算以及阻抗值的计算均为本领域的公知技术。
如图1所示,一种STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,包括如下步骤:
(1)读取原始时间序列数据,获取单个测点的多次独立观测时间序列原始结果;
(2)用FFT变换方法,实现由时间域到频率域的转化,得到每个频点多次独立观测的频谱值;
(3)采用最大似然估计(M估计)的方法,对同一频点的多个独立观测结果给出不同的权重;
(4)进行原始频点谱的叠加,实现由FFT变换得到的原始频点到阻抗频点的叠加;
(5)由每个阻抗频点叠加频谱值计算不同测量场(电、磁场)的互功率谱;
(6)通过融合标定使测量场转化为地球物理单位意义的电磁场数值,从而计算出阻抗值。
如图2所示,最大似然估计(M估计)的方法包括如下步骤:
①用LS方法(最小二乘法)计算同一频点的多个独立观测结果的初始阻抗;
②计算初始尺度估计参数;
③计算权函数w(r),并计算迭代更新电场值数据;
④采用LS方法(最小二乘法)计算同一频点的多个独立观测结果的新的阻抗张量与误差尺度;
⑤判断残差是否达到设定值以及是否达到设定的迭代次数;如果否,则返回步骤③,如果是则进行后续的阻抗计算。
在最大似然估计(M估计)过程中主要函数选取如下:
●损失函数
损失函数采用Huber(1981)的分段函数,残差r的表达式如下式所示:
其中,Ei为实测电场,Hi为实测磁场,为阻抗张量,σ为尺度估计参数。
阻抗张量与σ尺度估计参数初始值由最小二乘方法对N次独立观测给出,分别为则,对于第i次观测值而言,其预测电场值及其与实测电场值之间的残差ri0为:
预测电场值 其中,HT i为估算阻抗张量。
权函数w(r)形式为损失函数的导数与残差r的比值。
则以权函数w(r)的权,迭代后的电场值可以表示为:
●尺度参数选取
尺度参数的选取与残差r的计算紧密联系,选取残差的均方根值计算方法,具体公式如下:
其中,mjj磁场帽子矩阵对角元素,n表示独立观测次数,p表示自由度。
图3为对存在干扰的STRATAGEM系统的时序数据采用本发明的处理方法得到的阻抗曲线,图4为STRATAGEM系统配置的IMAGEM处理方法得到的阻抗曲线。从两幅图的比较可以看出本发明提供的估算方法更稳健,估算的可靠性更高。
Claims (1)
1.一种STRATAGEM大地电磁测量系统阻抗张量的稳健估算方法,包括以下步骤:
(1)读取原始时间序列数据,获取单个测点的多次独立观测时间序列原始结果;
(2)用FFT变换方法,实现由时间域到频率域的转化,得到每个频点多次独立观测的频谱值;
(3)采用最大似然估计的方法,对同一频点的多个独立观测结果给出不同的权重,包括如下步骤:
①用最小二乘法计算同一频点的多个独立观测结果的初始阻抗;
②计算初始尺度估计参数;
③计算权函数w(r),并计算迭代更新电场值数据,以权函数w(r)的权,迭代后的电场值表示为:
其中,预测电场值HT i为估算阻抗张量,ri0为预测电场值及其与实测电场值之间的残差,分别为阻抗张量与尺度估计参数σ初始值由最小二乘方法对N次独立观测所得的结果;
④采用最小二乘法计算同一频点的多个独立观测结果的新的阻抗张量与误差尺度;
⑤判断残差是否达到设定值以及是否达到设定的迭代次数;如果否,则返回步骤③,如果是则进行后续的阻抗计算;残差r的表达式采用如下形式:
其中,Ei为实测电场,Hi为实测磁场,为阻抗张量,σ为尺度估计参数,σ具体公式为:
其中,mjj磁场帽子矩阵对角元素,n表示独立观测次数,p表示自由度;
(4)进行原始频点谱的叠加,实现由FFT变换得到的原始频点到阻抗频点的叠加;
(5)由每个阻抗频点叠加频谱值计算不同测量场的互功率谱;
(6)通过融合标定使测量场转化为地球物理单位意义的电磁场数值,从而计算出阻抗。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN200910216980.8A CN102116872B (zh) | 2009-12-31 | 2009-12-31 | 一种stratagem大地电磁测量系统阻抗张量的稳健估算方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN200910216980.8A CN102116872B (zh) | 2009-12-31 | 2009-12-31 | 一种stratagem大地电磁测量系统阻抗张量的稳健估算方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN102116872A CN102116872A (zh) | 2011-07-06 |
CN102116872B true CN102116872B (zh) | 2014-10-01 |
Family
ID=44215726
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN200910216980.8A Active CN102116872B (zh) | 2009-12-31 | 2009-12-31 | 一种stratagem大地电磁测量系统阻抗张量的稳健估算方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN102116872B (zh) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104570131B (zh) * | 2014-12-10 | 2017-03-08 | 中国船舶重工集团公司第七二二研究所 | 一种估计大地电磁参数的方法和装置 |
TWI684781B (zh) * | 2018-08-09 | 2020-02-11 | 國立中央大學 | 大地電磁時序資料處理方法及其系統 |
CN109188542B (zh) * | 2018-11-12 | 2020-04-14 | 国科(重庆)仪器有限公司 | 一种波区相关性检测的远参考大地电磁阻抗计算方法 |
CN111273367B (zh) * | 2020-03-11 | 2021-01-08 | 中南大学 | 一种大地电磁阻抗的估算方法 |
CN113341373B (zh) * | 2021-05-31 | 2024-05-14 | 中国电子科技集团公司第三十六研究所 | 一种定位方法、装置和电子设备 |
Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN100480734C (zh) * | 2007-03-08 | 2009-04-22 | 刘俊昌 | 一种高分辨率去静态频率域大地电磁法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH06274344A (ja) * | 1993-03-23 | 1994-09-30 | Chugoku Nippon Denki Software Kk | 推論戦略自動生成方式 |
-
2009
- 2009-12-31 CN CN200910216980.8A patent/CN102116872B/zh active Active
Patent Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN100480734C (zh) * | 2007-03-08 | 2009-04-22 | 刘俊昌 | 一种高分辨率去静态频率域大地电磁法 |
Non-Patent Citations (3)
Title |
---|
EH-4系统的数据二次处理技术及应用;化希瑞等;《地球物理学进展》;20080831;第23卷(第4期);1261~1268 * |
JP特开平6-274344A 1994.09.30 |
化希瑞等.EH-4系统的数据二次处理技术及应用.《地球物理学进展》.2008,第23卷(第4期),1261~1268. |
Also Published As
Publication number | Publication date |
---|---|
CN102116872A (zh) | 2011-07-06 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN102116872B (zh) | 一种stratagem大地电磁测量系统阻抗张量的稳健估算方法 | |
CN103248043B (zh) | 一种基于同步相角测量装置的电力系统多区域分布式状态估计方法 | |
CN103326358B (zh) | 基于同步相角测量装置的电力系统动态状态估计方法 | |
CN101776710B (zh) | 一种高压直流输电线路雷电绕击电流波形反演恢复方法 | |
CN104330471B (zh) | 航空结构损伤的Lamb波时变概率模型监测方法 | |
CN104155580B (zh) | 一种关联分析与电力计算相结合的电压暂降源定位方法 | |
CN102163844B (zh) | 基于相量测量装置的电力系统状态检测方法 | |
CN102749521A (zh) | 一种电力系统谐波阻抗计算方法 | |
CN103222151A (zh) | 在大规模电系统应用中的数据对准 | |
CN105548718A (zh) | 一种基于混合整体最小二乘法的系统谐波阻抗计算方法 | |
CN102636693A (zh) | 一种结合fft与非线性最小二乘的谐波分析算法 | |
CN104331846A (zh) | 一种窃电行为多源建模与协同分析方法 | |
CN103198184A (zh) | 一种电力系统低频振荡特征类噪声辨识方法 | |
CN104408295A (zh) | 一种大跨桥梁下部结构风-浪耦合作用荷载数值模拟方法 | |
CN102023010A (zh) | 基于mems的小波域多传感器信息融合系统及融合方法 | |
CN103336866A (zh) | 一种电磁暂态仿真中含负电阻支路的处理方法 | |
CN102510060B (zh) | 一种电力系统频率特性系数的计算方法 | |
CN104897958A (zh) | 一种输电线路雷击类型的辨识方法 | |
CN102333052A (zh) | 一种适用于浅海低频条件的水声信号盲解卷方法 | |
CN104122594A (zh) | 时间域激电全波形采样的多参数提取数据处理方法 | |
CN107894552A (zh) | 一种故障行波检测方法 | |
CN102944901B (zh) | 一种大地电磁阻抗估计方法 | |
CN104777356A (zh) | 一种基于神经网络的实时高精度谐波检测方法 | |
Couch et al. | Large-scale physical response of the tidal system to energy extraction and its significance for informing environmental and ecological impact assessment | |
CN103560509B (zh) | 一种基于小波分析的电压下陷检测装置及该装置的控制方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant |