JP4065946B2 - 多成分ae波形の初動検出方法 - Google Patents
多成分ae波形の初動検出方法 Download PDFInfo
- Publication number
- JP4065946B2 JP4065946B2 JP2003077595A JP2003077595A JP4065946B2 JP 4065946 B2 JP4065946 B2 JP 4065946B2 JP 2003077595 A JP2003077595 A JP 2003077595A JP 2003077595 A JP2003077595 A JP 2003077595A JP 4065946 B2 JP4065946 B2 JP 4065946B2
- Authority
- JP
- Japan
- Prior art keywords
- point
- wave
- waveform
- initial motion
- reading
- 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.)
- Expired - Lifetime
Links
Images
Classifications
-
- Y02C10/14—
Landscapes
- Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)
- Geophysics And Detection Of Objects (AREA)
Description
【発明の属する技術分野】
本発明は、多成分AE(アコースティックエミッション)波形の初動検出方法に関するものであり、微小地震の波形、また、例えば、地熱貯留層水圧破砕時などに発生する多成分AE波形の初動を効率的に検出するための方法に関するものである。
【0002】
【従来の技術】
地震観測を中心に、観測波形から地震波のP波およびS波の自動読み取り法が開発されている。しかしながら、より複雑な岩盤破砕に伴う多成分AE波形を対象にしたものは少なく、アルゴリズムも人間の目視作業の流れに則ったものは存在しなかった。
【0003】
地震観測の技術分野において、これまでに、発表されている論文としては、例えば、次に説明する非特許文献1、非特許文献2のようなものがある。
【0004】
非特許文献1において、発表されている内容は、現在も広く地震観測の分野で用いられている局所定常AR(Autoregressive)モデルを用いた手法に関するものであり、この論文においては、何種類かのモデルの用い方が検討されて、P波およびS波は共に基本的に一つのアルゴリズムに基づいて自動読み取りを行っているものが示されている。この手法では、S波について振動状態の利用は無く、岩盤破砕時の多成分AE波形への適用は行われていない。
【0005】
また、非特許文献2に示されるように、岩盤破砕に伴う多成分(3成分)AE波形を対象にして1観測点でのAE源位置標定法・3軸ホドグラム法を前提にした自動読み取り法が検討されているが、ここにおいても、基本的に、S波振動状態に関する情報は利用せずに、振幅情報にのみ立脚した信号処理を行う方法が取られている。従来においても、地震波のP波およびS波の両方を対象に検討が行われたが、現在では、このような論文の内容が、多成分AE波形の現場解析で利用されることは無い。
【非特許文献1】
横田 崇・周 勝奎・溝上 恵・中村 功(1981):地震波データの自動検測方式とオンライン処理システムにおける稼動実験, Bull. Earthq. Res. Inst., 55, 449-484.
【非特許文献2】
Nagano, K., Niitsuma, H., and Chubachi, N., 1989, Automatic algorithm for triaxial hodogram source location in downhole acoustic emission measurement. Geophysics, 54, No. 4, p. 508-513.
【0006】
【発明が解決しようとする課題】
地熱開発などにおいて利用される岩盤破砕時に発生するアコースティックエミッション(AE、ここでは微小地震とも呼ばれるもの)波形は、地震波に比べて波形が複雑であることが多いうえ、大量に連続的に発生する。微弱なAE波形の観測手法としては、少数観測点においてAEセンサとして坑井内多成分弾性波検出器を設置して観測する方法が想定されるが、その場合には、より複雑なS波の読み取りも一般に必要である。
【0007】
また、従来において、地震観測などで開発された単純あるいは単一の手法に依存するP波およびS波の初動の自動読み取り法は、これまで、実際のAE波形の読み取り方法として適用されることはほとんど無かった。このため、結果的に、目視による波形観察と3次元粒子運動軌跡観察により、手作業でP波およびS波の初動読み取りが行われることが、一般的であった。
【0008】
例えば、発明者らも参加した2000年ソルツ高温岩体テストサイト(フランス)で行われた水圧破砕時のAE(微小地震)観測では、現場でのAE発生位置の速報を2人交代24時間体制の目視作業で行うことが試みられた。しかしながら、AE発生総数が膨大であるため、多大なる労力に比して、全記録数の僅か2%程度しか解析結果を得ることができなかった。
【0009】
このように、これまでのAE波形の読み取り方法において、S波の読み取り方法としては、P波読み取りと同様に振幅変化や統計モデルを利用した手法が試みられているが、S波の性状とは無関係である手法に基づいており、誤差の原因となっていた。特に、S波の振動方向に関する情報は、目視作業では頻繁に用いられていたが、自動読み取りには適用されてこなかった。
【0010】
また、これまでのほとんどのAE波形の自動読み取りの場合、計算処理の時間短縮を重視して、ひとつの信号処理法が支配的であることが多く、複雑な多成分AE波形を読みとる場合に対する適用性は高くは無かった。
【0011】
一方、近年の超解像AEマッピング法の開発や、そのAEを用いる解析手法の応用分野の広がり、適用事例の広がりに伴い、対象分野においては初動読み取り誤差を極力少なくすることが求められている。初動読み取り結果は、後々の解析手法に影響を与えることから、数ポイント程度のズレを問題にすることも頻繁に起こってきている。このため、地震観測で広く実用的に用いられている局所定常ARモデルを用いる手法も、波の不連続点に読み取りポイントが必ずしも一致しないという面が問題になっている。
【0012】
本発明は、上記のような問題点を解決するためになさたれものであり、本発明の目的は、多成分AE波形の初動検出方法を提供することにある。具体的に、本発明の目的は、第1に、複雑な多成分AE波形を対象にした目視作業と遜色無い実用的な自動読み取り手法を提供することができ、また、第2に、S波性状を反映した初動の読み取り方法を提供することができ、第3に、より高度な解析手法に適用可能である不連続点を初動到来点として正確に読み取る手法を提供することができる多成分AE波形の初動検出方法を提供することにある。
【0013】
【課題を解決するための手段】
上記のような目的を達成するため、本発明による多成分AE波形の初動検出方法は、基本原理として、波形振幅変化や局所定常ARモデルを用いる手法に加えて、S波振動方向の評価と波の不連続点検出を行い、複数の要素を組み合わせて目視作業による解析の流れに則った方法で、高精度に多成分AE波形のP波およびS波の初動の自動検出を行う。
【0014】
具体的には、本発明による多成分AE波形の初動検出方法は、複数の観測点における多成分AEセンサにより検出したAE信号を取得し、取得したAE信号の信号波形の波形記録を行うと共に、検出した信号波形における信号対雑音比をチャンネルごとに判定して、信号対雑音比が高いチャンネルを選択し、選択したチャンネルからAE信号の信号波形のエネルギー変化の解析により、多成分AE波形から初動を検出する方法であって、P波とS波の検出対象時間範囲をそれぞれ設定するステップと、前記検出対象時間範囲においてP波初動検出処理を行うステップと、前記検出対象時間範囲においてS波初動検出処理を行うステップと、前記P波初動検出処理およびS波初動検出処理に基づき、P波初動読み取りポイントの観測点の間の相互関係から矛盾がないようにS波読み取りポイントの不整合除去を行うステップとを備える。
【0015】
また、この多成分AE波形の初動検出方法において、P波初動検出処理を行うステップは、1つの観測点について、検出対象時間範囲における所定の基準点の周辺において、波形のエネルギー変化が顕著である点と局所定常ARモデルにより決定される変化点とを算出して1次P波到来ポイント候補点とする第1ステップと、前記波形のエネルギー変化が顕著である点および前記変化点の前後で波形エネルギーの変化率を計算し、変化が大なものを2次P波到来ポイント候補点とする第2ステップと、前記2次P波到来ポイント候補点の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出し、前記2次P波到来ポイント候補点に最も近い不連続点を、最終P波到来ポイント候補点とする第3ステップと、前記最終P波到来ポイント候補点の前後で有意なエネルギー変化が確認された場合に、当該最終P波到来ポイント候補点をP波読み取りポイントとする第4ステップと、前記第1ステップ〜第4ステップの処理を全観測点に対して行い、各観測点におけるP波読み取りポイントを得る第5ステップとを備えることを特徴とするものである。
【0016】
また、この多成分AE波形の初動検出方法において、S波初動検出処理を行うステップは、1つの観測点について、検出対象時間範囲における所定の基準点の周辺において、前記P波読み取りポイント付近で、時間−周波数領域での振動方向解析を行い、P波到来方向を抽出する第6ステップと、前記抽出したP波到来方向がP軸に一致するように、AE波形の3成分波形をP−SV−SH系に座標変換する第7ステップと、波形のエネルギー変化が顕著である点、局所定常ARモデルにより決定される変化点、さらに時間−周波数領域での振動方向解析によりP波到来方向と振動方向の直交度が最大になる直交度最大点を、1次S波到来ポイント候補点として算出すると共に、S波エネルギーの最大点を算出する第8ステップと、直交度最大点との差が所定のしきい値以下のとき、波形のエネルギー変化が顕著である点、局所定常ARモデルにより決定される変化点、直交度最大点、S波エネルギーの最大点の前後での波形エネルギーの変化率をそれぞれ計算し、変化が最大なものを2次S波到来ポイント候補点とする第9ステップと、2次S波到来ポイント候補点の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出すると共に、2次S波到来ポイント候補点に最も近い不連続点を、最終S波到来ポイント候補点とする第10ステップと、最終S波到来ポイント候補点の前後で有意なエネルギー変化が確認された場合に、当該最終S波到来ポイント候補点をS波読み取りポイントとする第11ステップと、前記第6ステップ〜第11ステップの処理を全観測点に対して行い、各観測点におけるS波読み取りポイントを得る第12ステップとを備えることを特徴とするものである。
【0017】
また、この多成分AE波形の初動検出方法において、S波読み取りポイントの不整合除去を行うステップは、各観測点におけるS波読み取りポイントの観測点間の大小関係が、各観測点におけるP波読み取りポイントにおける観測点間の大小関係と矛盾がないことを条件に、S波読み取りポイントから読み取り異常のものを削除することを特徴とするものである。
【0018】
【発明の実施の形態】
以下、本発明を実施する場合の形態について、図面を参照して具体的に説明する。図5は、本発明の多成分AE波形の初動検出方法を一形態で実施する場合のシステムの構成の概略を説明する図である。図において、10はAE信号伝送のためのネットワーク、11〜13は複数の各観測点の多成分AEセンサ、20はデータ処理装置、21はAE信号波形記録装置、22は波形解析処理装置である。
【0019】
複数の各観測点に設けられて、AE信号を検出する多成分AEセンサ11〜13は、AE信号伝送のためのネットワーク10を介して、AE信号波形に対するデータ処理を行うデータ処理装置20に接続されている。データ処理装置20には、データ処理を行ったデータを記憶する共に、多成分AEセンサ11〜13からのAE信号波形を記録するためのAE信号波形記録装置21、AE波形に対する波形解析などの特別なデータ処理を行うための波形解析処理装置22が接続されている。なお、この波形解析処理装置22が行うデータ処理を、データ処理装置20において行われるようにシステムを構成しても良い。また、図示されないが、データ処理装置20には、AE波形解析のデータ処理のためのコマンドを入力する入力装置、データ処理した結果を所望の表示形態で画面表示するための表示装置や、AE波形などを印刷文書データとして出力するための印刷装置が接続されている。
【0020】
多成分AEセンサ11〜13は、地熱開発、石油開発、炭酸ガスの地中貯留などでの地下計測を行う場合には、AE信号を3成分以上の多成分で検出するための坑井内弾性波検出器が用いられる。また、多成分AEセンサ11〜13からの信号は直接にデータ処理装置20に入力されるようにしても良い。
【0021】
本発明による多成分AE波形の初動検出方法の理解を容易とするため、目視作業によるAE波形の初動の判定方法について、作業工程順に説明する。
【0022】
《目視作業での流れ》
工程1:観測波形全体を眺め、大雑把にP波およびS波の受信されている範囲を定める。
工程2:大雑把なP波付近をチャンネル毎に拡大表示させて、エネルギーの変化、受信状態(周波数や信号分散)の変化に注目し、さらに詳しくP波到来候補点を決める。
工程3:P波到来候補点付近を極端に拡大し、波形記録に不連続が生じているポイントにP波到来候補点を一致させて、P波初動読み取り点とする。
工程4:P波読み取り点付近の3次元粒子運動軌跡を表示し、P波到来方向を表示する。
工程5:大雑把に定めたS波受信範囲において、エネルギー変化、受信状態(周波数や信号分散)の変化からS波到来候補点を定める。
工程6:S波到来候補点付近で3次元粒子運動軌跡を観察し、P波到来方向との直交度が急増するポイントを見つける。
工程7:直交度の急増ポイントの前で、波形記録に不連続が生じているポイントをみつけて、S波初動の読み取り点とする。
工程8:工程2〜工程7を使用する観測点のすべてについて行い、最後に観測点間相互関係についてP波初動とS波初動で矛盾の生じているものをP波を優先して削除する。
【0023】
このような目視作業によるAE波形の初動の判定方法に対応して、本発明による多成分AE波形の初動検出方法においては、その各作業工程を自動化して、機械的に処理するため、前述したようなシステム構成を用いて、データ処理を行って、多成分AE波形の初動検出を行う。次に、自動化して機械的に処理するための多成分AE波形の初動検出方法について説明する。この場合についても、理解を容易とするため、上述した目視作業によるAE波形の初動の判定方法に対応して説明する。
【0024】
《自動初動検出の流れ》
<前処理>
工程1:波形記録(例:3成分×3観測点+1成分×2観測点=11チャンネルなど)のうち、チャンネル毎にSN比を評価し、最も良好のものを選択する。この波形のエネルギー変化の解析により、大雑把なP波とS波の検出対象時間範囲(P波初動解析のための基準点をbaseP、S波初動解析のための基準点をbaseSとする)をそれぞれ設定する。
【0025】
<P波検出部分の処理>
工程2:1つの観測点Aについて、基準点basePの周辺において、波形のエネルギー変化が顕著である点(Pc1a)、局所定常ARモデルにより決定される変化点(Pc1b)を1次P波到来ポイント候補点として算出する。
工程3:点(Pc1a)および変化点(Pc1b)の前後での波形エネルギーの変化率を計算し、変化が大なものを2次P波到来ポイント候補点(Pc2)とする。
工程4:2次P波到来ポイント候補点(Pc2)の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出する(1)。2次P波到来ポイント候補点(Pc2)に最も近い不連続点を、最終P波到来ポイント候補点(Pc)とする。
工程5:最終P波到来ポイント候補点(Pc)の前後において有意なエネルギー変化が確認されたなら、最終P波到来ポイント候補点(Pc)をP波読み取りポイント(Pc-A)として採用する。
工程6:工程2〜工程5を全観測点に適用して処理を行い、各観測点A〜EにおけるP波読み取りポイント(Pc-A)〜(Pc-E)を得る。
【0026】
<S波検出部分の処理>
工程7:1つの観測点Aについて、基準点baseSの周辺において、P波の読み取りポイント(Pc-A)付近で、時間−周波数領域での振動方向解析(2)を行って、P波到来方向(θA)を抽出する。
工程8:P波到来方向(θA)の方向がP軸に一致するように、3成分波形(X−Y−Z)をP−SV−SH系に座標変換する。
工程9:波形のエネルギー変化が顕著である点(Sc1a)、局所定常ARモデルにより決定される変化点(Sc1b)、時間−周波数領域での振動方向解析(2)によりP波到来方向(θA)の方向と振動方向の直交度が最大になる直交度最大点(Sc1c)を1次S波到来ポイント候補点として算出すると共に、S波エネルギーの最大点を求めて(Sc-max)とする。
工程10:直交度最大点(Sc1c)との差が所定のしきい置T1以下のとき、先に算出した波形のエネルギー変化が顕著である点(Sc1a)、変化点(Sc1b)、直交度最大点(Sc1c)、S波エネルギーの最大点(Sc-max)の前後での波形エネルギーの変化率をそれぞれ計算し、変化が最大なものを2次S波到来ポイント候補点(Sc2)とする。
工程11:2次S波到来ポイント候補点(Sc2)の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出する。2次S波到来ポイント候補点(Sc2)に最も近い不連続点を、最終S波到来ポイント候補点(Sc)とする。
工程12:最終S波到来ポイント候補点(Sc)の前後で有意なエネルギー変化が確認されたなら、最終S波到来ポイント候補点(Sc)をS波読み取りポイント(Sc-A)として採用する。
工程13:工程7〜工程12を全観測点に適用し、各観測点A〜EにおけるS波読み取りポイント(Sc-A)〜(Sc-E)を得る。
工程14:各観測点A〜EにおけるS波読み取りポイント(Sc-A)〜(Sc-E)の観測点間の大小関係が各観測点A〜EにおけるP波読み取りポイント(Pc-A)〜(Pc-E)における観測点間の大小関係と矛盾がないことを条件に、各観測点A〜EにおけるS波読み取りポイント(Sc-A)〜(Sc-E)から読み取り異常のものを削除する。
【0027】
<最終表示の処理>
工程15:以上の操作を全対象波形データについて実施し、PおよびS波の初動の自動読み取りが完了する。
工程16:読み取りポイントを元にして、AE源位置標定などの解析を行い、所望の表示形態で表示出力する。
【0028】
以上のような各工程のデータ処理を行うことにより、目視作業による詳細解析の流れを完全に踏襲しながら、自動で多成分AE波形からのP波およびS波の初動の検出が可能となる。
【0029】
これらの処理においては、一部において、グラフィックユーザインタフェース処理を介しての操作者によるコマンド入力を受け付けながらデータ処理装置による処理によって実行されるが、各観測点からAE信号を取り込みながらAE波形に対するデータ処理が自動化されるようにシステム構成されても良い。
【0030】
図1は、本発明によるAE波形の自動初動検出処理の流れを説明するフローチャートである。この処理は、前述したシステム構成におけるデータ処理装置において実行される。この処理では、前処理として、複数の観測点における多成分AEセンサにより検出したAE信号を取得し、取得したAE信号の信号波形の波形記録を行うと共に、検出した信号波形における信号対雑音比をチャンネルごとに判定して、信号対雑音比が高いチャンネルを選択する。そして、選択したチャンネルからAE信号の信号波形のエネルギー変化の解析により、P波とS波の検出対象時間範囲をそれぞれ設定する。これらの処理は、各観測点毎のAE信号について行う。
【0031】
このため、初動検出処理を自動化して処理する場合は、処理の流れの制御として、各観測点毎のAE信号について行うため、未処理イベントがあるか否かを判定する(ステップ100)。そして、未処理イベントがある場合に、検出対象範囲決定の処理を行う(ステップ101)。次に、P波未検出の観測点があるか否かを判定し、前記検出対象時間範囲においてP波未検出の観測点がある場合にP波初動検出処理を行い(ステップ102,103)、次に、S波未検出の観測点があるか否かを判定し、前記検出対象時間範囲においてS波未検出の観測点がある場合にS波初動検出処理を行う(ステップ104,105)。そして、P波初動検出処理およびS波初動検出処理に基づき、P波初動読み取りポイントの観測点の間の相互関係から矛盾がないようにS波読み取りポイントの不整合除去を行う(ステップ106)、最終表示を行って、次の観測点における処理を行うために、ステップ100に戻り、処理を続行する。
【0032】
なお、このS波読み取りポイントの不整合除去を行う場合には、各観測点におけるS波読み取りポイントの観測点間の大小関係が、各観測点におけるP波読み取りポイントにおける観測点間の大小関係と矛盾がないことを条件に、S波読み取りポイントから読み取り異常のものを削除する。
【0033】
図2は、検出対象範囲決定処理を説明するフローチャートである。検出対象範囲の決定の処理では、各観測点から伝送されて記録されているAE波形の記録から基本波形記録の選択を行い、例えば、観測している波形記録のうち最大SN比の波形トレースを選択し(ステップ201)、波形エネルギー変化解析を行って(ステップ202)、検出対象範囲を決定する(ステップ203)。この検出対象範囲の決定では、基準点として、P波初動解析のための基準点baseP、S波初動解析のための基準点baseSを決定して、また、検出範囲幅を設定する(ステップ203)。
【0034】
図3は、P波初動検出処理を説明するフローチャートである。多成分AE波形の初動検出方法において、まず、P波初動検出処理が、各観測点について行われる。図3に示すように、未処理観測点があるか否かを判定し、未処理観測点がある場合に処理を行う(ステップ300)。1つの観測点について、波形エネルギー変化解析の処理を行い、検出対象時間範囲における所定の基準点の周辺において、波形のエネルギー変化が顕著である点と局所定常ARモデルにより決定される変化点とを算出して1次P波到来ポイント候補点とする(ステップ301)。次に、局所定常ARモデル解析を行い(ステップ302)、前記波形のエネルギー変化が顕著である点および前記変化点での波形エネルギーの変化率を計算し、変化が大なものを2次P波到来ポイント候補点とする。そして、2次候補点の選択の処理(ステップ303)を行い、2次P波到来ポイント候補点の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出し、前記2次P波到来ポイント候補点に最も近い不連続点を、最終P波到来ポイント候補点とする(ステップ304)。最終P波到来ポイント候補点の前後で有意なエネルギー変化が確認された場合に、当該最終P波到来ポイント候補点をP波読み取りポイントとする(ステップ305)。そして、これらの処理を全観測点に対して行い、各観測点におけるP波読み取りポイントを得る。
【0035】
図4は、S波初動検出処理を説明するフローチャートである。S波初動検出処理についても各観測点について行われる。図4に示すように、未処理観測点があるか否かを判定し、未処理観測点がある場合に処理を行う(ステップ400)。1つの観測点について検出対象時間範囲における所定の基準点の周辺において、前記P波読み取りポイント付近で、時間−周波数領域での振動方向解析を行い、P波到来方向を抽出する(ステップ401)。次に、抽出したP波到来方向がP軸に一致するように、AE波形の3成分波形をP−SV−SH系に座標変換する(ステップ402)。次に、波形エネルギー変化解析の処理、局所定常ARモデル解析の処理を行い(ステップ403、404)、波形のエネルギー変化が顕著である点、局所定常ARモデルにより決定される変化点、さらに時間−周波数領域での振動方向解析によりP波到来方向と振動方向の直交度が最大になる直交度最大点を、1次S波到来ポイント候補点として算出すると共に、S波エネルギーの最大点を算出する(ステップ405、406)。そして、直交度最大点との差が所定のしきい値以下のとき、波形のエネルギー変化が顕著である点、局所定常ARモデルにより決定される変化点、直交度最大点、S波エネルギーの最大点の前後での波形エネルギーの変化率をそれぞれ計算し、変化が最大なものを2次S波到来ポイント候補点とする(ステップ407、408)。次に、2次S波到来ポイント候補点の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出すると共に、2次S波到来ポイント候補点に最も近い不連続点を、最終S波到来ポイント候補点として(ステップ409)、最終S波到来ポイント候補点の前後で有意なエネルギー変化が確認された場合に、当該最終S波到来ポイント候補点をS波読み取りポイントとする(ステップ410)。これらの処理を全観測点に対して行い、各観測点におけるS波読み取りポイントを得る。
【0036】
次に、前述した各処理において用いられるウェーブレット変換による不連続点検出について説明する。ウェーブレット変換は積分変換の一種であり、連続ウェーブレット変換は
【数1】
で与えられる。a,bは拡大倍率と時刻に関する位置を表すパラメータでスケールとシフトと呼ばれる。h(t)はアナライジングウェーブレットと呼ばれる関数で、*は複素共益を示す。アナライジングウェーブレットとして、時間と周波数に関する不確定性が最小であるガボール関数を用いる。ガボール関数は、
【数2】
で与えられる。ウェーブレット変換を用いて波の不連続点を検出するため、位相指標を用いる。これはウェーブレット係数の位相を積分したもので、以下の関数の実部で表される。
【数3】
w(a)は正規化のための関数であり、次式で表される。
【数4】
ウェーブレット変換の特徴の一つとしては、「シフト−スケール(時間−周波数)平面上でウェーブレット係数の位相を表示したとき等位相線が不連続点の存在する時刻に向かって収束する」ということがあり、このため、不連続点近傍で位相指標は高い値を示し、それ以外では、位相が一定でないために小さい値を示す。しきい値を与えて位相指標を評価することにより、波の不連続点の検出が可能である。
【0037】
次に、時間−周波数領域での振動方向解析について説明する。連続ウェーブレット変換の定義式において、以下のようなMayer’s waveletと呼ばれるアナライジングウェーブレットを用いることにより、正規直交基底を有するウェーブレット変換を行う。メイヤーズウェーブレット(Mayer's wavelet)は、次のように定義される。
【数5】
【0038】
正規直交基底を有するため、3成分信号にウェーブレット変換を適用した場合、成分間のエネルギーの相互関係は完全に保持されるので正確な3次元粒子運動の解析に応用可能である。一方、ウェーブレット変換は信号の時間−周波数表現の一種であるとみなせるので、これを用いてパワースペクトルおよびクロススペクトルを算出し、時間−周波数領域での分散共分散行列すなわちスペクトル行列が、以下のように定義できる。
【数6】
この行列は、ウェーブレット変換で言うところのスケールとシフトの関数であるが、これらはそれぞれ周波数と時刻に相当するため、完全に時間−周波数領域の情報として取り扱うことができ、振動方向の時間−周波数領域での解析は、行列の第一固有値と第一固有ベクトルを用いて行う。固有値、固有ベクトルは、通常の共分散行列と同様に導出できる。
【0039】
P波到来方向はP波初動時刻付近の第一固有ベクトル方向Dir(Pp,f)から求められる。これは周波数の関数であるが、通常ノイズの周波数を避けるとともに、3次元粒子運動がもっとも直線的になる周波数の値を用いる。3次元粒子運動の直線性は例えば第一主成分寄与率から評価できる。
【0040】
また、S波の受信状態は、P波到来方向Dir(Pp,f)と時々刻々と変化する波の振動方向(第一固有ベクトル方向)Dir(t,f)が成す角Ang(t,f)を算出することにより評価される。ここでは、評価対象の時間範囲内で、Ang(t,f)が最大(<90°)となる時刻tを、S波が最も優勢に受信されている時刻としている。
【0041】
【実施例】
本発明によるAE波形の初動検出方法を評価するため、具体的に、2000年にソルツ高温岩体テストサイト(フランス)で記録された11チャンネル(3成分×3観測点、1成分×2観測点)波形データを用いて、本発明による初動検出方法によるデータ処理についての検討を行った。これについて説明する。
【0042】
図6は、ウェーブレット変換による不連続点検出結果を説明する図である。図6においては、上部側に2つの正弦波(25kHzおよび50kHz)の重ね合わせによるシミュレーション波形を示し、下部側に位相指標による不連続点検出結果を示している。上部側の図において点線部分が不連続点である。下部側の図において色の濃い部分が、推定された不連続を表している。このように、ウェーブレット変換による不連続点検出のシミュレーション例のように、微妙な波の不連続点を抽出することが可能である。
【0043】
図7は、ウェーブレット変換を用いた時間−周波数領域での振動方向解析の例を説明する図である。図7においては、P波到来方向と波が成す角度を時間−周波数領域で表示しており、S波到来付近ではその角度が大きくなっていることが分かる。図7の上部側の図は、ソルツ2000年の観測波形例を示しており、下部側の図において、P波振動方向に対する各時刻の波の振動方向が成す角の時間−周波数領域での変化の例を示している。図中の点線で囲んだ付近がS波の到来に対応しているものである。
【0044】
図8は、P波の自動読み取り結果と目視および従来法(ARモデル)との比較を説明するための図であり、図9は、S波の自動読み取り結果と目視および従来法(ARモデル)との比較を説明するための図である。従来法では、波の不連続点と無関係な、やや後ろの点が読まれているが、本方法では目視と同じ点を検出している。本方法では振動方向の情報を利用しているため、目視に近い点が検出されている。
【0045】
図8および図9に示すように、ソルツ地域の観測波形について、本発明による初動検出方法により、自動処理によってP波およびS波の初動の読み取りを行った例では、S波を含めて、高精度に初動の自動での読み取りが行えていることが分かる。比較のために、一般的な局所定常ARモデルのみを用いた読み取り結果を同時に示している。このように、本発明による初動検出方法によると、従来法に比べて目視作業に近い読み取りが自動で行えている。
【0046】
このように、P波の読み取りが目視と同様に不連続点に一致するようになり、また、S波の目視とのズレも改善される傾向となっている。従来法でも地震観測に対しては十分な精度であったと思われるが、読み取りポイントが波の不連続点に一致するなど、本発明による初動検出方法では、岩盤破砕時のAEの超解像マッピングに求められる、目視作業と遜色が無い正確さが達成できている。
【0047】
図10は、目視によるAE波形のP波およびS波の初動検出値により計算されたAE源分布(左:3472点)を示す図であり、図11は、自動読み取りにより計算されたAE波形のP波およびS波の初動検出値によるAE源分布(右:3322点)を示す図である。図10および図11を比較して参照すると、自動処理においても、地下構造の解釈に繋がるAE源分布の主要構造は、目視による結果とほぼ一致している。
【0048】
初動読み取りの処理効率は、目視作業に比べて大幅に改善され、約10倍の処理速度の改善が見られている。図12に示している表は、作業効率の改善について示している。この表により理解されることは、慎重な目視作業による初動の読み取り結果が、現場でのデータ取得後数ヵ月後に公表されたが、本発明による初動検出方法により、それを類似した正確さで数日間で達成することが可能になることが分かる。また、前述したように、これを元にするAE源分布(図10および図11)についても、目視作業によるものとほとんど変わらない正確さのものが、AEの観測期間内に導出できることが分かる。
【0049】
このように、本発明の初動検出方法によると、岩盤破砕時に観測される3成分AE波形を対象とする複数の信号処理手法を組み合わせた人間による目視作業の流れを模擬した高精度なP波およびS波初動の自動読み取り法が提供され、時間−周波数領域の振動方向解析を用いたP波振動方向との直交性評価によるS波受信状態の評価法に用いることができ、また、ウェーブレット変換による不連続点検出を利用した初動時刻の正確な決定法が提供される。
【0050】
【発明の効果】
本発明の初動検出方法によると、従来法の適用が困難であった岩盤破砕時の多成分AE波形に対して、P波およびS波の自動初動読み取りが適用可能になる。従来法では、S波の正常な読み取りが困難であったり、波の不連続点と無関係な点を読み取ってしまったりしていたが、本発明の方法によって、目視作業と同等の正確さで自動読み取りが可能になる。目視作業に対して10倍程度改善される処理効率についても、正確さを考慮すれば、再読み取り作業が不要になるので、十分有用なレベルである。
【0051】
また、本発明の初動検出方法により、岩盤破砕が発生しているその場での、AE源発生位置推定が実施できるようになり、地熱開発だけでなく、CO2地中貯留などの際の安全監視的な用途にも利用可能である。また、本発明によって、複雑な多成分AE解析技術も現場で直ちに利用できることが想定でき、幅広い分野に適用できる。
【図面の簡単な説明】
【図1】本発明によるAE波形の自動初動検出処理の流れを説明するフローチャートである。
【図2】自動初動検出処理における検出対象範囲決定処理を説明するフローチャートである。
【図3】自動初動検出処理のP波初動検出処理を説明するフローチャートである。
【図4】自動初動検出処理のS波初動検出処理を説明するフローチャートである。
【図5】本発明の多成分AE波形の初動検出方法を一形態で実施する場合のシステムの構成の概略を説明する図である。
【図6】ウェーブレット変換による不連続点検出結果を説明する図である。
【図7】ウェーブレット変換を用いた時間−周波数領域での振動方向解析の例を説明する図である。
【図8】P波の自動読み取り結果と目視および従来法(ARモデル)との比較を説明するための図である。
【図9】S波の自動読み取り結果と目視および従来法(ARモデル)との比較を説明するための図である。
【図10】目視によるAE波形のP波およびS波の初動検出値により計算されたAE源分布(左:3472点)を示す図である。
【図11】自動読み取りにより計算されたAE波形のP波およびS波の初動検出値によるAE源分布(右:3322点)を示す図である。
【図12】自動処理による作業効率を説明する図である。
【符号の説明】
10…AE信号伝送のためのネットワーク、
11〜13…各観測点の多成分AEセンサ、
20…データ処理装置、
21…AE信号波形記録装置、
22…波形解析処理装置。
Claims (3)
- 複数の観測点における多成分AEセンサにより検出したAE信号を取得し、取得したAE信号の信号波形の波形記録を行うと共に、検出した信号波形における信号対雑音比をチャンネルごとに判定して、信号対雑音比が高いチャンネルを選択し、選択したチャンネルからAE信号の信号波形のエネルギー変化の解析により、多成分AE波形から初動を検出する方法であって、
P波とS波の検出対象時間範囲をそれぞれ設定するステップと、
前記検出対象時間範囲においてP波初動検出処理を行うステップと、
前記検出対象時間範囲においてS波初動検出処理を行うステップと、
前記P波初動検出処理およびS波初動検出処理に基づき、P波初動読み取りポイントの観測点の間の相互関係から矛盾がないようにS波読み取りポイントの不整合除去を行うステップと、
を備える多成分AE波形の初動検出方法において、
P波初動検出処理を行うステップは、
1つの観測点について、検出対象時間範囲における所定の基準点の周辺において、波形のエネルギー変化が顕著である点と局所定常ARモデルにより決定される変化点とを算出して1次P波到来ポイント候補点とする第1ステップと、
前記波形のエネルギー変化が顕著である点および前記変化点の前後での波形エネルギーの変化率を計算し、変化が大なものを2次P波到来ポイント候補点とする第2ステップと、
前記2次P波到来ポイント候補点の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出し、前記2次P波到来ポイント候補点に最も近い不連続点を、最終P波到来ポイント候補点とする第3ステップと、
前記最終P波到来ポイント候補点の前後で有意なエネルギー変化が確認された場合に、当該最終P波到来ポイント候補点をP波読み取りポイントとする第4ステップと、
前記第1ステップ〜第4ステップの処理を全観測点に対して行い、各観測点におけるP波読み取りポイントを得る第5ステップと
を備えることを特徴とする多成分AE波形の初動検出方法。 - 前記請求項1に記載の多成分AE波形の初動検出方法において、
S波初動検出処理を行うステップは、
1つの観測点について、検出対象時間範囲における所定の基準点の周辺において、前記P波読み取りポイント付近で、時間−周波数領域での振動方向解析を行い、P波到来方向を抽出する第6ステップと、
前記抽出したP波到来方向がP軸に一致するように、AE波形の3成分波形をP−SV−SH系に座標変換する第7ステップと、
波形のエネルギー変化が顕著である点、局所定常ARモデルにより決定される変化点、さらに時間−周波数領域での振動方向解析によりP波到来方向と振動方向の直交度が最大になる直交度最大点を、1次S波到来ポイント候補点として算出すると共に、S波エネルギーの最大点を算出する第8ステップと、
直交度最大点との差が所定のしきい値以下のとき、波形のエネルギー変化が顕著である点、局所定常ARモデルにより決定される変化点、直交度最大点、S波エネルギーの最大点の前後での波形エネルギーの変化率をそれぞれ計算し、変化が最大なものを2次S波到来ポイント候補点とする第9ステップと、
2次S波到来ポイント候補点の周辺で、ガボール関数によるウェーブレット変換を計算し、位相指標から不連続点を算出すると共に、2次S波到来ポイント候補点に最も近い不連続点を、最終S波到来ポイント候補点とする第10ステップと、
最終S波到来ポイント候補点の前後で有意なエネルギー変化が確認された場合に、当該最終S波到来ポイント候補点をS波読み取りポイントとする第11ステップと、
前記第6ステップ〜第11ステップの処理を全観測点に対して行い、各観測点におけるS波読み取りポイントを得る第12ステップと
を備えることを特徴とする多成分AE波形の初動検出方法。 - 前記請求項2に記載の多成分AE波形の初動検出方法において、
S波読み取りポイントの不整合除去を行うステップは、
各観測点におけるS波読み取りポイントの観測点間の大小関係が、各観測点におけるP波読み取りポイントにおける観測点間の大小関係と矛盾がないことを条件に、S波読み取りポイントから読み取り異常のものを削除する
ことを特徴とする多成分AE波形の初動検出方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2003077595A JP4065946B2 (ja) | 2003-03-20 | 2003-03-20 | 多成分ae波形の初動検出方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2003077595A JP4065946B2 (ja) | 2003-03-20 | 2003-03-20 | 多成分ae波形の初動検出方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
JP2004286527A JP2004286527A (ja) | 2004-10-14 |
JP4065946B2 true JP4065946B2 (ja) | 2008-03-26 |
Family
ID=33292310
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2003077595A Expired - Lifetime JP4065946B2 (ja) | 2003-03-20 | 2003-03-20 | 多成分ae波形の初動検出方法 |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP4065946B2 (ja) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105137485A (zh) * | 2015-08-31 | 2015-12-09 | 中国石油天然气股份有限公司 | 一种地震剖面位图中解释信息的自动拾取转化方法及装置 |
Families Citing this family (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2008046077A (ja) * | 2006-08-21 | 2008-02-28 | Taiheiyo Cement Corp | 粒度測定方法 |
EP2116921B1 (en) * | 2008-05-07 | 2011-09-07 | Tyco Electronics Services GmbH | Method for determining the location of an impact on a surface of an object |
JP6131098B2 (ja) * | 2013-05-16 | 2017-05-17 | 大成建設株式会社 | 弾性波波群解析方法 |
CN103336063B (zh) * | 2013-06-20 | 2016-01-20 | 江苏大学 | 一种声发射信号初至点检测方法 |
KR101570085B1 (ko) | 2013-11-01 | 2015-11-19 | 한국원자력연구원 | P파 도달 시간 검출 장치 및 검출 방법 |
JP6945917B2 (ja) * | 2017-12-22 | 2021-10-06 | 株式会社奥村組 | トンネル切羽前方探査方法 |
CN113237957B (zh) * | 2021-05-31 | 2023-09-29 | 郑州大学 | 一种基于声发射的平行钢丝拉索损伤空间定位算法 |
Family Cites Families (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH03261810A (ja) * | 1990-03-12 | 1991-11-21 | Fujitsu Ltd | 外観検査装置 |
JP3077739B2 (ja) * | 1995-04-11 | 2000-08-14 | 東京電力株式会社 | 部分放電測定方法 |
JPH10153638A (ja) * | 1996-11-26 | 1998-06-09 | Fujikura Ltd | 電力ケーブル気中終端部部分放電測定方法 |
JPH11118939A (ja) * | 1997-10-17 | 1999-04-30 | Toshiba Corp | 地震観測システム |
-
2003
- 2003-03-20 JP JP2003077595A patent/JP4065946B2/ja not_active Expired - Lifetime
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105137485A (zh) * | 2015-08-31 | 2015-12-09 | 中国石油天然气股份有限公司 | 一种地震剖面位图中解释信息的自动拾取转化方法及装置 |
CN105137485B (zh) * | 2015-08-31 | 2017-12-05 | 中国石油天然气股份有限公司 | 一种地震剖面位图中解释信息的自动拾取转化方法及装置 |
Also Published As
Publication number | Publication date |
---|---|
JP2004286527A (ja) | 2004-10-14 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Jin et al. | Surface wave phase-velocity tomography based on multichannel cross-correlation | |
US7636275B2 (en) | Direct time lapse inversion of seismic data | |
CN106052848A (zh) | 双测量面近场声全息噪声源识别系统 | |
Xiantai et al. | Adaptive picking of microseismic event arrival using a power spectrum envelope | |
US10877175B2 (en) | Seismic acquisition geometry full-waveform inversion | |
CN110361792B (zh) | 一种地球物理数据融合及成像方法、介质与设备 | |
NO20121471A1 (no) | Fremgangsmåte og system for presentasjon av seismisk informasjon | |
CN108873063B (zh) | 一种微地震矩张量反演的方法及装置 | |
JP4065946B2 (ja) | 多成分ae波形の初動検出方法 | |
Liu et al. | Investigating DEM error patterns by directional variograms and Fourier analysis | |
CN109425894A (zh) | 一种地震异常道检测方法及装置 | |
CN113050189A (zh) | 一种测井曲线的重构方法、装置、设备及存储介质 | |
Klosterman et al. | Modal survey activity via frequency response functions | |
CN111694053A (zh) | 初至拾取方法及装置 | |
Bonekemper et al. | Automated set-up parameter estimation and result evaluation for SSI-Cov-OMA | |
De Matteis et al. | BISTROP: Bayesian inversion of spectral‐level ratios and P‐wave polarities for focal mechanism determination | |
NO343122B1 (no) | Analyse av geologiske strukturer basert på seismiske attributter | |
CN110109184B (zh) | 一种基于多日变点的被动场源类三维电场勘探方法 | |
CN112649848B (zh) | 利用波动方程求解地震波阻抗的方法和装置 | |
D’Angelo et al. | A new picking algorithm based on the variance piecewise constant models | |
NO20120514A1 (no) | Linje- og kantdeteksjon og forbedring | |
Moore‐Driskell et al. | Integration of arrival‐time datasets for consistent quality control: A case study of amphibious experiments along the Middle America Trench | |
Stepnov et al. | Earthworm-based automatic system for real-time calculation of local earthquake source parameters | |
Yin et al. | Research on interference signal recognition in p wave pickup and magnitude estimation | |
CN110095814B (zh) | 一种地壳应力状态检测方法、系统及相关组件 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20041124 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20070410 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20070530 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20070904 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20071105 |
|
TRDD | Decision of grant or rejection written | ||
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20071211 |
|
R150 | Certificate of patent or registration of utility model |
Ref document number: 4065946 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
S533 | Written request for registration of change of name |
Free format text: JAPANESE INTERMEDIATE CODE: R313533 |
|
R350 | Written notification of registration of transfer |
Free format text: JAPANESE INTERMEDIATE CODE: R350 |
|
EXPY | Cancellation because of completion of term |