JP2024057456A - 震度推定装置、震度推定プログラム及び震度推定方法 - Google Patents

Abstract

【課題】フィルター処理に要する時間を短縮する。【解決手段】震度推定装置１は、加速度観測記録ａ（ｔ）の加速度に対して、気象庁告示第４号による計測震度計算用のフィルター特性を近似する特性をもつＩＩＲフィルターＧ（ｆ）による演算を行う。ＩＩＲフィルターＧ（ｆ）は、低周波数帯域での特性を再現する１つのハイパスフィルターＧＨ（ｆ）と、低周波数帯域より高い中間周波数帯域及び高周波数帯域での特性を再現する１つのローパスフィルターＧＬ（ｆ）とを備える。ハイパスフィルターＧＨ（ｆ）は、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成される。ローパスフィルターＧＬ（ｆ）は、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである。【選択図】図１

Description

本開示は、震度を推定する震度推定装置、震度推定プログラム及び震度推定方法に関する。

地震が頻繁に発生する我が国では、地震が発生すると気象庁等が地震情報を発信し、ラジオ番組の放送中やテレビ番組の放映中であっても地震速報として、音声や字幕で各地での震度、震源及び地震規模（マグニチュード）が発表されている。
しかしながら、その地震情報は震度３以上の場合に最も早く発信される震度速報であっても通常２分程度要している。これは、地震による揺れの強さに関する指標の１つである震度が、加速度センサーによって測定される加速度記録の時刻歴データを一定時間蓄積した上でフーリエ変換し、周波数領域において気象庁が定めるフィルター処理を行った後、それを逆フーリエ変換して得られた時刻歴データから算出されているためである。すなわち、即時的に震度が計測されるのではなく、一旦地震の加速度波形データを一定時間蓄積した上で解析されているため、時間遅れが生じるのである。

上述の震度計測に関しては、官報第１８３１号で告示された気象庁震度階級表（気象庁告示第４号）に記載されている。加速度観測記録から計測震度を計算するときに用いるフィルターの総合的周波数特性をＪ（ｆ）とすると、式（１），（２），（３），（４）のように３つのフィルターＪ_１（ｆ），Ｊ_２（ｆ），Ｊ_３（ｆ）の重ね合わせとして表される。

ここで式（３）の中のｃ_ｋは、式（５）で示すように、定数である。

フィルターＪ_１（ｆ），Ｊ_２（ｆ），Ｊ_３（ｆ）の周波数特性、及び総合特性Ｊ（ｆ）フィルター特性を図１５に示す。フィルターＪ_１（ｆ），Ｊ_２（ｆ），Ｊ_３（ｆ）は、入力信号と出力信号の振幅のみが異なるゼロ位相推移フィルターである。図１５に示すように、フィルターＪ_１（ｆ）は、減衰定数が１以下、遮断周波数が０．５Ｈｚ、遮断周波数より低周波数帯域での減衰率が１．５のハイパスフィルターである。フィルターＪ_２（ｆ）は、全帯域で減衰率が０．５のフィルターである。フィルターＪ_３（ｆ）は、減衰定数が１以下、遮断周波数が１０Ｈｚ、遮断数波数より高周波数帯域での減衰率が６．０のローパスフィルターである。これら３つのフィルターを重ね合わせた総合的フィルターＪ（ｆ）は、低周波数帯域では遮断周波数が０．５Ｈｚ、減衰定数が１以下、減衰率が１．０のハイパスフィルター、遮断周波数０．５Ｈｚと遮断周波数１０Ｈｚとの間の中間周波数帯域では減衰率が０．５のフィルター、高周波数帯域では遮断周波数が１０Ｈｚ、減衰定数が１以下、減衰率が６．５のローパスフィルターとして機能する。低周波数、中間周波数、高周波数帯域、及び全周波数帯域でのフィルターの周波数特性を図１６に示す。

上記特性をもつフィルター処理は、６０秒間の加速度時刻歴データに対して周波数領域で演算するものとされており、フィルター処理され３成分ベクトル合成された時刻歴について積算継続時間が０．３秒となる振幅値を採用することとされている。ここで、積算継続時間が０．３秒となる振幅値とは、時刻歴の振幅がその振幅値以上の値をとる時間を積算した長さが０．３秒となる振幅値のことである。例えば１秒間のサンプル数が１００である時刻歴に対しては値の大きいほうから３０番目となる振幅値となる。

これまでの問題点を整理し課題解決の検証に資するために、１９９６～２０１９年の期間に日本周辺で発生した地震の観測記録を用いた解析により、既往特許文献によるフィルターの問題点、既往の震度換算方法の改良すべき点を抽出し、本開示によるフィルターの性能と震度換算方法について定量的検証を行った。観測記録は防災科学研究所が収集し提供しているＫ－ＮＥＴ及びＫｉＫｎｅｔの３成分加速度データである。解析に用いた地震総数は９３６、観測点数は１，７１８、記録数は９１，０５８である。解析に用いた地震の震央と観測点位置を図１７に示す。

地震が発生したときに即座に有効な警報を発出するために、震度計算のためのフィルター処理を周波数領域ではなく時間領域で行い即時性をもたせる工夫や、フィルター処理され３成分ベクトル合成された時刻歴から震度への変換の際に対象とする区間長を定義された６０秒間ではなく１０秒とする事例や１秒間とする工夫がなされている。

特許文献１では、「計測震度概算装置、それを用いた計測震度概算システム及び計測震度概算方法」という発明が開示されており、周波数領域でのフィルター処理を行う代わりに時間領域でのフィルター処理が提案されている。特許文献１においては、計測震度計算のためのフィルターは、（１）低周波数帯域でのハイパス特性を再現するための１つの１次ＩＩＲフィルター、（２）中間周波数帯域で周波数の平方根に逆比例するゲイン特性を再現するための３つの１次ＩＩＲフィルター、（３）高周波数帯域でのローパス特性をもつ１つの２次ＩＩＲフィルターで構成されている。ＩＩＲは、Infinite Impulse Responseの略である。図１８に、特許文献１によるフィルターの周波数特性と震度計算値の誤差頻度分布を示す。特許文献１によるフィルターでは、（１）低周波数帯域におけるハイパス特性の再現に１次ハイパスＩＩＲフィルターを採用しているため、低周波数帯域でのゲインピーク付近の特性の再現精度は不十分であり誤差が大きいことが問題点であった。また、誤差の頻度分布は概ね誤差０にピークをもつが、誤差の分布が正負非対称で負の誤差分布が正の誤差分布に比べて広がっている。震度換算の際のデータ区間長は定義された６０秒から１０秒に変更する事例が示されている。

特許文献２では、「リアルタイム震度計とそれを用いた震度等の予知方法」という発明が開示されており、周波数領域でのフィルター処理を行う代わりに時間領域でのフィルター処理が提案され、震度変換の際のデータ区間長は定義の６０秒から１秒に変更されている。計測震度計算のためのフィルターは、（１）低周波数帯域でのハイパス特性を再現するための１つの１次ＩＩＲフィルター、（２）中間周波数帯域で周波数の平方根に逆比例するゲイン特性を再現するための１つの１次ローパスＩＩＲフィルター、（３）高周波数帯域でのローパス特性をもつ１つの２次ＩＩＲフィルターで構成されている。この発明におけるフィルターの周波数特性と震度計算値の誤差頻度分布を図１９に示す。特許文献２によるフィルターでは（１）低周波数帯域におけるハイパス特性の再現に１次ハイパスＩＩＲフィルターを採用しているため、低周波数帯域でのゲインピーク付近の特性の再現精度は不十分であり誤差が大きいことが問題点であった。誤差の頻度分布は概ね誤差０にピークをもつが、誤差の分布が正負非対称で負の誤差分布が正の誤差分布に比べて広がっている点は特許文献１によるフィルターと同等である。区間長は１秒に短縮しているが、積算継続時間は定義の０．３秒のままでありこれに起因するＰ波検出には遅延が生じていた。

特許文献３では、「リアルタイム震度計測装置とその方法」という発明が開示されている。周波数領域でのフィルター処理を行う代わりに時間領域でのフィルター処理が提案され、低周波数帯域でのゲインピーク付近の特性再現の改良が行われている。この発明においては、計測震度計算のためのフィルターは、（１）低周波数帯域でのハイパス特性を再現するための１つの２次／１次のＩＩＲフィルター、（２）中間周波数帯域で周波数の平方根に逆比例するゲイン特性を再現するための３つの１次のＩＩＲフィルター、（３）高周波数帯域でのローパス特性をもつ１つの２次ＩＩＲフィルターで構成されている。特許文献３におけるフィルターの周波数特性と震度計算値の誤差頻度分布を図２０に示す。周波数特性の近似精度は十分であり誤差は周波数帯域全般で小さく抑えられているが、（１）低周波数帯域、（２）中間周波数帯域、（３）高周波数帯域のそれぞれで個別に特性の再現が試みられているため、ＩＩＲフィルターの構成は多段にわたる短所がある。また、震度変換の区間長や積算継続時間は特許文献２から変更されていない。

特許文献４では、「計測震度概算装置、それを用いた計測震度概算システム」という発明が開示されており、周波数領域でのフィルター処理を行う代わりに時間領域でのフィルター処理が提案され、低周波数帯域でのゲインピーク付近の特性再現の改良が行われている。この発明においては、計測震度計算のためのフィルターは、（１）低周波数帯域でのハイパス特性と（２）中間周波数帯域で周波数の平方根に逆比例するゲイン特性を再現するために２つの２次ＩＩＲフィルターと１つの２次／２次ＩＩＲフィルター、（３）高周波数帯域でのローパス特性をもつ３つの２次ＩＩＲフィルターで構成されている。ここでは、（１）低周波数帯域での１つの１次ＩＩＲフィルターと（２）中間周波数帯域での３つの１次のＩＩＲフィルターを２つの２次ＩＩＲフィルターで置き換えているが本質は同等である。特許文献４におけるフィルターの周波数特性と震度計算値の誤差頻度分布を図２１に示す。周波数特性の再現精度は十分であり誤差は周波数帯域全般で小さく抑えられているが、（１）低周波数帯域、（２）中間周波数帯域、（３）高周波数帯域のそれぞれで個別に特性の再現が試みられており、（１）低周波数帯域及び（３）高周波数帯域において多段のＩＩＲフィルター構成となっているため、ＩＩＲフィルターの全体構成は多段にわたる短所がある。また、震度換算におけるデータ区間長や積算継続時間は変更されていない。

これまでに提案されたこれらの改良においては、時間領域でのフィルター処理に関して、（１）低周波数帯域で遮断周波数０．５Ｈｚ、減衰率１．０、減衰定数１．０以下のハイパス特性をもつＩＩＲフィルター、（２）中間周波数帯域で周波数の平方根に逆比例するゲイン特性をもつＩＩＲフィルター、（３）高周波数帯域で遮断周波数１０Ｈｚ、減衰率６．５、減衰定数１．０以下のローパス特性をもつＩＩＲフィルター、これらそれぞれに対して個別に再現する方式をとっている。特許文献３による時間領域フィルターはその一例である。低周波数、中間周波数、高周波数帯域それぞれでのフィルター特性の再現の様子を図２２に示す。図２２において、実線は気象庁告示によるフィルター特性、破線は特許文献３のＩＩＲフィルター特性を示す。各周波数帯域毎に特性を近似するフィルターを設定するため、震度計算のためのフィルター特性を再現する時間領域フィルターの構成は精度を保証するためには多段にわたるため、演算時間がかかるとともに大きな計算負荷がかかる。

計測震度の計算において、加速度時刻歴のフィルター演算に加えて、フィルター処理された３成分時刻歴をベクトル合成して震度に変換する手順が示されている。気象庁告示では対象区間長は６０秒で積算継続時間は０．３秒となっている。即時性を実現するために区間長を６０秒から１０秒に変更する事例や１秒に変更する工夫は既に提案されているが、地震による揺れの消長を適切に追跡する立場から適切な区間長についての十分な議論はなされていない。また、積算継続時間は定義された０．３秒を変更することなく用いているため、これに起因して入力地震動時刻歴に対してフィルター出力は０．３秒の遅延が生じ。地震の揺れの即時的検知について問題であった。

特許第４２２９３３７号公報 特許第４４７２７６９号公報 特許第５８１９２２２号公報 特許第５９４６０６７号公報

特許文献１或いは特許文献２に開示された発明においては、周波数領域での演算の代わりに時間領域フィルターを用いることでリアルタイム性を改善することが可能となったものの、フィルター特性の再現精度に改善の余地があり、震度換算の際の区間長や積算継続時間の問題は深く議論されていない。特許文献３或いは特許文献４に開示された発明においては、フィルター特性の再現精度は改善されたものの、フィルター構成が多段にわたる課題があり、また、震度換算の際の区間長や積算継続時間については議論されていない。Ｐ波到達を的確に検知し速やかに揺れの大きさを予測し有効な警報を即時的に発出する観点から、未だ解決されていない課題である。

内陸で地震が発生すると、震源直上の地域ではＰ波到達時刻から最大の揺れが発現するまでの時間は、短い場合わずか２秒程度である。いわゆる直下地震が発生したときに、特に、震源近傍地域に対して即時的に有効な警報を発出するためには、フィルター処理や震度計算等を実行するのにかかる時間を極力短縮する工夫や揺れの消長を適切に追跡する区間長の設定は不可欠なのである。フィルターの多段構成によるフィルター処理時間や震度計算における積算継続時間に起因する遅延等は１秒以内であっても侮れないものである。警報発出の緊急性を踏まえれば、震度計算の際の時間領域でのフィルター処理における再現精度を保証しつつ、構成段数を少なくすることにより計算手順の簡潔化を実現するフィルターの最適化と、フィルター処理され３成分ベクトル合成された時刻歴から震度への変換計算において対象データ区間長と積算継続時間を調整することにより揺れの適切な追跡と遅延時間の解消を実現する震度換算方法の最適化が解決すべき具体的課題である。

本開示は、フィルター処理に要する時間を短縮することを目的とする。

本開示の一態様は、加速度観測記録作成部と、フィルター演算部と、震度算出部とを備える震度推定装置である。
加速度観測記録作成部は、地震動により生じる加速度を検出する地震計から取得した加速度検出信号に基づいて、加速度の時間変化を示す加速度観測記録を作成するように構成される。

フィルター演算部は、加速度観測記録の加速度に対して、震度計算のために予め設定されたフィルター特性を有するＩＩＲフィルターによる演算を行うように構成される。
震度算出部は、ＩＩＲフィルターによる演算が行われた加速度観測記録の加速度を用いて、地震動を発生させた地震の震度を算出するように構成される。

ＩＩＲフィルターは、第１周波数帯域での特性を再現する１つのハイパスフィルターと、第１周波数帯域より高い第２周波数帯域での特性を再現する１つのローパスフィルターとを備える。

１つのハイパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成される。

１つのローパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである。
このように構成された本開示の震度推定装置は、第１周波数帯域における特性を再現する２次／１次フィルターと、第１周波数帯域より高い第２周波数帯域での特性を再現する２次フィルターの２つのみを用いている。このため、本開示の震度推定装置は、震度計算のために用いられるＩＩＲフィルターの次数を従来の７或いは１２から４に低減することができ、フィルター処理に要する時間を短縮することができる。

本開示の別の態様は、コンピュータを、加速度観測記録作成部、フィルター演算部、及び、震度算出部として機能させるための震度推定プログラムである。
ＩＩＲフィルターは、第１周波数帯域での特性を再現する１つのハイパスフィルターと、第１周波数帯域より高い第２周波数帯域での特性を再現する１つのローパスフィルターとを備える。１つのハイパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成される。１つのローパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである。

本開示の震度推定プログラムによって制御されるコンピュータは、本開示の震度推定装置の一部を構成することができ、本開示の震度推定装置と同様の効果を得ることができる。

本開示の更に別の態様は、地震の震度を推定する震度推定装置が実行する震度推定方法である。
本開示の震度推定方法では、震度推定装置が、地震動により生じる加速度を検出する地震計から取得した加速度検出信号に基づいて、加速度の時間変化を示す加速度観測記録を作成する。

本開示の震度推定方法では、震度推定装置が、加速度観測記録の加速度に対して、震度計算のために予め設定されたフィルター特性を有するＩＩＲフィルターによる演算を行う。

本開示の震度推定方法では、ＩＩＲフィルターによる演算が行われた加速度観測記録の加速度を用いて、地震動を発生させた地震の震度を算出する。
ＩＩＲフィルターは、第１周波数帯域での特性を再現する１つのハイパスフィルターと、第１周波数帯域より高い第２周波数帯域での特性を再現する１つのローパスフィルターとを備える。１つのハイパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成される。１つのローパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである。

本開示の震度推定方法は、本開示の震度推定装置にて実行される方法であり、当該方法を実行することで、本開示の震度推定装置と同様の効果を得ることができる。

震度推定装置の構成を示すブロック図である。 本開示によるフィルターの周波数特性を示すグラフである。 本開示による各周波数帯域でのフィルター特性を示すグラフである。 加速度記録とフィルター処理された時刻歴とを示すグラフである。 震度波形の計算例を示すグラフである。 本開示によるフィルターを用いた震度計算の精度を示す図である。 本開示によるフィルターを用いた震度計算の誤差頻度分布を示すヒストグラムである。 本開示によるフィルターを用いた震度計算の誤差の分布を示す図である。 気象庁による計測震度計検定用波形に対する精度を示す図である。 区間長による震度誤差頻度分布の相違を示す図である。 積算継続時間による震度誤差頻度分布の相違を示す図である。 係数調整後の震度誤差頻度分布を示すヒストグラムである。 積算継続時間による即時震度時刻歴の相違を示す図である。 震度推定処理を示すフローチャートである。 フィルターＪ_１（ｆ），Ｊ_２（ｆ），Ｊ_３（ｆ），Ｊ（ｆ）の周波数特性を示すグラフである。 気象庁告示第４号の計測震度計算用フィルターの特性を示すグラフである。 フィルター性能の検証に用いた地震の震央位置と観測点位置とを示す図である。 特許文献１によるフィルターの周波数特性と震度計算値の精度とを示す図である。 特許文献２によるフィルターの周波数特性と震度計算値の精度とを示す図である。 特許文献３によるフィルターの周波数特性と震度計算値の精度とを示す図である。 特許文献４によるフィルターの周波数特性と震度計算値の精度とを示す図である。 特許文献３による各周波数帯域でのフィルター特性を示すグラフである。

以下に本開示の実施形態を図面とともに説明する。
［１］震度推定装置の構成
本実施形態の震度推定装置１は、図１に示すように、処理部２と、入出力部３とを備える。

処理部２は、ＣＰＵ２ａ、ＲＯＭ２ｂ及びＲＡＭ２ｃ等を備えたマイクロコンピュータを中心に構成された電子制御装置である。マイクロコンピュータの各種機能は、ＣＰＵ２ａが非遷移的実体的記録媒体に格納されたプログラムを実行することにより実現される。この例では、ＲＯＭ２ｂが、プログラムを格納した非遷移的実体的記録媒体に該当する。また、このプログラムの実行により、プログラムに対応する方法が実行される。なお、ＣＰＵ２ａが実行する機能の一部または全部を、一つあるいは複数のＩＣ等によりハードウェア的に構成してもよい。また、処理部２を構成するマイクロコンピュータの数は１つでも複数でもよい。

入出力部３は、震度推定装置１の外部と処理部２との間でアナログ信号またはデジタル信号の入出力を行わせるための回路である。入出力部３には、地震計４が接続される。
地震計４は、例えば地表面または地中に設置され、地震計４が設置された観測点における地震動により生じる加速度を、南北方向、東西方向及び上下方向の３成分で検出し、検出した加速度の時間変化を示す３成分加速度検出信号を出力する。

［２］フィルターの最適化
気象庁告示第４号に示されるフィルターの周波数特性を時間領域フィルターで近似する際に、従来は、フィルターＪ_１（ｆ），Ｊ_２（ｆ），Ｊ_３（ｆ）のそれぞれについて個別に近似フィルターを適用してきた。これは、気象庁告示にしたがって、上記（１），（２），（３）の特性それぞれを個別に再現する考え方に依っている。そのため、近似の精度を保つためにフィルター構成は多段にわたり、構成が複雑になり計算時間がかかる問題点があった。

本開示はかかる問題を解決するために行われたものであり、上記（１），（２），（３）それぞれの特性を個別に再現するのではなく、総合的特性の再現という考え方にしたがって、（１）低周波数帯域での特性を再現する１つのハイパスフィルター（Ａ）、及び（２）中間周波数帯域と（３）高周波数帯域での特性を合わせて再現する１つのローパスフィルター（Ｂ）のみにより、震度計算のためのフィルター処理を高精度で行える演算手段を提供する。

周波数領域でのフィルター処理を時間領域での演算に置き換えるとき、フィルターの総合的特性を十分な精度で再現するためには、（ａ）低周波数帯域でのハイパス特性を実現するためのＩＩＲハイパスフィルター、及び（ｂ）高周波数帯域でのローパス特性を実現するためのＩＩＲローパスフィルターは個別に必要となることは明白である。本開示によるフィルター構成は、これ以上は簡略不可能な究極の構成となっている。

（１）低周波数帯域での特性を再現するＩＩＲハイパスフィルターの減衰率は１．０であるが、遮断周波数付近の帯域における特性を忠実に再現するために２次ＩＩＲフィルターと１次ＩＩＲフィルターの比を用いている。また、（２）中間周波数帯域と（３）高周波数帯域での特性を合わせて再現するＩＩＲローパスフィルターでは、減衰定数が１．０以上の２次ＩＩＲフィルターを用いているところが本開示の特長であり、最も少ないフィルター要素で十分な精度を保証することを特徴としている。

フィルターＪ（ｆ）の周波数特性を近似するフィルターを設計するに際して、フィルターＪ（ｆ）の利得は図１５に示されるとおり周波数０．０１Ｈｚ以下、或いは２０Ｈｚ以上の帯域では十分に小さくなること、更に、観測される地震動の加速度スペクトル振幅も周波数０．０１Ｈｚ以下、或いは２０Ｈｚ以上の帯域では小さいことを考慮して、周波数０．０１～２０Ｈｚの帯域でフィルターＪ（ｆ）に最もよく一致するように、フィルターの係数を最小二乗法で決定した。本開示によるフィルターの周波数特性を図２に示す。また、低周波数、及び中間周波数・高周波数帯域でのフィルター特性を図３に示す。図２において、実線は気象庁定義のフィルター特性、破線は本開示のＩＩＲフィルター特性を示す。図３において、実線は気象庁告示によるフィルター特性、破線は本開示のＩＩＲフィルター特性を示す。

図２２に示された特許文献３によるフィルターは低周波数、中間周波数、及び高周波数ごとに特性を近似する方式をとっているが、本開示では低周波数及び高周波数でのフィルターのみで同様の近似精度を得ていることが確認される。再現の精度は特許文献３や特許文献４と同等であるが、必要最小限の２つのフィルターによる構成である本開示によるＩＩＲフィルターは合計次数が４であり、合計次数７或いは１２の特許文献３或いは特許文献４のフィルターより計算時間の短縮とアルゴリズムの簡潔さが実現されている。なお、特許文献３では、低周波数帯域においては２次／１次フィルターを１つ、中間周波数帯域においては１次フィルターを３つ、高周波数帯域においては２次フィルターを１つ用いており、次数の合計は７である。特許文献４では、低周波数帯域においては１次フィルターを１つと２次／１次フィルターを１つ、中間周波数帯域においては１次フィルターを３つ、高周波数帯域においては２次フィルターを３つ用いており、次数の合計は１２である。

気象庁が定義するフィルターＪ（ｆ）＝Ｊ_１（ｆ）Ｊ_２（ｆ）Ｊ_３（ｆ）の周波数特性は、式（６）に示すように、低周波数帯域におけるハイパスフィルターＧ_Ｈ（ｆ）と高周波数帯域におけるローパスフィルターＧ_Ｌ（ｆ）を組み合わせたＩＩＲフィルターＧ（ｆ）により近似される。

低周波数帯域での特性を近似するフィルターＧ_Ｈ（ｆ）は、周波数特性が極大となるおよそ０．６Ｈｚ付近の特性を充分な精度で近似するように、式（７）に示すように、２次ハイパスフィルターと１次ハイパスフィルターの比で構成される。

ここで、Ｈ_１（ｆ；ｆ_Ｈ１，ｇ_Ｈ１）は遮断周波数ｆ_Ｈ１、利得ｇ_Ｈ１の１次ハイパスフィルター、Ｈ_２（ｆ；ｆ_Ｈ２，ｈ_Ｈ２，ｇ_Ｈ２）は遮断周波数ｆ_Ｈ２、減衰定数ｈ_Ｈ２、利得ｇ_Ｈ２の２次ハイパスフィルターであり、式（８），（９）のように表わされる。

中間周波数及び高周波数帯域での特性を近似するフィルターＧ_Ｌ（ｆ）は、式（１０）で表され、減衰定数１．０以上の２次ローパスフィルターとした。

ここで、Ｌ_２（ｆ；ｆ_Ｌ２，ｈ_Ｌ２，ｇ_Ｌ２）は遮断周波数ｆ_Ｌ２、減衰定数ｈ_Ｌ２、利得ｇ_Ｌ２の２次ローパスフィルターで、式（１１）のように表わされる。

ここで、周波数帯域での総合特性が合致するように利得ｇ_Ｌ２を調整した。気象庁が定義するフィルターＪ（ｆ）は高周波数帯域では周波数の６．５乗で減衰するが、ＩＩＲフィルターにおいては高周波数帯域では減衰が実際の次数より過大になることを利用して、ここでは２次ローパスフィルターによる近似で十分であることが確認された。

地震動の時刻歴をｘ（ｔ）、ＩＩＲフィルターＧ_Ｌ（ｆ）の出力時刻歴をｙ_Ｌ（ｔ）とすると、入出力の関係は、式（１１）を時間領域で表現すると式（１２）のように表せる。

次に、フィルターＧ_Ｌ（ｆ）の出力時刻歴ｙ_Ｌ（ｔ）をフィルターＧ_Ｈ（ｆ）の入力時刻歴として、フィルターＧ（ｆ）＝Ｇ_Ｌ（ｆ）Ｇ_Ｈ（ｆ）の出力時刻歴をｙ（ｔ）とすると、ｙ（ｔ）は式（１３）のように表される。

以上により、地震動ｘ（ｔ）から、震度フィルター出力ｙ（ｔ）が得られる。フィルターの合計次数は４である。
図４に、本開示によるフィルターで処理した南北成分時刻歴と気象庁告示による周波数領域のフィルターで処理した南北成分時刻歴を、フィルター処理する前の観測記録とともに示す。図４において、破線は加速度記録の南北成分、実線はフィルター処理された時刻歴を示す。観測記録は、２０００年鳥取県西部地震（Ｍ７．３）のＫｉＫ－ｎｅｔ日野観測点（ＴＴＲＨ０２）における南北成分である。震央距離はおよそ７ｋｍ、計測震度は６．６である。気象庁告示による周波数領域のフィルターはゼロ位相推移フィルターであるため、観測記録とピーク時刻はほぼ一致しているが、本開示による時間領域フィルターは観測記録に対して位相が遅れている。時間領域でフィルター特性を実現する近似フィルターは、周波数領域でスペクトル振幅のみを変更するゼロ位相フィルターとは異なり因果律を満たすフィルターであり、近似フィルターの振幅特性誤差とは別の要因による相違を生じることは注意すべき点である。この相違を除けば本開示による時間領域フィルター処理と気象庁告示による周波数領域フィルターの処理はほぼ一致していることがわかる。

フィルターの精度を検証するために、区間長を６０秒、積算継続時間を０．３秒として、本開示によるフィルターで処理した震度時刻歴と気象庁告示による周波数領域のフィルターで処理した震度時刻歴との比較を２０００年鳥取県西部地震のＫｉＫ－ｎｅｔ日野観測点における例について図５に示す。実線は時間領域のＩＩＲフィルターによる結果、破線は気象庁告示で定義される周波数領域のフィルターによる結果、縦軸のゼロ付近の値をとる曲線は両者の差である。ＩＩＲフィルターを用いて処理された時刻歴から計算される震度時刻歴は気象庁告示に定義される計測震度時刻歴に充分な精度で一致していることが確認される。両者はよく一致している。

本開示によるフィルターは気象庁告示による周波数特性を精度よく再現しているが、計算震度の精度を確認するために、図６に計測震度と計算震度との関係、図７に計算震度の誤差頻度分布を示す。本開示によるフィルターを用いて計算される震度は、広い範囲にわたって計測震度と良い一致を示しており、計算震度の計測震度からの誤差の平均は－０．００１、標準偏差は０．０２７であった。ここで認められる計測震度相当値と計測震度との差は、フィルターの振幅特性の誤差よりも、ゼロ位相フィルターと因果律を満たすフィルターとの相違によるものであり、ＩＩＲフィルターによりフィルター特性を近似する限りにおいて、これ以上誤差は減少しない。特許文献１によるフィルターを用いた計算震度の計測震度からの誤差の平均は０．００８、標準偏差は０．０４７、特許文献２の場合は平均０．００４、標準偏差０．０４７、特許文献３の場合は平均－０．００１、標準偏差０．０２６、特許文献４の場合は平均０．００３、標準偏差０．０２８である。周波数特性の近似精度に加えて、計算震度の再現性からも、本開示によるフィルターは震度を高精度で即時的に計算するための要件を全て備えていることがわかる。

また、計測震度相当値と計測震度との差について、図８に、地震のマグニチュードとの関係、震源の深さとの関係、観測点の震央距離との関係、及び、計測震度との関係を示す。計算震度と計測震度との差は、地震のマグニチュード、震源の深さ、観測点の震央距離、あるいは計測震度の広い範囲にわたって、これらの値とは無相関に平均値ゼロの周りに±０．１程度の範囲内に収まっている。本開示によるフィルターによる処理は、気象庁告示によるフィルターの近似フィルターとして充分な精度と実用性を有することが確認される。

気象庁では計測震度の精度を評価するために検定を行っている。気象庁による検定の際に用いられる入力波形に対して、本開示によるフィルターを用いた計算震度と計測震度とを比較した。図９に結果を示す。ここで試験した波形の全てに対して差は０．１未満であった。このうち他と比較して誤差が大きい１つは、振幅と周波数が頻繁に変化する複雑な波形に対するもので、気象庁告示に定義される周波数領域でのゼロ位相推移フィルターと時間領域での位相推移ＩＩＲフィルターとの差が現れた例である。

［３］震度換算方法の変更
フィルター処理され３成分ベクトル合成された時刻歴を震度に変換する際に、気象庁告示では対象区間長は６０秒で積算継続時間は０．３秒となっている。

直下で地震が発生した場合、地震の揺れの検出は極めて緊急性を要し、積算継続時間による検知の遅延を解消することが望まれる。特に、Ｐ波到達から最大の揺れまでの時間が２秒程度と短い震源近傍の地域に対しては、Ｐ波到達から実際にＰ波が検知されるまでの時間が０．３秒以上短縮される効果は大きい。また、地震による揺れの消長を即時的に追跡することは、揺れの最大値と発現時刻を予測するうえで必要な情報である。

本開示では、揺れの消長を適切に追跡することを可能とするために区間長を０．５秒と設定し、震度計算の遅延時間を回避するために積算遅延時間を撤廃し区間内での最大値を採用し、そのために生じる震度誤差の偏差を解消するために震度変換に用いる係数を新たな係数に置き換えることにより、高精度で遅延の無い地震の揺れの検知と最適な揺れの消長の追跡を実現する手段を提供する。区間長と積算継続時間を変更することにより生じる計算震度の計測震度からの偏差を解消するために、多くのデータ解析に基づいて震度換算における係数を最適に再決定している。

図５に示された震度時刻歴の例からわかるように、震度変換の際のデータ区間長を６０秒とした場合、入力時刻歴としての地震動振幅がこの区間で徐々に小さくなっている場合でも、計算震度の値はその最大値に達してからほぼ６０秒の間その最大値に維持される。このため、震度時刻歴は、揺れの強さの時間推移が単調ではなく時間とともに増減する場合は、地震動の揺れの強さの時間的推移を反映するとこができない。揺れの時間的推移を見るには震度変換の際のデータ区間長を短くする必要があるが、同時に、平滑化と同様の効果を確保するためには一定のデータ区間長が必要である。ここでは、区間長を１０秒、２秒、１秒、及び０．５秒とした場合を検討した。

また、積算継続時間は出力震度時刻歴が入力時刻歴に対してその時間だけ遅れる効果をもたらす。時間遅れを減少あるいは解消するためには積算継続時間の変更が必要である。変更による影響を調べる目的で、積算継続時間を計測震度計算の場合と同様０．３秒とした場合、０．２秒、０．１秒、０．０秒（区間内での最大値を取ることに相当）とした場合を検討した。

震度変換の際の積算継続時間を０．３秒とし、データ区間長を１０秒、２秒、１秒、及び０．５秒とした場合の計算震度の計測震度からの誤差の頻度分布を図１０に示す。用いた地震及び観測点は先に用いた図１７に示されたものと同じである。区間長が短かくなるにしたがって、即時震度と計測震度との誤差の平均は負方向にずれ、誤差の分散はおおきくなる傾向が認められる。

区間長を０．５秒とし、積算継続時間を０．３秒、０．２秒、０．１秒、及び０．０秒とした場合の計算震度の計測震度からの誤差の頻度分布を図１１に示す。用いた地震及び観測点は先に用いた図１７に示されたものと同じである。区間長を０．５秒とした場合、積算継続時間を短くするにしたがって、即時震度と計測震度との誤差の平均は正方向にずれ、誤差の分散はあまり変わらないが、積算継続時間が０．０秒では最も大きい。積算継続時間を０．０秒とした場合、平均値は＋０．３１５である。

震度変換におけるデータ区間長を０．５秒、積算計測時間を０．０秒としたときの計算震度の計測震度からの偏差は、震度変換に用いる定数０．９４を０．６２５と変更すれば解消される。図１２に、区間長を０．５秒、積算継続時間を０．０秒、震度変換時の定数を０．６２５としたときの計算震度と計測震度との誤差の頻度分布を示す。誤差の平均値は０となっており系統的偏差は解消されている。

図１３に、震度変換の際のデータ区間長を０．５秒とし、積算継続時間を０．３秒、或いは０．０秒としたときの即時震度時刻歴の３成分合成時刻歴からの遅延時間を比較した。図に明らかなように、当然ながら、即時震度時刻歴は３成分合成時刻歴に対して、Ｐ波到達直後には積算継続時間に等しい時間遅延が認められる。震源近傍地点において、即時的に揺れの最大値を予測し警報を発出する際には、この遅延時間は望ましいものではない。図１３における積算継続時間０．３秒と積算継続時間０．０秒との比較から、計算震度の計測震度からの誤差は積算継続時間を０．０秒（区間内の最大値）とした場合は積算継続時間を０．３秒のままとした場合より大きいが、遅延時間の解消は、誤差の増大を考慮しても即時的警報にためには有利である。

ここで検討した震度変換の際のデータ区間長の範囲においては、即時震度と計測震度との差の標準偏差は積算継続時間が０．１秒以上では０．１以下、積算継続時間が０．０秒のときはおよそ０．３であり、データ区間長を０．５秒、積算継続時間を０．０秒としても、十分な精度と実用性を有することがわかる。

ここで採用した即時的評価手法による計算震度は、最大速度値、最大加速度値、あるいはそれらの積を用いる推定値よりばらつきが小さく、地震規模、震源深さ、震源からの距離、計測震度によらない安定した値である。演算に用いた時間領域のフィルターは、気象庁告示によるフィルター特性を必要最小限の段数で精度よく近似することが確かめられる。

［４］処理部において実行される処理
震度推定装置１の処理部２が実行する震度推定処理の手順を説明する。震度推定処理は、処理部２の動作中において繰り返し実行される処理である。

震度推定処理が実行されると、処理部２のＣＰＵ２ａは、図１４に示すように、まずＳ１０にて、地震計４から取得した３成分加速度検出信号に基づいて、南北方向加速度成分、東西方向加速度成分及び上下方向加速度成分それぞれの時間変化を示す加速度観測記録ａ（ｔ）を作成する。

次にＣＰＵ２ａは、Ｓ２０にて、加速度観測記録ａ（ｔ）の各成分ａ_ｃ（ｔ）に対して、上記の式（６）に示すＩＩＲフィルターＧ（ｆ）による演算を行う。
この演算による効果は、ａ_ｃ（ｔ）のフーリエ変換Ａ_ｃ（ｆ）に対してＧ（ｆ）を乗じてＦ_ｃ（ｆ）を得てＦ_ｃ（ｆ）を逆フーリエ変換してｆ_ｃ（ｔ）を得ることと同値であり、式（１４），（１５），（１６）のように表すことができる。

次にＣＰＵ２ａは、Ｓ３０にて、式（１７）に示すように、ｆ_ｃ（ｔ）の３成分ベクトル合成波ｖ（ｔ）を算出する。

次にＣＰＵ２ａは、Ｓ４０にて、式（１８）に示すように、ｖ（ｔ）の常用対数を２倍して定数項Ｃを加え関数ｗ（ｔ）を計算する。

次にＣＰＵ２ａは、Ｓ５０にて、式（１９）に示すように、ｗ（ｔ）に対して、区間長０．５秒について積算継続時間を０（区間内の最大値を取ることに相当）とした演算Ｓを実施して即時震度時刻歴Ｉ（ｔ）を算出し、算出した即時震度時刻歴Ｉ（ｔ）を震度推定値とし、震度推定処理を終了する。例えば、「区間長が１秒であり、１秒間のサンプル数が１００である場合において、積算継続時間が０．３秒となる値」とは、区間長内における関数ｗ（Ｒ，ｔ）の値を大きい順に並べて３０番目となる値である。従って、例えば、「区間長が１秒であり、１秒間のサンプル数が１００である場合において、積算継続時間が０秒となる値」とは、区間長内における関数ｗ（Ｒ，ｔ）の値を大きい順に並べて１番目となる値、すなわち最大値である。

［５］効果
このように構成された震度推定装置１は、地震動により生じる加速度を検出する地震計４から取得した３成分加速度検出信号に基づいて、加速度の時間変化を示す加速度観測記録ａ（ｔ）を作成する。

震度推定装置１は、加速度観測記録ａ（ｔ）の加速度に対して、気象庁告示第４号による計測震度計算用のフィルター特性を近似する特性をもつＩＩＲフィルターＧ（ｆ）による演算を行う。

震度推定装置１は、ＩＩＲフィルターＧ（ｆ）による演算が行われた加速度観測記録ａ（ｔ）の加速度を用いて、地震動を発生させた地震の震度を算出する。
ＩＩＲフィルターＧ（ｆ）は、低周波数帯域での特性を再現する１つのハイパスフィルターＧ_Ｈ（ｆ）と、低周波数帯域より高い中間周波数帯域及び高周波数帯域での特性を再現する１つのローパスフィルターＧ_Ｌ（ｆ）とを備える。

１つのハイパスフィルターＧ_Ｈ（ｆ）は、遮断周波数ｆ_Ｈ２及び減衰定数ｈ_Ｈ２の２つの自由度を有する２次ハイパスフィルターＨ_２（ｆ；ｆ_Ｈ２，ｈ_Ｈ２，ｇ_Ｈ２）と、遮断周波数ｆ_Ｈ１の１つの自由度を有する１次ハイパスフィルターＨ_１（ｆ；ｆ_Ｈ１，ｇ_Ｈ１）との比で構成される。

１つのローパスフィルターＧ_Ｌ（ｆ）は、遮断周波数ｆ_Ｌ２及び減衰定数ｈ_Ｌ２の２つの自由度を有する２次ローパスフィルターＬ_２（ｆ；ｆ_Ｌ２，ｈ_Ｌ２，ｇ_Ｌ２）である。
このような震度推定装置１は、低周波数帯域における特性を再現する２次／１次フィルターと、低周波数帯域より高い中間周波数帯域及び高周波数帯域での特性を再現する２次フィルターの２つのみを用いている。このため、震度推定装置１は、震度計算のために用いられるＩＩＲフィルターの次数を従来の７或いは１２から４に低減することができ、フィルター処理に要する時間を短縮することができる。

また震度推定装置１は、ＩＩＲフィルターの次数の低減によって、フィルター処理のアルゴリズムを簡略化することができる。
また地震計４は、地震動により生じる加速度について、３次元直交座標系を構成する南北方向軸、東西方向軸及び上下方向軸のそれぞれの方向に沿った南北方向加速度成分、東西方向加速度成分及び上下方向加速度成分を検出する。

加速度観測記録ａ（ｔ）は、南北方向加速度成分、東西方向加速度成分及び上下方向加速度成分それぞれの時間変化を示す。
震度推定装置１は、ＩＩＲフィルターによる演算が行われた南北方向加速度成分、東西方向加速度成分及び上下方向加速度成分をベクトル合成して３成分ベクトル合成波ｖ（ｔ）を算出する。

震度推定装置１は、３成分ベクトル合成波ｖ（ｔ）の常用対数を含む関数ｗ（ｔ）を算出する。
震度推定装置１は、関数ｗ（ｔ）に対して、０．５秒に設定された区間長について積算継続時間を０秒とする演算Ｓを実行して即時震度時刻歴Ｉ（ｔ）を算出する。

震度推定装置１は、即時震度時刻歴Ｉ（ｔ）が示す震度を、地震の震度として算出する。
このような震度推定装置１は、区間長を６０秒から０．５秒に短くしたことにより、地震動振幅が増加する様子に加えて振幅が増減する様子を即時的に監視できるようになり、最大震度の予測に有効な即時震度時刻歴を提供することができる。また震度推定装置１は、積算継続時間を０秒としたことにより、計測震度換算に伴う時間遅延がなくなりＰ波初動を即時に検知することができる。

これにより、震度推定装置１は、緊急性を要する震源近傍地域に対して、最大の揺れが襲来する前に揺れの強さを定量的に予測できるようになり、即時的で有効な警報を発出することを可能とする。

以上より、震度推定装置１は、地震防災措置を講じ緊急警報を発出する必要がある施設及び機関において有意義であると同時に、特に緊急性が要請される震源近傍地域における警報発出に対して有効である。

以上説明した実施形態において、Ｓ１０は加速度観測記録作成部としての処理に相当し、Ｓ２０はフィルター演算部としての処理に相当し、Ｓ３０～Ｓ５０は震度算出部としての処理に相当し、低周波数帯域は第１周波数帯域に相当し、中間周波数帯域及び高周波数帯域は第２周波数帯域に相当する。

また、Ｓ３０はベクトル合成波算出部としての処理に相当し、Ｓ４０は対数関数算出部としての処理に相当し、Ｓ５０は即時震度時刻歴算出部としての処理に相当する。
また、南北方向軸は第１軸に相当し、東西方向軸は第２軸に相当し、上下方向軸は第３軸に相当し、南北方向加速度成分は第１加速度成分に相当し、東西方向加速度成分は第２加速度成分に相当し、上下方向加速度成分は第３加速度成分に相当し、３成分加速度検出信号は加速度検出信号に相当する。

また、３成分ベクトル合成波ｖ（ｔ）はベクトル合成波に相当し、関数ｗ（ｔ）は対数関数に相当する。
以上、本開示の一実施形態について説明したが、本開示は上記実施形態に限定されるものではなく、種々変形して実施することができる。

［変形例１］
例えば上記実施形態では、１つの観測点に設置された地震計４から取得した加速度検出信号に基づいて震度を推定する形態を示した。しかし、複数の観測点のそれぞれに設置された複数の地震計４から取得した加速度検出信号に基づいて震度を推定するようにしてもよい。

上記実施形態における１つの構成要素が有する複数の機能を、複数の構成要素によって実現したり、１つの構成要素が有する１つの機能を、複数の構成要素によって実現したりしてもよい。また、複数の構成要素が有する複数の機能を、１つの構成要素によって実現したり、複数の構成要素によって実現される１つの機能を、１つの構成要素によって実現したりしてもよい。また、上記実施形態の構成の一部を省略してもよい。また、上記実施形態の構成の少なくとも一部を、他の上記実施形態の構成に対して付加または置換してもよい。

上述した震度推定装置１の他、当該震度推定装置１を構成要素とするシステム、当該震度推定装置１としてコンピュータを機能させるためのプログラム、このプログラムを記録した半導体メモリ等の非遷移的実体的記録媒体、震度推定方法など、種々の形態で本開示を実現することもできる。

１…震度推定装置、４…地震計

Claims (5)

1. 地震動により生じる加速度を検出する地震計から取得した加速度検出信号に基づいて、前記加速度の時間変化を示す加速度観測記録を作成するように構成された加速度観測記録作成部と、
   前記加速度観測記録の前記加速度に対して、震度計算のために予め設定されたフィルター特性を有するＩＩＲフィルターによる演算を行うように構成されたフィルター演算部と、
   前記ＩＩＲフィルターによる演算が行われた前記加速度観測記録の前記加速度を用いて、前記地震動を発生させた地震の震度を算出するように構成された震度算出部と
   を備え、
   前記ＩＩＲフィルターは、
   第１周波数帯域での特性を再現する１つのハイパスフィルターと、
   前記第１周波数帯域より高い第２周波数帯域での特性を再現する１つのローパスフィルターとを備え、
   前記１つのハイパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成され、
   前記１つのローパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである震度推定装置。
2. 請求項１に記載の震度推定装置であって、
   前記地震計は、前記加速度について、３次元直交座標系を構成する第１軸、第２軸及び第３軸のそれぞれの方向に沿った第１加速度成分、第２加速度成分及び第３加速度成分を検出し、
   前記加速度観測記録は、前記第１加速度成分、前記第２加速度成分及び前記第３加速度成分それぞれの時間変化を示し、
   前記震度推定装置は、
   前記ＩＩＲフィルターによる演算が行われた前記第１加速度成分、前記第２加速度成分及び前記第３加速度成分をベクトル合成してベクトル合成波を算出するように構成されたベクトル合成波算出部と、
   前記ベクトル合成波の対数を含む対数関数を算出するように構成された対数関数算出部と、
   前記対数関数に対して、予め設定された区間長について積算継続時間を所定値とする演算を実行して即時震度時刻歴を算出するように構成された即時震度時刻歴算出部とを備え、
   前記震度算出部は、前記即時震度時刻歴が示す震度を、前記地震の震度として算出し、
   前記区間長は、０．５秒である震度推定装置。
3. 請求項１または請求項２に記載の震度推定装置であって、
   前記地震計は、前記加速度について、３次元直交座標系を構成する第１軸、第２軸及び第３軸のそれぞれの方向に沿った第１加速度成分、第２加速度成分及び第３加速度成分を検出し、
   前記加速度観測記録は、前記第１加速度成分、前記第２加速度成分及び前記第３加速度成分それぞれの時間変化を示し、
   前記震度推定装置は、
   前記ＩＩＲフィルターによる演算が行われた前記第１加速度成分、前記第２加速度成分及び前記第３加速度成分をベクトル合成してベクトル合成波を算出するように構成されたベクトル合成波算出部と、
   前記ベクトル合成波の対数を含む対数関数を算出するように構成された対数関数算出部と、
   前記対数関数に対して、予め設定された区間長について積算継続時間を所定値とする演算を実行して即時震度時刻歴を算出するように構成された即時震度時刻歴算出部とを備え、
   前記震度算出部は、前記即時震度時刻歴が示す震度を、前記地震の震度として算出し、
   前記所定値は、０秒である震度推定装置。
4. コンピュータを、
   地震動により生じる加速度を検出する地震計から取得した加速度検出信号に基づいて、前記加速度の時間変化を示す加速度観測記録を作成するように構成された加速度観測記録作成部、
   前記加速度観測記録の前記加速度に対して、震度計算のために予め設定されたフィルター特性を有するＩＩＲフィルターによる演算を行うように構成されたフィルター演算部、及び、
   前記ＩＩＲフィルターによる演算が行われた前記加速度観測記録の前記加速度を用いて、前記地震動を発生させた地震の震度を算出するように構成された震度算出部
   として機能させるための震度推定プログラムであって、
   前記ＩＩＲフィルターは、
   第１周波数帯域での特性を再現する１つのハイパスフィルターと、
   前記第１周波数帯域より高い第２周波数帯域での特性を再現する１つのローパスフィルターとを備え、
   前記１つのハイパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成され、
   前記１つのローパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである震度推定プログラム。
5. 地震の震度を推定する震度推定装置が実行する震度推定方法であって、
   前記震度推定装置が、地震動により生じる加速度を検出する地震計から取得した加速度検出信号に基づいて、前記加速度の時間変化を示す加速度観測記録を作成し、
   前記震度推定装置が、前記加速度観測記録の前記加速度に対して、震度計算のために予め設定されたフィルター特性を有するＩＩＲフィルターによる演算を行い、
   前記震度推定装置が、前記ＩＩＲフィルターによる演算が行われた前記加速度観測記録の前記加速度を用いて、前記地震動を発生させた地震の震度を算出し、
   前記ＩＩＲフィルターは、
   第１周波数帯域での特性を再現する１つのハイパスフィルターと、
   前記第１周波数帯域より高い第２周波数帯域での特性を再現する１つのローパスフィルターとを備え、
   前記１つのハイパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ハイパスフィルターと、遮断周波数の１つの自由度を有する１次ハイパスフィルターとの比で構成され、
   前記１つのローパスフィルターは、遮断周波数及び減衰定数の２つの自由度を有する２次ローパスフィルターである震度推定方法。