JP2011242231A - 地下密度構造の推定方法 - Google Patents

Abstract

【課題】深さ方向を含めて地下の密度構造を一意的に求めることができる地下密度構造の推定方法を提供する。
【解決手段】探査機器位置情報、地形情報、反射断面の情報、地震波屈折法の情報、総合化測定データ群の情報とから、地中の地震波速度の水平・深さ方向の分布を示す地震波速度構造モデルの第一段階と、密度と地震波速度の間の最適関係式を用いて初期密度構造モデルの第二段階と、計算重力ベクトル値と、観測重力ベクトル値とを比較し、修正密度構造モデルの第三段階と、計算重力テンソル値と観測から得られる観測重力テンソル値とを比較し、最適密度構造モデルの第四段階と、前記最適密度構造モデルを地震波速度構造モデルに変換し、前記第一段階から第四段階を繰り返し、理論地震波屈折波、理論地震波広角反射波、理論地震波反射断面を作成し、所定の収束条件を満たす最適密度構造を推定密度構造として出力する第五段階と、を有する。
【選択図】図１

Description

本発明は、地下や海底下の鉱物資源、石油、天然ガス、メタンハイドレート、地下水脈、断層帯、地下埋設構造等によって生じる密度異常による地下密度構造の推定方法に関する。

特に、本発明は探査機器の経度、緯度、高度や水深などの位置情報と、探査対象域の地表や海底の地形情報と、地震波反射法による反射断面の情報と、地震波屈折法データの情報と、地震波広角反射データの情報(総合化測定データ群と呼ぶ)に総合的に解釈を加え、地震波速度の水平方向及び深さ方向の分布を示す地震波速度構造を作成し、次に密度と地震波速度の間の最適関係式を用いて前記地震波速度構造から密度構造モデルを作成し、前記密度構造モデルから算出される計算重力ベクトル値と観測による観測重力ベクトル値の差、および、密度構造モデルから算出される計算重力テンソル値と観測から得られる観測重力テンソル値の差から、観測データと計算値との差の2乗和を最小にする最適な密度構造を推定密度構造として出力する地下密度構造の推定方法に関する。

地下あるいは海底下の密度の異常体分布を推定する方法は従来から種々提案されていた。

まず、「地震波反射法探査」という方法がある。地震波反射法探査は、地震波が地下に垂直に入射したときに、地下の地層境界で反射して地表や海底面に戻ってくる地震波の反射に要する時間や振幅を用いて、地下の地震波速度や密度分布の状態を知る方法である。

この地震波反射法探査は、地下構造が作る「音響インピーダンス」すなわち密度と地震波速度の積が地層境界を挟む境界上下でコントラストがある場合、それを利用し、反射地震波の強さの分布（「地震波の反射強度」、反射面の形）、および、反射面間の平均的な地震波速度あるいは反射面間の地震波速度（「区間地震波速度」分布）を求める。

つぎに、地震波が地下で屈折して伝わる部分を用いて地下構造を求める「地震波屈折法」と、入射波と反射波の作る角度がたとえば３０°以上の広角な「地震波広角反射法」を組み合わせた「広角反射法」がある。

「地震波屈折法」と「広角反射法」は地下の「地震波速度の分布（構造）」を求めるのに使用される。この場合、地震波の震源と受震点間を伝搬するのに要する時間（「走時」）を用いて地下の地震波速度の分布（構造）を求める。これを地震波の「走時逆解析（走時インバージョン或いは走時トモグラフィ）法」。波形を利用する場合を「波形逆解析（地震波波形インバージョン）法」と呼ぶ。

さらに、「重力解析」の方法がある。

「重力解析」は重力異常から地下の密度分布（構造）を求める方法であり、重力の異常を観測し、地下の密度異常体を検知する方法である。

特開２００８−０２０１９１号公報 特開２００４−１８４２６８号公報

しかし、上記「地震波反射法探査」は、地震波の反射強度と地下の密度構造が一対一に対応しないこと、浅部の強反射層によりそこより地下深部に地震波が到達しにくいこと、或いは密度変化が緩やかな場合、などにより、その解析によって地下の密度の分布を求めることは困難な場合があった。

たとえば、地震波速度や密度が地中深くなるにしたがい徐々に早くなったり、徐々に遅くなったりするような地下構造は、地震波の「反射面」がないため、地層の境界を反射断面とする地震波反射法を用いることができない。

図５は、従来の地震波反射法探査による地震波反射断面の例を示す図であるが、図から明らかなように、従来の地震波反射法探査による地震波反射断面は、地層境界が途切れ、地層境界を推定できないことがあった。

また、「初動走時インバージョン」は、地震波の走時からはデータの持つ物理的性質のため、極狭い場所ですなわち「局所的」に観測値を説明できるような解（「局所的な解」）が無数に存在し、唯一の地下地震波速度分布（正しい地下構造を示す「唯一の解」）を求めることができなかった。また、浅くて軽い（密度が小さく、地震波速度が小さい）地層と、深くて重い（密度が大きく、地震波速度が大きい）地層が存在する場合、深くまで伝搬した地震波の方が先に到達し、浅いところを通過する地震波が遅れて到達し、これらを区別することが困難であった。また、軽い地層を通過する地震波の振幅が小さいため、軽い地層に対する地震波速度の分布を求めるのは困難であった。

最後に「重力解析」による方法は、海面や航空機などを用いて観測した重力の値を用い、地下の密度構造を推定する「重力逆解法（重力異常のインバージョン法）」が用いられるが、深くて大きな密度が地表や海底面上で密度の異常値を与えるために、唯一の解を求めることが困難であった。「重力解析」による方法では深さ方向の絶対位置を一意的に求めることは困難であった。

そこで、本願発明が解決しようとする課題は、入手可能な地震探査データと重力データをすべて合理的に満足するような、深さ方向を含めて地下の最適密度構造を求めることができる地下密度構造の推定方法を提供することにある。

本発明による地下密度構造の推定方法は、
探査機器の経度、緯度、高度を含む位置情報と、探査対象域の地表または海底の地形情報と、地震波反射法による反射断面の情報と、地震波屈折法のデータの情報（総合化測定データ群）と、総合化測定データ群全体或いは一部を用い得られる地中の地震波速度の水平方向及び深さ方向の分布を示す地震波速度構造モデルを作成する第一段階と、
密度と地震波速度の間の最適関係式を用いて前記地震波速度構造から初期密度構造モデルを作成する第二段階と、
前記初期密度構造モデルから算出される計算重力ベクトル値と、観測による観測重力ベクトル値とを比較し、それらの差である重力異常分布を最小化するように修正密度構造モデルを作成する第三段階と、
前記修正密度構造モデルから算出される計算重力テンソル値と観測から得られる観測重力テンソル値とを比較し、それらの差が最小化するように最適密度構造モデルを作成する第四段階と、
前記最適密度構造モデルを地震波速度構造モデルに変換し、前記第一段階から第四段階を繰り返し、所定の収束条件(総合化測定データ群、と観測重力ベクトルデータ、観測重力テンソルデータなどを所定の誤差範囲で満たす条件)を満たす最適密度構造を推定密度構造として出力する第五段階と、を有することを特徴とする。

本発明は、探査機器の経度、緯度、高度を含む位置情報と、探査対象域の地表または海底の地形情報と、地震波反射法による反射断面の情報と、地震波屈折法のデータの情報と、地震波広角反射のデータの情報(総合化測定データ群)から、地中の地震波速度の水平方向及び深さ方向の分布を示す地震波速度構造を作成する（第一段階）。

次に、密度と地震波速度の間の最適関係式を用いて前記地震波速度構造から初期密度構造モデルを作成する（第二段階）。

その後に、前記初期密度構造モデルから算出される重力のベクトルとテンソルと、観測による重力のベクトルとテンソルを比較し、その差が最小となるように、密度構造モデルに修正を加えるようにする。

このように本発明は、地震波速度構造と重力による密度構造モデルを一体として評価している。

一般に、地震波速度構造は重力データより深部の情報と地層の境界を捉える特性を有しており、一方「重力解析」は密度異常体を深さｘ密度として捉える特性を有している。深い大きな密度体と浅いやや小さな密度体が同じ重力値を与える。本発明では、地震波による地震波速度構造を密度構造モデルに転換し、該密度構造モデルを重力ベクトルと重力テンソルによって評価することにより、両計測データの相補的特徴を生かし、測定データ群全体を合理的かつ最適に説明可能な地下密度構造モデルを求めることができる。

本発明の一実施形態による地下密度構造の推定方法のフローチャート。 固有値解析において、らプラ思案固有関数の固有値のパターン例を示す図。 重力ベクトルと重力テンソルによる重力異常の観測値の面的分布を示す図（右図：一般的なフリーエアー重力異常分布と左図：重力テンソルのGzz成分の分布）。 重力テンソル観測値により地下の密度分布を推定した結果を示す図。 従来の地震波反射法探査による反射法地震波断面図。

以下に添付図面を参照して、本発明の実施形態を説明する。

地下や海底下に、鉱物資源、石油、天然ガス、メタンハイドレート、地下水脈、断層帯、地下埋設構造が存在すると、それらを構成する物質の密度が周辺の物質と異なることにより、密度差の分布を生じる。地下の密度3次元分布によって生じる重力ベクトル分布と重力テンソル分布と地震波の総合化測定データ群を用いてすべての測定データを最適に説明できる様な地下密度構造を推定することにより、上記鉱物資源、石油等の存在場所と埋蔵量を推定できる可能性がある。

図１は、本発明の一実施形態による地下密度構造の推定方法のフローチャートを示している。

図１に示すように、本実施形態による地下密度構造の推定方法は、最初に探査機器の経度、緯度、高度を含む位置情報と、探査対象域の地表または海底の地形情報（図１中の「標高（深度）詳細データ」）と、地震波反射法による反射断面の情報（図１中の「地震波反射法探査」）と、地震波屈折法のデータの情報と、地震波広角反射のデータの情報（図１中の「地震波屈折法探査」）を含む総合化測定データ群を作成し（ステップ１００）、地中の地震波速度の水平方向及び深さ方向の分布を示す地震波速度構造モデルを作成する（ステップ１１０）第一段階を有する。

次に、本実施形態による地下密度構造の推定方法は、密度と地震波速度の間の最適関係式を用いて前記地震波速度構造から初期密度構造モデルを作成する第二段階（ステップ１２０）を有する。

ここで、密度構造は、密度と地震波速度の間にある関係式、たとえばFrancis Birchの実験式などを用い、地震波速度（縦波速度、横波速度）から密度構造に変換する。

また、密度構造に変換するに際しては、対象域情報である位置情報、地形の凹凸と総合測定データ群から得られる、地下密度構造と堆積層の厚さ変化、密度推測値、岩石の密度の取り得べき範囲、を考慮する。

初期密度構造モデルは総合測定データの解釈により真の地下密度構造モデルに近いモデルから解析を始めることにより、局所的解に陥らないようにすることができる。

次に、初期密度構造モデルから算出される計算重力ベクトル値と、観測による観測重力ベクトル値とを比較し、それらの差である重力異常分布を最小化するように修正密度構造モデル（図１中の「地下密度構造修正モデル１」）を作成する第三段階（ステップ１３０−１６０）を有する。

ここで、初期密度構造モデルから算出される計算重力ベクトル値と、観測による観測重力ベクトル値とを比較し、重力異常分布を最小化するように修正密度構造モデルを作成するには、２つの方法がある。

一つの方法（方法１）は線形逆問題（線形インバージョン）解法である（ステップ１４０）。

ここで「線形」とは仮定した地下の構造から推定される重力の理論値が地下の密度分布と一次式の関係にあることをいう。すなわち、
重力異常∝地下密度分布１（ρ_０）
と仮定し、重力異常（観測値と推定値の差）を満たす補正すべき地下密度分布（Δρ）を求める。

修正地下密度分布２（ρ_ｍ）＝地下密度分布１（ρ_０）＋補正地下密度分布（Δρ）
により、新たな修正地下密度分布を求める。

この解析では、重力ベクトルの３成分の観測値がある場合は、これらの観測値を同時に満たす密度分布（構造）を同時に求めることができる（重力の同時解析（重力のジョイントインバージョン法）という）。

なお、上記線形インバージョンにおいて、上下方向、水平方向に近接するモデル間の密度の値に関連があるという仮定をおき、モデルの未知数を減少させる。これを「超（ハイパ）パラメータ」の導入と呼ぶ。導入した密度構造モデルの未知数の数に対し十分な観測値がないときは、ハイパパラメータの導入により、解くべき未知数の数を減らすことができる。

初期密度構造モデルから算出される計算重力ベクトル値と、観測による観測重力ベクトル値とを比較して修正密度構造モデルを作成する２つ目の方法（方法２）は、固有関数による解法である。

固有関数による解法は、２次元あるいは３次元の密度分布を求めるため、地下密度分布を「大局的な変化パターン」から「局所的な変化パターン」までパターン分けする方法において、特殊な固有関数、たとえばラプラシアン固有関数を用いる。

図２にプラシアン固有関数の固有値１から６までのパターンを示す。図２において、♯１は大まかな構造、♯６は細かい構造を示し、固有値が大きいものから小さいものへ並べ、三次元重力ベクトルをラプラシアン固有関数で近似する。

［観測値−現在の密度構造モデルから推定される推定値］^２ ≦ ε
である値εを所定の値、たとえば０．００１ｍ／ｓ^２以下にする。

上記方法１，２のいずれかにより、修正密度構造モデル（図１中の「地下密度構造修正モデル１」）を作成する。

次に、本実施形態は、修正密度構造モデル（図１中の「地下密度構造修正モデル１」）から算出される計算重力テンソル値と観測から得られる観測重力テンソル値とを比較し、それらの差が最小化するように最適密度構造モデルを作成する第四段階（ステップ１７０−２００）を有する。

第三段階と同様に、この第四段階でも２つの方法がある。

一つの方法（方法１）は逆解析の解法である（ステップ１８０）。

上記ステップ１６０で算出された修正密度構造モデル（図１中の「地下密度構造修正モデル１」）から、ステップ１４０で行った方法と同様に線形インバージョンを行う。ただし、重力ベクトルの代わりに、重力ベクトルの３成分を方向微分してできる重力テンソルを使用する。

地下密度構造修正モデル１から算出される重力テンソルと、観測から得られる重力テンソルを満足する最適密度構造モデル（図１中の「地下密度構造修正モデル２」）を求める（ステップ２００）。

修正密度構造モデル（図１中の「地下密度構造修正モデル１」）から算出される計算重力テンソル値と観測から得られる観測重力テンソル値とを比較して最適密度構造モデルを作成する２つ目の方法（方法２）は、固有関数による解法である（ステップ１９０）。

この固有関数による解法は、重力ベクトルの場合と同様であるが、重力ベクトルの代わり重力テンソルを使用する。

固有関数による解法は、２次元あるいは３次元の密度分布を求めるため、地下密度分布を「大局的な変化パターン」から「局所的な変化パターン」までパターン分けする方法において、たとえばラプラシアン固有関数を用いて、三次元重力テンソルをラプラシアン固有関数で近似する。

上記最適密度構造モデル（図１中の「地下密度構造修正モデル２」）が求められると、修正密度構造モデル（図１中の「地下密度構造修正モデル１」）を「改良地下密度構造修正モデル２」に改良する（ステップ２１０）。

次に、上記「改良地下密度構造修正モデル２」を密度から地震波速度へ変換し（ステップ２２０）、また地震波速度分布から、理論地震波屈折波、理論地震波広角反射波、理論地震波反射断面を作成する。以降第一段階から第四段階を繰り返し、所定の収束条件を満たす最適密度構造（図１中の「最終地下密度構造」）を推定密度構造として出力する（第五段階、ステップ２３０）。なお、所定の収束条件とは、総合化測定データ群と観測重力ベクトルデータ、観測重力テンソルデータなどを所定の誤差範囲で満たす条件をいう。

図３は、重力ベクトルと重力テンソルによる重力異常の観測値の面的分布を示す図である。

図に示すように、重力ベクトルに加えて重力テンソルによる解析を行うことにより、重力異常の観測値の面的分布がより詳細に分かることが明らかである。

図４は、重力テンソル観測値により地下の密度分布を推定した結果を示す図である。

図に示すように、重力テンソルにより地下の密度分布を推定した結果、深さ方向に広がった密度構造が得られることが分かる。

Claims (1)

1. 探査機器の経度、緯度、高度を含む位置情報と、探査対象域の地表または海底の地形情報と、地震波反射法による反射断面の情報と、地震波屈折法のデータの情報と、地震波広角反射のデータの情報からなる総合化計測データ群から、地中の地震波速度の水平方向及び深さ方向の分布を示す地震波速度構造モデルを作成する第一段階と、
   密度と地震波速度の間の最適関係式を用いて前記地震波速度構造から初期密度構造モデルを作成する第二段階と、
   前記初期密度構造モデルから算出される計算重力ベクトル値と、観測による観測重力ベクトル値とを比較し、それらの差である重力異常分布を最小化するように修正密度構造モデルを作成する第三段階と、
   前記修正密度構造モデルから算出される計算重力テンソル値と観測から得られる観測重力テンソル値とを比較し、それらの差が最小化するように最適密度構造モデルを作成する第四段階と、
   前記最適密度構造モデルを地震波速度構造モデルに変換し、前記第一段階から第四段階を繰り返し、また地震波速度分布から、理論地震波屈折波、理論地震波広角反射波、理論地震波反射断面を作成し、所定の収束条件を満たす最適密度構造を推定密度構造として出力する第五段階と、を有する地下密度構造の推定方法。