JP4867016B2

JP4867016B2 - Earthquake prediction system, earthquake prediction method, earthquake prediction program, and recording medium recording earthquake prediction program

Info

Publication number: JP4867016B2
Application number: JP2007510376A
Authority: JP
Inventors: 正士早川
Original assignee: THE UNIVERSITY OF ELECTRO-COMUNICATINS
Current assignee: THE UNIVERSITY OF ELECTRO-COMUNICATINS
Priority date: 2005-03-29
Filing date: 2006-03-15
Publication date: 2012-02-01
Anticipated expiration: 2026-03-15
Also published as: WO2006103929A1; JPWO2006103929A1

Description

【技術分野】
【０００１】
本発明は、地震の予兆として発生する電離層擾乱が既存の超長波（ＶＬＦ）／長波（ＬＦ）に影響を及ぼすことを利用して創作した地震予知システム、地震予知方法、地震予知プログラム及び地震予知プログラムを記録した記録媒体に関する。
【背景技術】
【０００２】
従来における地震の予知技術は地殻の移動を検知して地盤内に蓄積された歪みエネルギーを推定することが基本となっていた。すなわち、三角測量法によって地盤の移動を計測して歪みエネルギーの分布状態を算出し、地盤が歪み応力に耐えられなくなって正断層もしくは逆断層を生じる可能性を推測していた。近年に至ってレーザー測量技術が適用されたり歪ゲージの性能が向上したりして地震予知技術が向上したが、本質的には地盤が歪み応力に耐えられなくなる状態を予測するものであるから、危険地域を指定することは出来ても、地震発生の期日については１０年単位の予告すら困難であった。すなわち、例えば相模湾トラフに地盤の歪みエネルギーが蓄積されていて断層を生じる危険が有るということは指摘できても、地盤の耐力が限界に達して断層を生じる時期については１０年以上の誤差を含んだ予告しかできない。このため有史以来の大地震発生の周期性までが参考に供されている。最近の例としては、阪神大震災に関して予告した者が無かった。
【特許文献１】
特開平８−３３４５６９号公報
【特許文献２】
特許第２８７５３９８号公報
【非特許文献１】
A.V. Shvets, M.Hayakawa, O.A.Molchanov, and Y.Ando, A study ofionospheric response to regional seismic activity by VLF radio sounding,Physics and Chemistry of the Earth, vol. 29, 627-638, 2004
【発明の開示】
【発明が解決しようとする課題】
【０００３】
地震発生に先立って岩盤の破壊が始まろうとするとき、地表から空気中へラドンなどが放出され、大気プラズマが形成されると言う現象が知られている。このことを具体的に図を用いて説明する。図１には、地球１０と、下部電離層１１と、震源１２と、擾乱域１３とが示されている。擾乱域１３は、下部電離層１１に局所的に発生した異常電離部分であり、震源の上空に発生するものである。そして、擾乱域１３の下部電離層１１が降下している。
【０００４】
図２に示されるように、送信点２０から受信点２１へ送信されたＶＬＦ波は、下部電離層１１と地球１０の間の導波管の内を伝播するが、この擾乱域１３により、下部電離層が降下しているので伝播距離が短くなり位相が進む。
【０００５】
この大気プラズマの本質や、それに伴う電離層擾乱の発生メカニズムについては、目下解明が進められているが、地震発生の数日前の夜間に発光が目視されたり、通信電波が障害を被ったりすることから、前述の大気プラズマの効果や待機の振動などが一原因であると考えられている。前記の電波障害の情況を解析することによって地震の予兆を検知し得るのではなかろうかという着想の下に研究が試みられたが、次に述べるような理由により、電波障害の情況を解析して地震の発生を予知することは出来なかった。すなわち、
ａ．大気プラズマによる電波障害は、通信電波の映像波形の中に混入するノイズとして、もしくは音声波形の中に混入するノイズとして捉えられていた。従って、ノイズの混入を受ける原波形自体が既に複雑な振幅変調もしくは複雑な周波数変調を受けているのでノイズの選別除去が容易でなく、その上、ノイズ選別除去前の原波形とノイズを含んだ受信波形との比較も困難である。
ｂ．通信電波の中に混入するノイズとしては、地震発生を予告する電離層プラズマに因るノイズの他に、イ．地球電離圏／磁気圏で発生する電波に因るノイズ（以下、宇宙由来のノイズと呼ぶ）、ロ．地球大気圏の気象現象である雷などの空電に由来するノイズ（以下、気象由来のノイズと呼ぶ）、および、ハ．もろもろの電気機器から発生する人工的なノイズが有るので、これらのノイズの中から地殻に由来する電離層プラズマに因るノイズを選別することが困難である。
【０００６】
本発明は上述の事情に鑑みて為されたものであって、地震発生の予兆である電離層プラズマに因る電波伝搬の変化を検出して高い確率で地震を予知し得る地震予知システム、地震予知方法、地震予知プログラム及び地震予知プログラムを記録した記録媒体を提供することを目的とする。
課題を解決するための手段
［０００７］
［０００８］
［０００９］
［００１０］
［００１１］
［００１２］
［００１３］
［００１４］
［００１５］
また、上記課題を解決するために、本発明は、第１の地震予知方法と、第２の地震予知方法とを併用した地震予知方法であって、前記第１の地震予知方法は、所定の送信局より送信された所定の振幅と位相を有するＶＬＦ／ＬＦ信号を受信する信号受信段階と、受信したＶＬＦ／ＬＦ信号と基準となるＶＬＦ／ＬＦ信号の位相差を求め、該位相差を記憶する位相差取得段階と、所定の時間帯に記憶した位相差の最小値を求める最小位相差取得段階と、前記最小値をとった時刻である最小位相差時刻を求める時刻取得段階と、前記最小位相差時刻の変動に基づき、地震発生を予知する第１の地震予知段階とを有し、前記第１の地震予知段階は、所定の期間における前記最小位相差時刻の平均値と標準偏差とを求める統計処理段階と、前記平均値と前記最小位相差時刻との差分の絶対値が前記標準偏差の２倍を超えるかどうか判断する異常判断段階とを有し、前記異常判断段階が前記標準偏差の２倍を超えたと判断すると、地震発生を予知し、前記第２の地震予知方法は、前記ＶＬＦ／ＬＦ信号の所定の時刻における振幅を求め、該振幅を記憶する振幅取得段階と、所定の期間の前記振幅の平均値と、
前記振幅に基づき、地震発生を予知する第２の地震予知段階とを有し、前記第２の地震予知段階は、前記平均値と、前記振幅との差分をとり、その差分の２乗を求める差分取得段階と、前記所定の期間における前記差分の２乗の平均値と標準偏差とを求める統計処理段階と、前記差分の２乗が、前記平均値に前記標準偏差の２倍を加えた値を超えるかどうか判断する異常判断段階とを有し、前記異常判断段階が前記平均値に前記標準偏差の２倍を加えた値を超えたと判断すると、地震発生を予知し、前記第１の地震予知方法と、前記第２の地震予知方法とを併用したことを特徴とする。
［００１６］
［００１７］
［００１８］
［００１９］
また、上記課題を解決するために、本発明は、前記所定の時間帯は、日出時刻または日没時刻を含むことを特徴とする。
［００２０］
［００２１］
［００２２］
また、上記課題を解決するために、本発明は、前記所定の時刻は、午後２１時から午前３時までの時刻であることを特徴とする。
［００２３］
また、上記課題を解決するために、本発明は、地震予知方法をコンピュータに実行させるための地震予知プログラム及び該地震予知プログラムを記録したコンピュータが読み取り可能な記録媒体として構成することができる。また、上記課題を解決するために、本発明は、第１の地震予知システムと、第２の地震予知システムとを併用した地震予知システムであって、前記第１の地震予知システムは、所定の送信局より送信された所定の振幅と位相を有するＶＬＦ／ＬＦ信号を受信する信号受信手段と、受信したＶＬＦ／ＬＦ信号と基準となるＶＬＦ／ＬＦ信号の位相差を求め、該位相差を記憶する位相差取得手段と、所定の時間帯に記憶した位相差の最小値を求める最小位相差取得手段と、前記最小値をとった時刻である最小位相差時刻を求める時刻取得手段と、前記最小位相差時刻の変動に基づき、地震発生を予知する第１の地震予知手段とを有し、前記第１の地震予知手段は、所定の期間における前記最小位相差時刻の平均値と標準偏差とを求める統計処理手段と、前記平均値と前記最小位相差時刻との差分の絶対値が前記標準偏差の２倍を超えるかどうか判断する異常判断手段とを有し、前記異常判断手段が前記標準偏差の２倍を超えたと判断すると、地震発生を予知し、前記第２の地震予知システムは、前記ＶＬＦ／ＬＦ信号の所定の時刻における振幅を求め、該振幅を記憶する振幅取得手段と、所定の期間の前記振幅の平均値と、前記振幅に基づき、地震発生を予知する第２の地震予知手段とを有し、前記第２の地震予知手段は、前記平均値と、前記振幅との差分をとり、その差分の２乗を求める差分取得手段と、前記所定の期間における前記差分の２乗の平均値と標準偏差とを求める統計処理手段と、前記差分の２乗が、前記平均値に前記標準偏差の２倍を加えた値を超えるかどうか判断する異常判断手段とを有し、前記異常判断手段が前記平均値に前記標準偏差の２倍を加えた値を超えたと判断すると、地震発生を予知するように構成することができる。また、上記課題を解決するために、本発明は、第１の地震予知方法と、第２の地震予知方法とを併用した地震予知方法であって、前記第１の地震予知方法は、所定の送信局より送信された所定の振幅と位相を有するＶＬＦ／ＬＦ信号を受信する信号受信段階と、受信したＶＬＦ／ＬＦ信号と基準となるＶＬＦ／ＬＦ信号の位相差を求め、該位相差を記憶する位相差取得段階と、日出時刻または日没時刻を含む所定の時間帯に記憶した位相差の最小値を求める最小位相差取得段階と、前記最小値をとった時刻である最小位相差時刻を求める時刻取得段階と、前記最小位相差時刻の変動に基づき、地震発生を予知する第１の地震予知段階とを有し、前記第１の地震予知段階は、所定の期間における前記最小位相差時刻の平均値と標準偏差とを求める統計処理段階と、前記平均値と前記最小位相差時刻との差分の絶対値が前記標準偏差の２倍を超えるかどうか判断する異常判断段階とを有し、前記異常判断段階が前記標準偏差の２倍を超えたと判断すると、地震発生を予知し、前記第２の地震予知方法は、前記ＶＬＦ／ＬＦ信号の午後２１時から午前３時までの所定の時刻における振幅を求め、該振幅を記憶する振幅取得段階と、所定の期間の前記振幅の平均値と、前記振幅に基づき、地震発生を予知する第２の地震予知段階とを有し、前記第２の地震予知段階は、前記平均値と、前記振幅との差分をとり、その差分の２乗を求める差分取得段階と、前記所定の期間における前記差分の２乗の平均値と標準偏差とを求める統計処理段階と、前記差分の２乗が、前記平均値に前記標準偏差の２倍を加えた値を超えるかどうか判断する異常判断段階とを有し、前記異常判断段階が前記平均値に前記標準偏差の２倍を加えた値を超えたと判断すると、地震発生を予知し、前記第１の地震予知方法と、前記第２の地震予知方法とを併用するように構成することができる。
また、上記課題を解決するために、本発明は、第１の地震予知システムと、２の地震予知システムとを併用した地震予知システムであって、前記第１の地震予知システムは、所定の送信局より送信された所定の振幅と位相を有するＶＬＦ／ＬＦ信号を受信する信号受信手段と、受信したＶＬＦ／ＬＦ信号と基準となるＶＬＦ／ＬＦ信号の位相差を求め、該位相差を記憶する位相差取得手段と、日出時刻または日没時刻を含む所定の時間帯に記憶した位相差の最小値を求める最小位相差取得手段と、前記最小値をとった時刻である最小位相差時刻を求める時刻取得手段と、前記最小位相差時刻の変動に基づき、地震発生を予知する第１の地震予知手段とを有し、前記第１の地震予知手段は、所定の期間における前記最小位相差時刻の平均値と標準偏差とを求める統計処理手段と、前記平均値と前記最小位相差時刻との差分の絶対値が前記標準偏差の２倍を超えるかどうか判断する異常判断手段とを有し、前記異常判断手段が前記標準偏差の２倍を超えたと判断すると、地震発生を予知し、前記第２の地震予知システムは、前記ＶＬＦ／ＬＦ信号の午後２１時から午前３時までの所定の時刻における振幅を求め、該振幅を記憶する振幅取得手段と、所定の期間の前記振幅の平均値と、前記振幅に基づき、地震発生を予知する第２の地震予知手段とを有し、前記第２の地震予知手段は、前記平均値と、前記振幅との差分をとり、その差分の２乗を求める差分取得手段と、前記所定の期間における前記差分の２乗の平均値と標準偏差とを求める統計処理手段と、前記差分の２乗が、前記平均値に前記標準偏差の２倍を加えた値を超えるかどうか判断する異常判断手段とを有し、前記異常判断手段が前記平均値に前記標準偏差の２倍を加えた値を超えたと判断すると、地震発生を予知するように構成することができる。
発明の効果
【００２４】
本発明によれば、地震発生の予兆である電離層プラズマに因るノイズを検出して高い確率で地震を予知し得る地震予知システム、地震予知方法、地震予知プログラム及び地震予知プログラムを記録した記録媒体を提供できる。
【図面の簡単な説明】
【００２５】
【図１】大気プラズマと電離層とを示す図である。
【図２】ＶＬＦ波が擾乱域と地球の間を伝播する様子を示す図である。
【図３】ＶＬＦ波が伝播する様子を示す図である。
【図４】オメガ発信基地の位置を示す図である。
【図５】対馬オメガ局と観測点（犬吠）とを示す図である。
【図６】犬吠でのＶＬＦデータを示す図である。
【図７】地震予知システムの構成図である。
【図８】クライアントＰＣとサーバＰＣのハードウェア構成を示す図である。
【図９】位相と振幅を測定する原理を示す図である。
【図１０】地震予知システムのソフトウェア構成を示す図である。
【図１１】ターミネータ・タイムの時系列変化を示す図である。
【図１２】ターミネータ・タイム法の処理を示すフローチャートである。
【図１３】振幅と、振幅の平均の時系列変化を示す図である。
【図１４】夜間ゆらぎ量の時系列データを示す図である。
【図１５】ゆらぎ量の時系列変化とマグニチュードを示す図である。
【図１６】ゆらぎ量を用いた方法の処理を示すフローチャートである。
【図１７】ゆらぎ量と地震発生との相互相関関係を時系列で示した図である。
【図１８】現在のＶＬＦ／ＬＦ観測点を示す図である（その１）。
【図１９】現在のＶＬＦ／ＬＦ観測点を示す図である（その２）。
【図２０】新潟中越地震で用いたパスを示す図である。
【図２１】紀伊半島南東沖地震と新潟中越地震の規模を示す図である。
【図２２】ターミネータ・タイムの時系列変化をプロットした図である。
【符号の説明】
【００２６】
１０地球
１１下部電離層
１２震源
１３擾乱域
２０送信点
２１受信点
３０垂直アンテナ
３１ＧＰＳアンテナ
３２前置増幅部
３３サービスユニット
３４オシレータ部
３５乗算部
３５ａ、３５ｂ乗算器
３６ローパスフィルタ部
３６ａ、３６ｂローパスフィルタ
３７Ａ／Ｄ変換部
３８クライアントＰＣ
３９サーバＰＣ
４０電話回線
４１データ伝送線
４３増幅部
４４９０°移相器
４７測定器
６０ＶＬＦ／ＬＦ送信局
６１電離層プラズマ
６２観測点
６３大地
１２０演算処理装置
１２１インタフェース装置
１２２入力装置
１２３表示装置
１２４ドライブ装置
１２５記録媒体
１２６補助記憶装置
１２７メモリ装置
２０１位相差取得部
２０２最小位相差取得部
２０３振幅取得部
２０４最小振幅取得部
２０５時刻取得部
２０６差分取得部
２０７統計処理部
２０８異常判断部
【発明を実施するための最良の形態】
【００２７】
本実施の形態では、一定周波数、一定振幅の長波電波の発信局と受信局とによって、両局を結ぶ監視線を設定する。ＶＬＦ／ＬＦ信号は、図３に示されるように、大出力のＶＬＦ／ＬＦ送信局６０から送信され、電離層プラズマ６１と大地３１との間の導波管の内を伝播し、観測点６２で受信される。
【００２８】
受信局は比較的小規模の機器で構成でき消費電力も小さいが、これに比して発信局は大規模設備であり消費電力も大きい。
【００２９】
そこで本実施の形態は電波航法援助システムとして既に設置されているオメガシステムの発信基地から送信される電波を利用する。上記のオメガ発信基地は図４に示されるように全世界の各地に配置され、１０．２ｋＨｚ，１１．３ｋＨｚ，１２．０５ｋＨｚ，または１３．６ｋＨｚの一定周波数・一定振幅の長波電波を放射しているから、受信局を設置すれば比較的容易かつ安価に、本発明に係る監視網を構成することができる。
【００３０】
なお、図４には数年前に停波したオメガ発信基地を全部記載してあるが、オメガシステムの受信基地は多数有るので、本実施の形態に関連する受信基地のみを記載した。
【００３１】
次に、世界で最も顕著な地震に伴う電離層プラズマ擾乱によるＶＬＦ伝播異常を紹介する。図５は、対馬オメガ局と観測点（犬吠）とを示す図である。対馬オメガ局と犬吠とを結ぶ大円の近傍に神戸地震（マグニチュード７．２、深さ２０ｋｍ）が１９９５年１月１７日に発生した。犬吠でのＶＬＦデータ（振幅も位相も同様の結果）を図６に示す。
【００３２】
図６に示されるグラフは、横軸が時刻を示し、左の縦軸が日付を示し、右の縦軸が位相差を示している。日出と日没時付近において位相差が最小となる時刻をターミネータ・タイムというが、図６中に示されているｔｍは、日出付近のターミネータ・タイムであり、ｔｅは、日没付近のターミネータ・タイムである。
【００３３】
この図において、ｔｍ、ｔｅにおける位相差が、影を受けた部分のようになって、通常とは異なる時間、異なる位相差を示し、異常が示されていることがわかる。この異常は、地震の数日前から地震日までとなっている。そして、このときのｔｍ、ｔｅと、所定の期間のｔｍ、ｔｅを統計値としたときの平均値との差分の絶対値は、標準偏差（以下、σと表現することがある）の２倍を超えている。
【００３４】
その後、多くの事例に対する発明者らの解析から、ターミネータ・タイムの変動に対して同じ基準（２σ）を用いて、マグニチュードが大きく（６．０以上）、浅い（深さ５０ｋｍ）、しかも送信局と観測点とを結ぶ大円の近くに発生する地震に対しては、約８０％の確率にて電離層擾乱が発生し、ＶＬＦ伝播異常が発生することがわかった。更にＶＬＦ伝播異常の数日前から約１週間程度の時間遅れがあることもわかっている。
【００３５】
そこで本実施の形態では、受信波形から宇宙由来のノイズと、気象由来のノイズとを選別除去する。ここ何年かの研究により宇宙由来のノイズと気象由来のノイズの諸特性（特に継続時間）が良く知られているので、濾波回路にて除去することができる。
【００３６】
従来のＶＬＦ受信器は、単一の周波数をターゲットとして受信するシステムであった。即ち、ある特定の周波数を送信する局電波は、その周波数に同調した受信器にて受信し、別の周波数を送信する局電波は別の受信器が受信する必要があった。そこで本地震予知システムでは、垂直アンテナにて広帯域（１０ｋＨｚ〜５０ｋＨｚ）にて受信し、前置増幅器を通し、サービスユニットへ導入される。
【００３７】
この地震予知システムの構成図を図７に示す。図７には、垂直アンテナ３０と、ＧＰＳアンテナ３１と、前置増幅器３２と、サービスユニット部３３と、クライアントＰＣ３８と、サーバＰＣ３９と、電話回線４０と、データ伝送線４１とが示されている。
【００３８】
このうち、ＧＰＳアンテナ３１は、位相測定のための時刻標準としてＧＰＳ信号を受信するためのものである。クライアントＰＣ３８は、サーバＰＣ３９に送信するデータを作成するとともに、そのデータを処理するＰＣである。サーバＰＣ３９は、各観測点で観測したデータを解析するＰＣである。
【００３９】
また、サービスユニット３３は、オシレータ部３４と、乗算部３５と、ローパスフィルタ部３６と、Ａ／Ｄ変換部３７とを含む。オシレータ３４は、受信したＧＰＳ信号に基づき、オメガ発信基地の周波数と同じ周波数の信号を発生する。乗算部３５は、オシレータ３４が発生する信号と受信した信号を乗算して同期検波を行う。ローパスフィルタ部３６は、同期検波した信号からノイズを除去するために用いられる。Ａ／Ｄ変換部３７は、クライアントＰＣ３８にデータをディジタル信号として送信するために、同期検波した信号をディジタル信号に変換する。
【００４０】
クライアントＰＣ３８とサーバＰＣ３９のハードウェア構成図を、図８を用いて説明する。図８に示されるように、ハードウェア構成は、それぞれバスＢで相互に接続されている入力装置１２２と、表示装置１２３と、ドライブ装置１２４と、記録媒体１２５と、補助記憶装置１２６と、メモリ装置１２７と、演算処理装置１２０と、インタフェース装置１２１とを含む。
【００４１】
入力装置１２２は、キーボード及びマウスなどで構成され、コンピュータシステムに各種操作信号を入力するために用いられる。表示装置１２３は、コンピュータシステムを操作するのに必要な各種ウィンドウやデータ等を表示する。インタフェース装置１２１は、コンピュータシステムをネットワークに接続する為のインタフェースや他の機器を接続するためのインタフェースであり、例えばＮＩＣ（Network Interface Card）やモデム、ＵＳＢ（Universal Serial Bus）等で構成される。
【００４２】
そして、上記コンピュータシステムを動作させるためのプログラムは、ＣＤ−ＲＯＭ等の記録媒体１２５によって提供されるか、ネットワークを通じてダウンロードされる。また、記録媒体１２５は、ドライブ装置１２４にセットされ、データやプログラムが記録媒体１２５からドライブ装置１２４を介して補助記憶装置１２６にインストールされる。
【００４３】
補助記憶装置１２６は、データやプログラムを格納されると共に、必要なファイル等を格納する。メモリ装置１２７は、コンピュータシステムの起動時に補助記憶装置１２６からプログラムを読み出して格納する。演算処理装置１２０は、メモリ装置１２７に読み出され格納されたプログラムに従って処理を実行する。
【００４４】
以上説明した地震予知システムで、位相と振幅を測定する原理を示したのが図９である。図７で説明した符号と同一の符号については説明を省略する。
【００４５】
図９には、アンテナ３０と、前置増幅器３２と、乗算器３５ａ、３５ｂと、オシレータ３４と、ローパスフィルタ３６ａ、３６ｂと、９０°移相器４４と、測定器４７とが示されている。
【００４６】
乗算器３５ａ、３５ｂならびに９０°移相器４４は、図７の乗算部３５に対応する。ローパスフィルタ３６ａ、３６ｂは、図７のローパスフィルタ部３６に対応する。測定器４７は、クライアントＰＣ３８に対応する。
【００４７】
この図９は、基準となるオシレータ３４からの信号（同相）とアンテナで受信したＶＬＦ／ＬＦ信号とを乗算し、更に、オシレータ３４からの信号の位相を９０度ずらしたもの（直交）とＶＬＦ／ＬＦ信号とを乗算して、受信した信号の同相成分と直交成分を検出する。測定器４７は、これらの同相成分と直交成分に基づいて、図に示すように、受信した信号の位相θと振幅Ａを得る。
【００４８】
次に本実施の形態における地震予知システムのソフトウェア構成について、図１０を用いて説明する。図１０には、位相差取得部２０１と、最小位相差取得部２０２と、振幅取得部２０３と、最小振幅取得部２０４と、時刻取得部２０５と、差分取得部２０６と、統計処理部２０７と、異常判断部２０８とが示されている。
【００４９】
位相差取得部２０１は、上記した位相θに基づいて、ＶＬＦ／ＬＦ信号の位相差を求め、その位相差を記憶する。最小位相差取得部２０２は、所定の時間帯に記憶した位相差の最小値を求める。振幅取得部２０３は、ＶＬＦ／ＬＦ信号の振幅Ａを求め、その振幅を記憶する。最小振幅取得部２０４は、所定の時間帯に記憶した振幅の最小値を求める。時刻取得部２０５は、位相差または振幅の最小値をとった時刻を求める。差分取得部２０６は、所定の期間の振幅の平均値と、ある振幅との差分をとり、その差分の２乗を求める。統計処理部２０７は、所定の期間の最小位相差時刻、所定の期間の最小振幅時刻、所定の期間の差分の２乗、の３つの統計値における平均値と標準偏差とを求める。異常判断部２０８は、所定の期間の最小位相差時刻の平均値と、取得した最小位相差時刻との差分の絶対値が標準偏差の２倍を超えるかどうか判断する。また、異常判断部２０８は、所定の期間の最小振幅時刻の平均値と、取得した最小振幅時刻との差分の絶対値が標準偏差の２倍を超えるかどうか判断する。さらに、異常判断部２０８は、取得した差分の２乗が、所定の期間の差分の２乗の平均値に標準偏差の２倍を加えた値を超えるかどうか判断する。
【００５０】
以上説明した地震予知システムにより、広帯域信号の中から所望のオメガ発信基地の送信信号（例えばオーストラリアＮＷＣ局（１９．８ｋＨｚ）、中国北京ＣＨＩ局（２２．２ｋＨｚ）、ハワイＮＰＭ局（２１．４ｋＨｚ）、）の振幅と位相を連続観測する。同様に、日本のＪＪＹ局（福島、４０ｋＨｚ）の信号の振幅と位相を測定するようにしてもよい。なお。位相測定のための時刻標準としてＧＰＳ信号を受信している。国内多点での観測点からデータは電話回線を通してマスターステーション（電通大）へ毎日伝送され、解析されるシステムになっている。
【００５１】
各観測点で観測されるＶＬＦ／ＬＦ信号（振幅と位相）は、上述したように毎日電通大のマスターステーションへ伝送される。その一日の変化パターン（図６参照）が振幅と位相に対して得られる。本実施の形態ではＶＬＦ／ＬＦ伝播異常を検出する信号処理法として二種類を採用する。その一つは、ターミネータ・タイムを用いるターミネータ・タイム法と、もう一つは夜間ゆらぎスペクトルを用いる方法である。
【００５２】
まず、ターミネータ・タイム法を用いる場合について説明する。図６のＶＬＦ信号（例えば振幅）の日変化パターンにおいてターミネータ・タイムを読み取る。ソフトウェアが自動的に読み取り、一日に一ヶ（即ち、ｔｍとｔｅ）を得る。例えば、ｔｅの毎日のデータを用い、当該日の過去１０日間のデータを用いて、ｔｅの平均値と分散（標準偏差）を計算する。
【００５３】
そして、図１１に示すように、平均値とｔｅとの差の絶対値が標準偏差σの２倍を超える日をＶＬＦ／ＬＦ伝播異常と定義する。この方法によるアルゴリズムを開発し、計算機に自動的に図１１を作成する手法を用いる。さらに異常が発生した時には生の元データ（日変化パターン）に立ち帰り、ターミネータ・タイムの異常を確認する。ターミネータ・タイムはｔｍ、ｔｅともに同様の処理を行う。また、上述した説明では、位相差を用いるものであったが、振幅を用いても良い。この場合も同様に、振幅が最小となる時刻がターミネータ・タイムである。
【００５４】
なお、図１１に示されるグラフは、横軸が日付を示し、縦軸がその日付におけるｔｅと平均値との差を示している。
【００５５】
このターミネータ・タイム法を用いる場合の処理を、図１２のフローチャートを用いて説明する。この処理は地震予知システム全体が実行する処理である。ステップＳ１０１で、ＶＬＦ／ＬＦ信号を受信する。ステップＳ１０２で位相差（又は振幅）を求め、記憶する。ステップＳ１０３で、所定の時間が経過したかどうか判断する。ステップＳ１０４で所定の時間内に記憶した位相差（又は振幅）の最小値を求める。ステップＳ１０５で、最小値をとった時刻を求め、それをターミネータ・タイムＴとする。ステップＳ１０６で、所定の期間におけるターミネータ・タイムの平均値と標準偏差σを求める。ステップＳ１０７で、平均値とＴの差の絶対値が２σを超えているか判断する。超えている場合、ステップＳ１０８で異常ありと判断し、地震発生を予知する。超えていない場合は、異常なしと判断する。
【００５６】
地震発生の予知は、例えばサーバＰＣの画面に表示するなどして、観測者に通知するようにする。
【００５７】
次に、夜間ゆらぎスペクトルを用いる方法について説明する。この信号処理法は、図６のような日変化パターンにおいて夜間の時刻（地方時２１時から３時まで）の変動に注目するものである。観測された日変化パターンに対して次式の操作を行う。
【００５８】
【数１】

但し、ｋは当該日を示し、Ｎは平均値を取る日数で、通常はＮ＝１０日とする。また、Ａ（ｔ）はある日の振幅の日変化パターンを示し、＜Ａ（ｔ）＞は前Ｎ日間の平均日変化パターンを示し、ｄＡ（ｔ）はその差分を示している。
【００５９】
図１３に示されるグラフは、＜Ａ（ｔ）＞とｄＡ（ｔ）の理解を助けるための説明用のグラフであって、＜Ａ（ｔ）＞とｄＡ（ｔ）を１時間毎にプロットしたグラフであり、縦軸が振幅を示し、横軸が時刻を示している。＜Ａ（ｔ）＞は、平均（例えば、１０時間の移動平均）である。したがって、グラフでプロットされている＜Ａ（ｔ）＞の各点は、現時点における値を含む、その１０時間前からの平均を示している。またｄＡ（ｔ）は、各時点における移動平均との差を示している。
【００６０】
これらのプロセスも計算機でのアルゴリズムが開発されており、毎日ｄＡ（ｔ）を計算し、２１時から３時までの時間のｄＡ（ｔ）の２乗を毎日計算し、一日のデータとする。すると、ターミネータ・タイムと同様の時系列データが得られる。平均と分散（標準偏差）の作成はターミネータ・タイムのときと同じ方法を用いる。ゆらぎ量が平均値＋２σを超えた時を異常とする。地震の前兆と考えられる。ゆらぎ量と地震の相関関係には発明者らの最新の論文である非特許文献１を参照されたい。
【００６１】
図１４は、夜間ゆらぎ量の時系列データを示すグラフである。このグラフは、縦軸がゆらぎ量を示し、横軸が日付を示している。このグラフにおいて、平均値＋２σを超えたものを異常と判断する。
【００６２】
次に図１５のグラフについて説明する。このグラフは、ゆらぎ量の変化とマグニチュードを示すもので、横軸が日付を示し、左の縦軸がゆらぎ量を示し、右の縦軸がマグニチュードを示す。太線は平均値＋２σを示し、細線はゆらぎ量を示し、棒線がマグニチュードを示している。このマグニチュードは、送信局と観測点とを結ぶパスの近傍で発生した地震のマグニチュードである（非特許文献１参照）。この図に示されるように、地震の前にはゆらぎの異常が出ていることがわかる。
【００６３】
このゆらぎ量を用いた処理を、図１６のフローチャートを用いて説明する。ステップＳ２０１で、時刻ｔにおける振幅Ａを求め、記憶される。ステップＳ２０２で、所定の期間におけるｔでの振幅の平均値と振幅Ａの差をとり、その差分をｄＡとし、その差分ｄＡの２乗を求める。ステップＳ２０３で、所定の期間におけるｔでの差分ｄＡ２乗の平均値と標準偏差σを求める。ステップＳ２０４で、差分ｄＡ２乗が、平均値＋２σを超えているかが判断される。超えている場合、ステップＳ２０５で異常ありと判断され、地震発生が予知される。超えていない場合は、異常なしと判断される。
【００６４】
地震発生の予知は、例えばサーバＰＣの画面に表示するなどして、観測者に通知するようにする。
【００６５】
次に、地震予知の確度の上昇について説明する。上述したターミネータ・タイムとゆらぎ量との併用により、異常の確度は高めることができる。即ちターミネータ・タイムの異常があり、更にゆらぎ量にも異常が発生している時には、地震前兆としての確率は著しく上昇する。従って、両者を用いることで地震予知の確度を上昇することができる。
【００６６】
では次に、いつ地震が発生するかという地震の発生時期について説明する。ＶＬＦ／ＬＦ伝播異常と地震との時間遅れに関しては、図１５からもわかる。しかし、より定量的にするため、相互相関関係を用いる。図１７はゆらぎ量と地震発生との相互相関関係を時系列にて示したものであり、縦軸が相関係数を示し、横軸が地震発生日を０日とした日付を示している。地震の数日前に相互相関のピークが出現しており、伝播異常の数日後に地震が発生することになる。ターミネータ・タイムの異常ともほぼ同様の時間遅れがある。この時間遅れは、数日から１週間程度である（非特許文献１参照）。
【００６７】
次に、地震が発生する場所の決定について説明する。先ず、以上よりＶＬＦ／ＬＦ伝播異常が検出された時には、いつ発生するかは予知できるが、どこで発生するかに関しては、これまでの解析では、送信局と観測点を結ぶ大円近傍（もう少し定量的に言うと、例えば第５フレネル帯）で発生することが予測される。そこで本実施の形態で用いる手法は多観測点での多送信局受信に基づく手法である。即ち、送信局と観測点とを結ぶ大円の廻りの第５フレネル帯の重なり合いを利用する方法である。
【００６８】
現在のＶＬＦ／ＬＦ観測点を図１８、１９を用いて説明する。図１８は、国内７観測点となるＶＬＦ送信局の受信施設を示している。図中の線は観測点とＶＬＦ送信局を結ぶ大円を示している。
【００６９】
また、図１９は、福島県東部の４０ｋＨｚ標準電波（ＪＪＹ）の送信施設とその観測点を示している。
【００７０】
このように、現在ＶＬＦ／ＬＦ観測点は、北から母子里（北海道）、調布、千葉館山、清水、春日井（名古屋）、舞鶴、高知の７地点で、各観測点ではＣＨＩ局（中国）、ＮＷＣ局（オーストラリア）、ＮＰＭ局（ハワイ）、ＪＪＹ局（福島）の送信局電波を受信している。
【００７１】
図１８及び図１９において、点線のパス上の地点に地震が発生した場合、点線のパス上の測定点で、受信の異常を検出することになる。
【００７２】
したがって、いろいろなパスでの異常の出現状況から、異常を引き起こす電離層擾乱の場所、即ち将来の地震発生位置を決定することができる。このような発明者らのＶＬＦ／ＬＦネットワークでは国内の地震に対しては、約１００ｋｍの精度がある。
【００７３】
また、このＶＬＦ／ＬＦネットワークによる将来の地震位置予報法は、本地震予知システムにて初めて提案されたものである。
【００７４】
以上説明した地震予知システムと送信局ならびに観測点による新潟中越地震に対する結果を示す。まず、このとき用いたパスを図２０に示す。図２０は、観測点と送信局とのパスを示すものである。この図において、実線は第５フレネル帯を示し、このフレネル帯に含まれる位置を震源とする地震には反応すると考えられる。
【００７５】
図２０に示されるように、母子里（ＨＯＫ）はＣＨＩ局（中国）を受信し、調布（ＣＨＯ）、千葉館山（ＣＢＡ）でもＣＨＩ局を受信し、高知（ＫＯＣ）もＣＨＩ局を受信している。更に、ＪＪＹ局を高知にて受信するパスを用いる。
【００７６】
２００４年後半には図２１に示される表の様に、二つの大きな地震が発生した。即ち、紀伊半島南東沖地震と新潟中越地震である。ともにマグニチュード６．０を超える地震であるが、前者は海の中という点が、後者と大きく異なる点である。従来の発明者らの研究で、陸上での浅い地震に対しては、電離層は極めて敏感であることがわかってきている。
【００７７】
図２２は、上述したすべてのパスに対するターミネータ・タイム（ｔｍ）の変化をプロットしたものである。このグラフは、横軸が日付を示し、左の縦軸がｔｅを示し、右の縦軸は、マグニチュードを示す。なお、左の縦軸は、プロットされた点に対してはｔｅを示しているが、棒グラフの場合は平均とｔｅの差分を示している。
【００７８】
また、左の縦軸に記載されているＸＸＸ−ＹＹＹは、ＸＸＸとＹＹＹ間のパスにおけるグラフであることを示している。例えばＣＨＩ−ＨＯＫは、ＣＨＩ局と母子里とのパスである。従って、図２２には、上のグラフからＣＨＩ局と母子里、ＣＨＩ局と調布、ＣＨＩ局と千葉館山、ＣＨＩ局と高知、ＪＪＹ局と高知のパスに対応するグラフが描かれている。
【００７９】
また、図２２の一番下のグラフは、太陽電波フラックスで太陽活動を示し、また地磁気活動指標を示しており、地球周辺での地磁気の荒れ具合を表すものである。このグラフにおいて、実線は太陽電波フラックスを示し、棒線表示は地磁気活動指標（Ａｐ）を示している。そして、右の縦軸は太陽電波フラックスに対応し、左の縦軸は棒線表示は地磁気活動指標に対応し、単位はそれぞれナノテラス、フラックスである。
【００８０】
新潟中越地震の前には、ＣＨＯ、ＣＢＡ観測点でのＣＨＩ局に伝播異常が明瞭に出ている。新潟中越地震両パスのすぐ近傍であるためである。北海道母子里−ＣＨＩパスには何の変化もない。伝播異常は地磁気活動が低い時に発生しており、地震との関係を示唆している。
【００８１】
このように、本実施の形態による地震予知システムは、パスの近傍で発生した地震による伝播異常を正確に把握しており、このようにすることで、地震発生を予知することができる。
【００８２】
本件国際出願は、２００５年３月２９日に出願した日本国特許出願２００５−９５１３７号に基づく優先権を主張するものであり、特許出願２００５−９５１３７号の全内容を本国際出願に援用する。【Technical field】
[0001]
The present invention relates to an earthquake prediction system, an earthquake prediction method, an earthquake prediction program, and an earthquake prediction that are created by utilizing the influence of ionospheric disturbances that occur as a sign of an earthquake on existing very long waves (VLF) / long waves (LF). The present invention relates to a recording medium on which a program is recorded.
[Background]
[0002]
Conventional earthquake prediction technology is based on detecting the movement of the crust and estimating the strain energy accumulated in the ground. In other words, the movement of the ground was measured by the triangulation method to calculate the strain energy distribution, and the possibility that the ground could not withstand the strain stress and would cause a normal fault or a reverse fault was estimated. In recent years, laser prediction technology has been applied and the performance of strain gauges has improved, so earthquake prediction technology has improved, but it is essentially a prediction of the situation where the ground cannot withstand strain stress. Although it was possible to specify the area, it was difficult to give a 10-year notice about the date of the earthquake. That is, for example, although it can be pointed out that there is a risk of generating a fault due to the accumulated strain energy of the ground in the Sagami Bay trough, there is an error of more than 10 years regarding the time when the ground strength reaches the limit and the fault is generated. You can only include an advance notice. For this reason, the periodicity of the occurrence of major earthquakes since its history has been used as a reference. As a recent example, no one has foretold the Great Hanshin Earthquake.
[Patent Document 1]
Japanese Patent Application Laid-Open No. 8-33469
[Patent Document 2]
Japanese Patent No. 2875398
[Non-Patent Document 1]
AV Shvets, M. Hayakawa, OAMolchanov, and Y. Ando, A study of ionospheric response to regional seismic activity by VLF radio sounding, Physics and Chemistry of the Earth, vol. 29, 627-638, 2004
DISCLOSURE OF THE INVENTION
[Problems to be solved by the invention]
[0003]
It is known that radon is released from the earth's surface into the air and atmospheric plasma is formed when the destruction of the rock mass begins before the occurrence of the earthquake. This will be specifically described with reference to the drawings. In FIG. 1, the earth 10, the lower ionosphere 11, the epicenter 12, and the disturbance region 13 are shown. The disturbance region 13 is an abnormal ionization portion locally generated in the lower ionosphere 11 and is generated above the epicenter. And the lower ionosphere 11 of the disturbance area 13 is falling.
[0004]
As shown in FIG. 2, the VLF wave transmitted from the transmission point 20 to the reception point 21 propagates in the waveguide between the lower ionosphere 11 and the earth 10, and this disturbance region 13 causes the lower ionosphere. , The propagation distance is shortened and the phase advances.
[0005]
The essence of this atmospheric plasma and the mechanism of ionospheric disturbances associated with it are currently being elucidated. However, light emission is visible at night several days before the occurrence of the earthquake, and communication radio waves are damaged. The above-mentioned atmospheric plasma effect and standby vibration are considered to be one cause. Research was attempted based on the idea that an earthquake sign could be detected by analyzing the situation of the radio interference described above, but for the reasons described below, the situation of the radio disturbance was analyzed. It was impossible to predict the occurrence of an earthquake. That is,
a. Radio interference due to atmospheric plasma has been perceived as noise mixed in video waveforms of communication radio waves or noise mixed in audio waveforms. Therefore, since the original waveform itself subject to noise has already undergone complex amplitude modulation or complex frequency modulation, it is not easy to select and remove noise. In addition, the original waveform and noise before noise selection and removal are included. Comparison with the received waveform is also difficult.
b. In addition to the noise caused by ionospheric plasma that predicts the occurrence of an earthquake, the noise mixed in communication radio waves is as follows. Noise caused by radio waves generated in the ionosphere / magnetosphere (hereinafter referred to as space-derived noise), b. Noise derived from aerodynamics such as lightning, which is a meteorological phenomenon in the Earth's atmosphere (hereinafter referred to as meteorological noise), and c. Since there are artificial noises generated from various electric devices, it is difficult to select noises caused by ionospheric plasma derived from the crust from these noises.
[0006]
The present invention has been made in view of the above circumstances, and an earthquake prediction system and an earthquake prediction capable of detecting a change in radio wave propagation caused by ionospheric plasma, which is a sign of earthquake occurrence, and predicting an earthquake with high probability It is an object to provide a method, an earthquake prediction program, and a recording medium on which the earthquake prediction program is recorded.
Means for solving the problem
[0007]
[0008]
[0009]
[0010]
[0011]
[0012]
[0013]
[0014]
[0015]
In order to solve the above problems, the present invention provides an earthquake prediction method using both the first earthquake prediction method and the second earthquake prediction method, wherein the first earthquake prediction method is a predetermined earthquake prediction method. A signal receiving stage for receiving a VLF / LF signal having a predetermined amplitude and phase transmitted from a transmitting station, and obtaining a phase difference between the received VLF / LF signal and a reference VLF / LF signal, and storing the phase difference A phase difference obtaining step, a minimum phase difference obtaining step for obtaining a minimum value of a phase difference stored in a predetermined time zone, a time obtaining step for obtaining a minimum phase difference time which is a time when the minimum value is taken, and the minimum A first earthquake prediction stage for predicting the occurrence of an earthquake based on fluctuations in the phase difference time, and the first earthquake prediction stage includes an average value and a standard deviation of the minimum phase difference time in a predetermined period. Statistical processing stage to be obtained and the average value An abnormality determination step of determining whether an absolute value of a difference from the minimum phase difference time exceeds twice the standard deviation, and determining that the abnormality determination step exceeds twice the standard deviation, The second earthquake prediction method calculates an amplitude of the VLF / LF signal at a predetermined time, stores the amplitude, an average value of the amplitude in a predetermined period,
A second earthquake prediction stage for predicting the occurrence of an earthquake based on the amplitude, wherein the second earthquake prediction stage takes a difference between the average value and the amplitude and obtains a square of the difference. A difference obtaining step, a statistical processing step for obtaining a mean value and a standard deviation of the square of the difference in the predetermined period, and a value obtained by adding the double of the standard deviation to the mean value as the square of the difference An abnormality determination step for determining whether or not the value exceeds, and if the abnormality determination step is determined to exceed a value obtained by adding twice the standard deviation to the average value, the occurrence of an earthquake is predicted, and the first earthquake The prediction method and the second earthquake prediction method are used in combination.
[0016]
[0017]
[0018]
[0019]
In order to solve the above problem, the present invention is characterized in that the predetermined time zone includes a sunrise time or a sunset time.
[0020]
[0021]
[0022]
In order to solve the above problem, the present invention is characterized in that the predetermined time is a time from 21:00 pm to 3:00 am.
[0023]
In order to solve the above-described problems, the present invention can be configured as an earthquake prediction program for causing a computer to execute an earthquake prediction method and a computer-readable recording medium recording the earthquake prediction program. In order to solve the above-described problem, the present invention provides an earthquake prediction system using both the first earthquake prediction system and the second earthquake prediction system, wherein the first earthquake prediction system includes a predetermined earthquake prediction system. A signal receiving means for receiving a VLF / LF signal having a predetermined amplitude and phase transmitted from a transmitting station, a phase difference between the received VLF / LF signal and a reference VLF / LF signal is obtained, and the phase difference is stored. A phase difference acquisition unit that performs the minimum phase difference acquisition unit that obtains the minimum value of the phase difference stored in a predetermined time zone, a time acquisition unit that obtains a minimum phase difference time that is a time when the minimum value is taken, and the minimum First earthquake prediction means for predicting the occurrence of an earthquake based on fluctuations in the phase difference time, and the first earthquake prediction means calculates an average value and a standard deviation of the minimum phase difference time in a predetermined period. Statistical processing hand And an abnormality determining means for determining whether or not an absolute value of a difference between the average value and the minimum phase difference time exceeds twice the standard deviation, and the abnormality determining means reduces twice the standard deviation. If it is determined that the frequency has exceeded, the occurrence of an earthquake is predicted, and the second earthquake prediction system obtains the amplitude of the VLF / LF signal at a predetermined time, stores the amplitude, and the amplitude for a predetermined period. And a second earthquake prediction means for predicting the occurrence of an earthquake based on the amplitude, and the second earthquake prediction means takes a difference between the average value and the amplitude, and the difference A difference obtaining means for obtaining the square of the difference, a statistical processing means for obtaining an average value and a standard deviation of the square of the difference in the predetermined period, and the square of the difference is equal to 2 of the standard deviation in the average value. Judge whether the value is doubled And a abnormality determination means determines that said abnormality determination means exceeds a value obtained by adding twice the standard deviation to the average value, it may be configured to predict earthquakes. In order to solve the above problems, the present invention provides an earthquake prediction method using both the first earthquake prediction method and the second earthquake prediction method, wherein the first earthquake prediction method is a predetermined earthquake prediction method. A signal receiving stage for receiving a VLF / LF signal having a predetermined amplitude and phase transmitted from a transmitting station, and obtaining a phase difference between the received VLF / LF signal and a reference VLF / LF signal, and storing the phase difference A phase difference obtaining step, a minimum phase difference obtaining step for obtaining a minimum value of phase differences stored in a predetermined time zone including a sunrise time or a sunset time, and a minimum phase difference time which is a time when the minimum value is taken And a first earthquake prediction stage for predicting the occurrence of an earthquake based on the fluctuation of the minimum phase difference time, wherein the first earthquake prediction stage includes the minimum phase difference in a predetermined period. Find time average and standard deviation A statistical processing step, and an abnormality determination step of determining whether an absolute value of a difference between the average value and the minimum phase difference time exceeds twice the standard deviation, and the abnormality determination step When it is determined that the frequency exceeds twice, the second earthquake prediction method predicts the amplitude of the VLF / LF signal at a predetermined time from 21:00 pm to 3:00 am and stores the amplitude. An amplitude acquisition step, an average value of the amplitude for a predetermined period, and a second earthquake prediction step for predicting the occurrence of an earthquake based on the amplitude, wherein the second earthquake prediction step includes the average value A difference obtaining step for obtaining a difference between the amplitude and obtaining the square of the difference, a statistical processing step for obtaining an average value and a standard deviation of the square of the difference in the predetermined period, and 2 of the difference The power of the standard deviation An abnormality determination stage that determines whether or not a value that is doubled is exceeded, and determining that the abnormality determination stage exceeds a value that is obtained by adding twice the standard deviation to the average value, predicts the occurrence of an earthquake, The first earthquake prediction method and the second earthquake prediction method can be used in combination.
In order to solve the above problem, the present invention provides an earthquake prediction system using both the first earthquake prediction system and the second earthquake prediction system, wherein the first earthquake prediction system has a predetermined transmission. A signal receiving means for receiving a VLF / LF signal having a predetermined amplitude and phase transmitted from the station, and obtaining a phase difference between the received VLF / LF signal and a reference VLF / LF signal, and storing the phase difference A phase difference acquisition means, a minimum phase difference acquisition means for obtaining a minimum value of a phase difference stored in a predetermined time zone including a sunrise time or a sunset time, and a minimum phase difference time which is a time when the minimum value is taken. A time acquisition means to be obtained; and a first earthquake prediction means for predicting the occurrence of an earthquake based on a change in the minimum phase difference time, wherein the first earthquake prediction means has the minimum phase difference time in a predetermined period. Mean and mark Statistical processing means for obtaining a deviation, and abnormality determination means for determining whether an absolute value of a difference between the average value and the minimum phase difference time exceeds twice the standard deviation, the abnormality determination means When it is determined that the standard deviation exceeds twice the standard deviation, an earthquake is predicted, and the second earthquake prediction system obtains an amplitude of the VLF / LF signal at a predetermined time from 21 pm to 3 am, Amplitude acquisition means for storing the amplitude, average value of the amplitude for a predetermined period, and second earthquake prediction means for predicting the occurrence of an earthquake based on the amplitude, wherein the second earthquake prediction means Difference obtaining means for taking a difference between the average value and the amplitude and obtaining a square of the difference; a statistical processing means for obtaining an average value and a standard deviation of the square of the difference in the predetermined period; The square of the difference is the average An abnormality determining means for determining whether or not a value obtained by adding twice the standard deviation to the average value, and determining that the abnormality determining means exceeds a value obtained by adding twice the standard deviation to the average value, It can be configured to predict the occurrence of an earthquake.
The invention's effect
[0024]
According to the present invention, an earthquake prediction system, an earthquake prediction method, an earthquake prediction program, and a recording medium in which an earthquake prediction program is recorded, which can detect an earthquake with high probability by detecting noise caused by ionospheric plasma that is a sign of an earthquake occurrence. Can provide.
[Brief description of the drawings]
[0025]
FIG. 1 shows atmospheric plasma and ionosphere.
FIG. 2 is a diagram showing how a VLF wave propagates between a disturbance area and the earth.
FIG. 3 is a diagram illustrating a state in which a VLF wave propagates.
FIG. 4 is a diagram showing the position of an omega transmission base.
FIG. 5 is a diagram showing a Tsushima omega station and an observation point (Inubo).
FIG. 6 is a diagram showing VLF data in an inugar.
FIG. 7 is a configuration diagram of an earthquake prediction system.
FIG. 8 is a diagram illustrating a hardware configuration of a client PC and a server PC.
FIG. 9 is a diagram showing the principle of measuring phase and amplitude.
FIG. 10 is a diagram showing a software configuration of the earthquake prediction system.
FIG. 11 is a diagram showing a time-series change in terminator time.
FIG. 12 is a flowchart showing processing of a terminator time method.
FIG. 13 is a diagram showing amplitude and time-series changes in average amplitude.
FIG. 14 is a diagram showing time-series data of nighttime fluctuation amount.
FIG. 15 is a diagram showing a time-series change and a magnitude of a fluctuation amount.
FIG. 16 is a flowchart showing processing of a method using a fluctuation amount.
FIG. 17 is a diagram showing a cross-correlation between fluctuation amount and occurrence of an earthquake in time series.
FIG. 18 is a diagram showing current VLF / LF observation points (No. 1);
FIG. 19 is a diagram showing current VLF / LF observation points (part 2);
FIG. 20 is a diagram showing a path used in the Niigata Chuetsu earthquake.
FIG. 21 is a diagram showing the magnitude of the Kii Peninsula Southeast Offshore Earthquake and the Niigata Chuetsu Earthquake.
FIG. 22 is a diagram plotting time-series changes in terminator time.
[Explanation of symbols]
[0026]
10 Earth
11 Lower ionosphere
12 epicenter
13 Disturbance area
20 transmission points
21 Receiving points
30 Vertical antenna
31 GPS antenna
32 Preamplifier
33 Service unit
34 Oscillator section
35 Multiplier
35a, 35b multiplier
36 Low-pass filter section
36a, 36b Low pass filter
37 A / D converter
38 Client PC
39 Server PC
40 Telephone line
41 Data transmission line
43 Amplifier
44 90 ° phase shifter
47 Measuring instrument
60 VLF / LF transmitter station
61 Ionospheric plasma
62 observation points
63 Earth
120 arithmetic processing unit
121 Interface device
122 Input device
123 Display device
124 drive unit
125 recording media
126 Auxiliary storage device
127 memory device
201 Phase difference acquisition unit
202 Minimum phase difference acquisition unit
203 Amplitude acquisition unit
204 Minimum amplitude acquisition unit
205 Time acquisition unit
206 Difference acquisition unit
207 Statistical processing section
208 Abnormality judgment part
BEST MODE FOR CARRYING OUT THE INVENTION
[0027]
In the present embodiment, a monitoring line connecting the two stations is set by a transmitting station and a receiving station of a long wave having a constant frequency and a constant amplitude. As shown in FIG. 3, the VLF / LF signal is transmitted from the high-power VLF / LF transmission station 60 and propagates in the waveguide between the ionospheric plasma 61 and the ground 31, and at the observation point 62. Received.
[0028]
The receiving station can be composed of a relatively small device and consumes little power. On the other hand, the transmitting station is a large-scale facility and consumes much power.
[0029]
Therefore, this embodiment uses radio waves transmitted from a transmission base of an omega system already installed as a radio navigation assistance system. As shown in FIG. 4, the above-mentioned omega transmission bases are arranged in various parts of the world, and emit long wave radio waves having a constant frequency and a constant amplitude of 10.2 kHz, 11.3 kHz, 12.05 kHz, or 13.6 kHz. Therefore, if a receiving station is installed, the monitoring network according to the present invention can be configured relatively easily and inexpensively.
[0030]
FIG. 4 shows all omega transmission bases stopped several years ago. However, since there are many omega system reception bases, only the reception bases related to the present embodiment are shown.
[0031]
Next, we will introduce VLF propagation anomalies due to ionospheric plasma disturbances associated with the world's most prominent earthquakes. FIG. 5 is a diagram showing a Tsushima omega station and an observation point (Inu-en). The Kobe earthquake (magnitude 7.2, depth 20 km) occurred on January 17, 1995 near the great circle connecting Tsushima Omega station and Inubo. FIG. 6 shows the VLF data (similar results for both amplitude and phase) in the dog cage.
[0032]
In the graph shown in FIG. 6, the horizontal axis indicates time, the left vertical axis indicates date, and the right vertical axis indicates phase difference. The time at which the phase difference between the sunrise and the sunset is the minimum is called the terminator time, but tm shown in FIG. 6 is the terminator time near the sunrise, and te is near the sunset. Terminator time.
[0033]
In this figure, it can be seen that the phase difference at tm and te becomes a shaded part, shows a different phase difference for a different time than usual, and shows an abnormality. This anomaly is from a few days before the earthquake until the earthquake date. The absolute value of the difference between tm and te at this time and the average value when tm and te for a predetermined period are statistical values is twice the standard deviation (hereinafter sometimes expressed as σ). Is over.
[0034]
After that, the inventors' analysis of many cases shows that the same criterion (2σ) is used for the terminator time variation, the magnitude is large (6.0 or more), shallow (depth 50 km), and the transmitting station It was found that ionospheric disturbances occur with an probability of about 80% and an anomaly in VLF propagation occurs for earthquakes that occur near the great circle connecting the observation point and the observation point. Furthermore, it is known that there is a time delay of about one week from several days before the VLF propagation abnormality.
[0035]
Therefore, in the present embodiment, the noise derived from the universe and the noise derived from the weather are selected and removed from the received waveform. The characteristics (especially duration) of space-derived noise and weather-derived noise are well known from research over the past several years, and can be removed by a filtering circuit.
[0036]
Conventional VLF receivers are systems that receive a single frequency as a target. That is, a local radio wave that transmits a specific frequency is received by a receiver tuned to that frequency, and a local radio wave that transmits another frequency needs to be received by another receiver. Therefore, in this earthquake prediction system, the signal is received in a wide band (10 kHz to 50 kHz) with a vertical antenna, and is introduced into a service unit through a preamplifier.
[0037]
A block diagram of this earthquake prediction system is shown in FIG. 7 shows a vertical antenna 30, a GPS antenna 31, a preamplifier 32, a service unit 33, a client PC 38, a server PC 39, a telephone line 40, and a data transmission line 41. .
[0038]
Among these, the GPS antenna 31 is for receiving a GPS signal as a time standard for phase measurement. The client PC 38 creates data to be transmitted to the server PC 39 and sends the data to the server PC 39. processing PC. The server PC 39 is a PC that analyzes data observed at each observation point.
[0039]
The service unit 33 includes an oscillator unit 34, a multiplication unit 35, a low-pass filter unit 36, and an A / D conversion unit 37. Based on the received GPS signal, the oscillator 34 generates a signal having the same frequency as the frequency of the omega transmission base. Multiplier 35 multiplies the signal generated by oscillator 34 and the received signal to perform synchronous detection. The low-pass filter unit 36 is used to remove noise from the synchronously detected signal. The A / D converter 37 converts the synchronously detected signal into a digital signal in order to transmit data as a digital signal to the client PC 38.
[0040]
A hardware configuration diagram of the client PC 38 and the server PC 39 will be described with reference to FIG. As shown in FIG. 8, the hardware configuration includes an input device 122, a display device 123, a drive device 124, a recording medium 125, an auxiliary storage device 126, and a memory connected to each other via a bus B. A device 127, an arithmetic processing device 120, and an interface device 121 are included.
[0041]
The input device 122 includes a keyboard and a mouse, and is used to input various operation signals to the computer system. The display device 123 displays various windows and data necessary for operating the computer system. The interface device 121 is an interface for connecting the computer system to a network or an interface for connecting other devices, and is composed of, for example, a NIC (Network Interface Card), a modem, a USB (Universal Serial Bus), or the like.
[0042]
A program for operating the computer system is provided by a recording medium 125 such as a CD-ROM or downloaded through a network. In addition, the recording medium 125 is set in the drive device 124, and data and programs are installed in the auxiliary storage device 126 from the recording medium 125 via the drive device 124.
[0043]
The auxiliary storage device 126 stores data and programs, and stores necessary files and the like. The memory device 127 reads the program from the auxiliary storage device 126 and stores it when the computer system is activated. The arithmetic processing unit 120 executes processing according to a program read and stored in the memory device 127.
[0044]
FIG. 9 shows the principle of measuring the phase and amplitude in the earthquake prediction system described above. The description of the same reference numerals as those described in FIG. 7 is omitted.
[0045]
FIG. 9 shows an antenna 30, a preamplifier 32,

multipliers

35a and 35b, an oscillator 34, low-pass filters 36a and 36b, a 90 ° phase shifter 44, and a measuring instrument 47. .
[0046]
The

multipliers

35a and 35b and the 90 ° phase shifter 44 correspond to the multiplication unit 35 in FIG. The low-pass filters 36a and 36b correspond to the low-pass filter unit 36 in FIG. The measuring device 47 corresponds to the client PC 38.
[0047]
In FIG. 9, the signal (in phase) from the reference oscillator 34 is multiplied by the VLF / LF signal received by the antenna, and the signal from the oscillator 34 is shifted by 90 degrees (orthogonal) and VLF. The / LF signal is multiplied to detect the in-phase component and quadrature component of the received signal. The measuring device 47 obtains the phase θ and the amplitude A of the received signal based on these in-phase component and quadrature component as shown in the figure.
[0048]
Next, the software configuration of the earthquake prediction system in the present embodiment will be described with reference to FIG. 10 includes a phase difference acquisition unit 201, a minimum phase difference acquisition unit 202, an amplitude acquisition unit 203, a minimum amplitude acquisition unit 204, a time acquisition unit 205, a difference acquisition unit 206, and a statistical processing unit 207. An abnormality determination unit 208 is shown.
[0049]
The phase difference acquisition unit 201 obtains the phase difference of the VLF / LF signal based on the above-described phase θ and stores the phase difference. The minimum phase difference acquisition unit 202 obtains the minimum value of the phase difference stored in a predetermined time zone. The amplitude acquisition unit 203 is the amplitude of the VLF / LF signal. A Obtain and store the amplitude. The minimum amplitude acquisition unit 204 obtains the minimum value of the amplitude stored in a predetermined time zone. The time acquisition unit 205 obtains the time when the minimum value of the phase difference or the amplitude is taken. The difference acquisition unit 206 calculates a difference between the average amplitude of a predetermined period and a certain amplitude, and obtains the square of the difference. The statistical processing unit 207 obtains an average value and a standard deviation in three statistical values of the minimum phase difference time in a predetermined period, the minimum amplitude time in the predetermined period, and the square of the difference in the predetermined period. The abnormality determination unit 208 determines whether or not the absolute value of the difference between the average value of the minimum phase difference times in a predetermined period and the acquired minimum phase difference time exceeds twice the standard deviation. In addition, the abnormality determination unit 208 determines whether or not the absolute value of the difference between the average value of the minimum amplitude time in a predetermined period and the acquired minimum amplitude time exceeds twice the standard deviation. Furthermore, the abnormality determination unit 208 determines whether or not the square of the acquired difference exceeds a value obtained by adding twice the standard deviation to the average value of the square of the difference during a predetermined period.
[0050]
With the earthquake prediction system described above, the transmission signal of a desired omega transmission base (for example, Australia NWC station (19.8 kHz), Beijing Beijing CHI station (22.2 kHz), Hawaii NPM station (21.4 kHz) from the broadband signal described above. )) Is continuously observed. Similarly, the amplitude and phase of a signal from a Japanese JJY station (Fukushima, 40 kHz) may be measured. Note that. A GPS signal is received as a time standard for phase measurement. Data is transmitted from the observation points at multiple points in Japan to the master station (UEC) through a telephone line every day for analysis.
[0051]
The VLF / LF signal (amplitude and phase) observed at each observation point is transmitted to the UEC master station every day as described above. The daily change pattern (see FIG. 6) is obtained for the amplitude and phase. In this embodiment, two types of signal processing methods for detecting a VLF / LF propagation abnormality are adopted. One is a terminator time method using a terminator time, and the other is a method using a nighttime fluctuation spectrum.
[0052]
First, the case where the terminator time method is used will be described. The terminator time is read in the daily change pattern of the VLF signal (for example, amplitude) in FIG. The software automatically reads and gets one a day (ie tm and te). For example, using the daily data of te, the average value and variance (standard deviation) of te are calculated using the data for the past 10 days of the day.
[0053]
And as shown in FIG. 11, the day when the absolute value of the difference between the average value and te exceeds twice the standard deviation σ is defined as VLF / LF propagation abnormality. An algorithm based on this method is developed, and a method of automatically creating FIG. 11 is used in a computer. Furthermore, when an abnormality occurs, return to the original raw data (daily change pattern) and confirm the abnormality of the terminator time. Terminator time is processed in the same manner for both tm and te. In the above description, the phase difference is used, but the amplitude may be used. In this case as well, the terminator time is the time when the amplitude is minimum.
[0054]
In the graph shown in FIG. 11, the horizontal axis indicates the date, and the vertical axis indicates the difference between te and the average value on that date.
[0055]
Processing when this terminator time method is used will be described with reference to the flowchart of FIG. This process is a process executed by the entire earthquake prediction system. In step S101, a VLF / LF signal is received. In step S102, the phase difference (or amplitude) is obtained and stored. In step S103, it is determined whether a predetermined time has elapsed. In step S104, the minimum value of the phase difference (or amplitude) stored within a predetermined time is obtained. In step S105, the time at which the minimum value is obtained is obtained, and is used as the terminator time T. In step S106, an average value and a standard deviation σ of terminator time in a predetermined period are obtained. In step S107, it is determined whether the absolute value of the difference between the average value and T exceeds 2σ. If exceeded, it is determined that there is an abnormality in step S108, and an earthquake is predicted. If not, it is determined that there is no abnormality.
[0056]
Foreseeing the occurrence of an earthquake is displayed on the screen of the server PC, for example, so as to notify the observer.
[0057]
Next, a method using a nighttime fluctuation spectrum will be described. This signal processing method pays attention to the fluctuation of the night time (from 21:00 to 3 o'clock) in the daily change pattern as shown in FIG. The following formula is used for the observed daily change pattern.
[0058]
[Expression 1]

However, k shows the said day and N is the number of days which takes an average value, and it is normally set as N = 10 days. Further, A (t) indicates a daily change pattern of the amplitude of a certain day, <A (t)> indicates an average daily change pattern for the previous N days, and dA (t) indicates the difference.
[0059]
The graph shown in FIG. 13 is an explanatory graph for helping understanding of <A (t)> and dA (t), and plots <A (t)> and dA (t) every hour. The vertical axis indicates the amplitude and the horizontal axis indicates the time. <A (t)> is an average (for example, a moving average of 10 hours). Therefore, each point of <A (t)> plotted in the graph represents an average from the previous 10 hours including the current value. DA (t) represents the difference from the moving average at each time point.
[0060]
As for these processes, a computer algorithm has been developed, and dA (t) is calculated every day, and the square of dA (t) of the time from 21:00 to 3 o'clock is calculated every day to obtain the data of one day. . Then, time series data similar to the terminator time is obtained. The average and variance (standard deviation) are created using the same method as for the terminator time. When the fluctuation amount exceeds the average value + 2σ, it is regarded as abnormal. It is considered a precursor to the earthquake. For the correlation between the fluctuation amount and the earthquake, see Non-Patent Document 1, which is the latest paper by the inventors.
[0061]
FIG. 14 is a graph showing time-series data of the nighttime fluctuation amount. In this graph, the vertical axis indicates the fluctuation amount, and the horizontal axis indicates the date. In this graph, those exceeding the average value + 2σ are determined to be abnormal.
[0062]
Next, the graph of FIG. 15 will be described. This graph shows the change and magnitude of fluctuation amount, the horizontal axis shows the date, the left vertical axis shows the fluctuation amount, and the right vertical axis shows the magnitude. The thick line indicates the average value + 2σ, the thin line indicates the fluctuation amount, and the bar indicates the magnitude. This magnitude is the magnitude of an earthquake that occurred near the path connecting the transmitting station and the observation point (see Non-Patent Document 1). As shown in this figure, it can be seen that fluctuation anomalies occurred before the earthquake.
[0063]
Processing using this fluctuation amount will be described with reference to the flowchart of FIG. In step S201, the amplitude A at time t is obtained and stored. In step S202, the difference between the average amplitude and the amplitude A at t in a predetermined period is taken, and the difference is set as dA, and the square of the difference dA is obtained. In step S203, an average value of the difference dA squared at t in a predetermined period and a standard deviation σ are obtained. In step S204, it is determined whether or not the difference dA square exceeds the average value + 2σ. If it exceeds, it is determined that there is an abnormality in step S205, and the occurrence of an earthquake is predicted. If not, it is determined that there is no abnormality.
[0064]
Foreseeing the occurrence of an earthquake is displayed on the screen of the server PC, for example, so as to notify the observer.
[0065]
Next, an increase in the accuracy of earthquake prediction will be described. By using the terminator time and the fluctuation amount described above in combination, the accuracy of abnormality can be increased. That is, when the terminator time is abnormal and the fluctuation amount is also abnormal, the probability as an earthquake precursor increases remarkably. Therefore, the accuracy of earthquake prediction can be increased by using both.
[0066]
Next, I will explain the time of the earthquake when it occurs. The time delay between the VLF / LF propagation abnormality and the earthquake can also be seen from FIG. However, to make it more quantitative, cross-correlation is used. FIG. 17 shows a cross-correlation between the fluctuation amount and the occurrence of an earthquake in a time series. The vertical axis indicates the correlation coefficient, and the horizontal axis indicates the date when the earthquake occurrence date is 0 day. A cross-correlation peak appears several days before the earthquake, and an earthquake occurs several days after the propagation abnormality. There is almost the same time delay as the terminator time. This time delay is about several days to one week (see Non-Patent Document 1).
[0067]
Next, the determination of the place where an earthquake occurs will be described. First of all, when a VLF / LF propagation abnormality is detected, it can be predicted when it will occur, but with regard to where it occurs, in the analysis so far, the vicinity of the great circle connecting the transmitting station and the observation point (a little more quantitative) Specifically speaking, it is predicted to occur in the fifth Fresnel zone, for example. Therefore, the method used in the present embodiment is a method based on reception of multiple transmitting stations at multiple observation points. In other words, this is a method using the overlap of the fifth Fresnel band around the great circle connecting the transmitting station and the observation point.
[0068]
The current VLF / LF observation point will be described with reference to FIGS. FIG. 18 shows a receiving facility of a VLF transmitting station which is 7 observation points in Japan. The line in the figure shows a great circle connecting the observation point and the VLF transmitting station.
[0069]
FIG. 19 shows a 40 kHz standard radio wave (JJY) transmission facility and its observation points in eastern Fukushima Prefecture.
[0070]
In this way, there are currently seven VLF / LF observation points from north to Maikori (Hokkaido), Chofu, Chiba Tateyama, Shimizu, Kasugai (Nagoya), Maizuru, and Kochi. The receiving station radio waves of the NWC station (Australia), NPM station (Hawaii), and JJY station (Fukushima) are received.
[0071]
In FIG. 18 and FIG. 19, when an earthquake occurs at a point on the dotted path, a reception abnormality is detected at a measurement point on the dotted path.
[0072]
Therefore, it is possible to determine the location of the ionospheric disturbance causing the abnormality, that is, the position of the future earthquake occurrence from the appearance of the abnormality in various paths. The inventors' VLF / LF network has an accuracy of about 100 km against domestic earthquakes.
[0073]
The future earthquake position prediction method using this VLF / LF network was first proposed in this earthquake prediction system.
[0074]
The results for the Niigata Chuetsu earthquake with the earthquake prediction system, transmitting station and observation point explained above are shown. First, the path used at this time is shown in FIG. FIG. 20 shows the path between the observation point and the transmitting station. In this figure, the solid line indicates the fifth Fresnel zone, which is considered to react to earthquakes whose location is included in this Fresnel zone.
[0075]
As shown in FIG. 20, Maikori (HOK) receives the CHI station (China), Chofu (CHO), Chiba Tateyama (CBA) also receives the CHI station, and Kochi (KOC) also receives the CHI station. ing. Further, a path for receiving the JJY station in Kochi is used.
[0076]
In the latter half of 2004, two major earthquakes occurred as shown in the table in FIG. That is, the Kii Peninsula Southeast Offshore Earthquake and the Niigata Chuetsu Earthquake. Both are earthquakes of magnitude over 6.0, but the former is significantly different from the latter in the sea. Previous inventors' research has shown that the ionosphere is extremely sensitive to shallow earthquakes on land.
[0077]
FIG. 22 is a plot of changes in terminator time (tm) for all the paths described above. In this graph, the horizontal axis indicates the date, the left vertical axis indicates te, and the right vertical axis indicates the magnitude. The left vertical axis shows te for the plotted points, but in the case of a bar graph, shows the difference between the average and te.
[0078]
Further, XXX-YYY written on the left vertical axis indicates a graph in a path between XXX and YYY. For example, CHI-HOK is a path between a CHI station and mother and child villages. Therefore, in FIG. 22, the graph corresponding to the paths of the CHI station and the mother-and-child village, the CHI station and Chofu, the CHI station and Chiba Tateyama, the CHI station and Kochi, and the JJY station and Kochi is drawn.
[0079]
The bottom graph of FIG. 22 shows solar activity by solar radio wave flux and a geomagnetic activity index, and represents the degree of geomagnetism around the earth. In this graph, the solid line indicates the solar radio wave flux, and the bar display indicates the geomagnetic activity index (Ap). The right vertical axis corresponds to the solar radio wave flux, the left vertical axis corresponds to the geomagnetic activity index, and the units are nano terrace and flux, respectively.
[0080]
Prior to the Niigata Chuetsu Earthquake, propagation anomalies clearly appeared at CHI stations at the CHO and CBA stations. This is because it is in the immediate vicinity of both Niigata Chuetsu earthquake paths. There is no change in the Hokkaido Maikori-CHI path. Propagation anomalies occur when geomagnetic activity is low, suggesting a relationship with earthquakes.
[0081]
As described above, the earthquake prediction system according to the present embodiment accurately grasps the propagation abnormality due to the earthquake that occurred in the vicinity of the path, and can thereby predict the occurrence of the earthquake.
[0082]
This international application claims priority based on Japanese Patent Application No. 2005-95137 filed on Mar. 29, 2005, and the entire contents of Patent Application No. 2005-95137 are incorporated herein by reference.

Claims

An earthquake prediction method using both the first earthquake prediction method and the second earthquake prediction method,
The first earthquake prediction method includes:
A signal receiving stage for receiving a VLF / LF signal having a predetermined amplitude and phase transmitted from a predetermined transmitting station, and a GPS signal ;
And VLF / LF signal received, obtains a phase difference between the reference and the ing signals generated from the GPS signal, and the phase difference obtaining step of storing the phase difference,
A minimum phase difference acquisition step for obtaining a minimum value of the phase difference stored in a predetermined time zone including the sunrise time or sunset time;
A time acquisition step for obtaining a minimum phase difference time which is a time when the minimum value is taken;
A first earthquake prediction step for predicting the occurrence of an earthquake based on the change in the minimum phase difference time;
The first earthquake prediction stage includes:
A statistical processing step of obtaining an average value and a standard deviation of the minimum phase difference time in a predetermined period;
An abnormality determination step of determining whether an absolute value of a difference between the average value and the minimum phase difference time exceeds twice the standard deviation;
The phase difference obtaining step includes an oscillator that oscillates at a frequency generated based on the GPS signal, a phase shifter that shifts the output of the oscillator by 90 degrees, an output of the oscillator, and a received VLF / LF signal. A first multiplier for multiplying the output, a second multiplier for multiplying the output of the phase shifter and the received VLF / LF signal, and a first low-pass filter provided on the output side of the first multiplier And a circuit having a second low-pass filter provided on the output side of the second multiplier detects the in-phase component and the quadrature component based on the output of the first low-pass filter and the output of the second low-pass filter. Get the phase shift of the received VLF / LF signal,
The first low-pass filter and the second low-pass filter selectively remove space-derived noise and weather-derived noise,
If it is determined that the abnormality determination step exceeds twice the standard deviation, an earthquake is predicted,
The second earthquake prediction method is:
Obtaining an amplitude of the VLF / LF signal at a predetermined time from 21:00 pm to 3:00 am, and storing the amplitude;
An average value of the amplitude over a predetermined period, and a second earthquake prediction stage for predicting the occurrence of an earthquake based on the amplitude,
The second earthquake prediction stage includes:
A difference acquisition step of taking a difference between the average value and the amplitude and obtaining a square of the difference;
A statistical processing step of obtaining a mean value and a standard deviation of the square of the difference in the predetermined period;
An abnormality determination step for determining whether the square of the difference exceeds a value obtained by adding twice the standard deviation to the average value,
If it is determined that the abnormality determination step exceeds a value obtained by adding twice the standard deviation to the average value, an earthquake occurrence is predicted,
An earthquake prediction method characterized by using both the first earthquake prediction method and the second earthquake prediction method.

An earthquake prediction system using both the first earthquake prediction system and the second earthquake prediction system,
The first earthquake prediction system includes:
Signal receiving means for receiving a VLF / LF signal having a predetermined amplitude and phase transmitted from a predetermined transmitting station, and a GPS signal ;
And VLF / LF signal received, obtains a phase difference between the reference and the ing signals generated from the GPS signal, the phase difference acquiring means for storing the phase difference,
Minimum phase difference acquisition means for obtaining a minimum value of phase difference stored in a predetermined time zone including sunrise time or sunset time;
Time acquisition means for obtaining a minimum phase difference time which is a time when the minimum value is taken;
First earthquake prediction means for predicting the occurrence of an earthquake based on the fluctuation of the minimum phase difference time;
The first earthquake prediction means includes:
Statistical processing means for obtaining an average value and a standard deviation of the minimum phase difference time in a predetermined period;
An abnormality determining means for determining whether an absolute value of a difference between the average value and the minimum phase difference time exceeds twice the standard deviation;
The phase difference acquisition means includes an oscillator that oscillates at a frequency generated based on the GPS signal, a phase shifter that shifts the output of the oscillator by 90 degrees, an output of the oscillator, and a received VLF / LF signal. A first multiplier that multiplies the output of the phase shifter and the received VLF / LF signal, and a first multiplier through a first low-pass filter. And a phase shift of the received VLF / LF signal by detecting an in-phase component and a quadrature component based on the output of the second multiplier and the output of the second multiplier through the second low-pass filter,
The first low-pass filter and the second low-pass filter selectively remove space-derived noise and weather-derived noise,
When the abnormality determining means determines that the standard deviation exceeds twice the standard deviation, it predicts the occurrence of an earthquake,
The second earthquake prediction system includes:
An amplitude acquisition means for obtaining an amplitude at a predetermined time from 21:00 pm to 3:00 am of the VLF / LF signal, and storing the amplitude;
An average value of the amplitude over a predetermined period, and second earthquake prediction means for predicting the occurrence of an earthquake based on the amplitude,
The second earthquake prediction means includes:
Difference obtaining means for taking a difference between the average value and the amplitude and obtaining a square of the difference;
Statistical processing means for obtaining an average value and a standard deviation of the square of the difference in the predetermined period;
An abnormality determining means for determining whether the square of the difference exceeds a value obtained by adding twice the standard deviation to the average value;
An earthquake prediction system characterized by predicting the occurrence of an earthquake when the abnormality determining means determines that the average value exceeds a value obtained by adding twice the standard deviation to the average value.

An earthquake prediction program for causing a computer to execute the earthquake prediction method according to claim 1 .

A computer-readable recording medium on which the earthquake prediction program according to claim 3 is recorded.