1.9MB - 高知工科大学

卒業研究報告
インフラサウンドN型波形イベント
自動検出ソフトウェアの開発
報 告 者
学籍番号:1160086
氏名:
指
山本
反町 玲聖
導
教
真行
員
教授
平成 28 年 2 月 10 日
高知工科大学
システム工学群
目次
第1章 序論・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・1
1.1 研究背景・・・・・・・・・・・・・・・・・・・・・・・・・・・・・1
1.2 研究の目的・・・・・・・・・・・・・・・・・・・・・・・・・・・・2
第2章 インフラサウンドのイベント検出・・・・・・・・・・・・・・・・・・・・・3
2.1 インフラサウンドとは・・・・・・・・・・・・・・・・・・・・・・・3
2.2 インフラサウンドによる地球物理学的イベントの計測・・・・・・・・・3
2.3 N 型波形の特徴・・・・・・・・・・・・・・・・・・・・・・・・・・4
2.4 N 型波形イベントの検出方法・・・・・・・・・・・・・・・・・・・・5
第3章 イベント検出ソフトウェア・・・・・・・・・・・・・・・・・・・・・・・・8
3.1 win フォーマットについて・・・・・・・・・・・・・・・・・・・・・8
3.2 ソフトウェアの仕様・・・・・・・・・・・・・・・・・・・・・・・10
第4章 アルゴリズムの検証・・・・・・・・・・・・・・・・・・・・・・・・・・14
4.1 検証方法・・・・・・・・・・・・・・・・・・・・・・・・・・・・14
4.2 研修結果・・・・・・・・・・・・・・・・・・・・・・・・・・・・・15
第5章 考察・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・19
5.1 N 型波形検出精度の評価・・・・・・・・・・・・・・・・・・・・・19
5.2 今後の展望・・・・・・・・・・・・・・・・・・・・・・・・・・・・21
第6章 結論・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・23
謝辞・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・24
参考文献・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・25
2
第1章
序論
1.1 研究背景
一般に人間の可聴域は周波数 20 Hz から 20 kHz とされており、20 Hz 以下の低周波域
圧力波をインフラサウンドと呼ぶ。インフラサウンドは大気中を伝搬する際に空気の粘性
による減衰を受けにくく長距離伝搬できる特性がある。また、火山の噴火や津波などの自然
現象や核実験等により発せられるものであり大気中のリモートセンシング技術のひとつと
して注目されている。
国内でのインフラサウンドに関する研究は、1980 年代に田平誠氏(現 愛知教育大学名誉
教授)による観測によって基礎が築かれた。近年、国際条約に基づく世界的な核実験検知網
として CTBTO (Comprehensive Nuclear Test Ban Treaty Organization)により全地球上
を 60 点でカバーする大規模観測プロジェクトが進行中であり、ホットな研究領域になりつ
つある。日本国内では、日本気象協会により CTBTO 観測網の 1 地点(IS30)が千葉県夷隅市
にて連続運用されているが、その他の観測点が非常に限られるのが現状である。本分野の研
究進展ならびに防災への活用のため、国内インフラサウンド観測網の充実が望まれるが、セ
ンサ 1 台の価格が高価であることもあり、理想的な密度での配置が出来ず、音波の到来方
向の推定や発生源の特定を行える環境が整っていない。
高知工科大学 システム工学群 宇宙地球探査システム研究室(山本研究室)では、先行研究
(西山、2007)においてピエゾ素子を用いた周波数 0.01 Hz から 100 Hz の音波が検出可能な
低コストセンサの開発が行われ、また山田(2009)において非接触光学型センサの開発に成功
している。鈴木(2009)において観測ロケットの打ち上げによって生じるインフラサウンドを
用いてインフラサウンドセンサの性能比較が行われ、その際到来方向を探知する試みが行
われた。センサ間の間隔が狭かったことから高精度で音波の到来方向探知は行えなかった
が、小松(2012)において桜島の噴火を対象とした音波源位置探査を多地点アレイ観測により
行った際は、桜島南岳の火口半径 15 km 以内に音波源位置を推定することに成功した。ま
た本研究室では南極で起こる氷震などの現象をターゲットとしデータを取得するために
2008 年より南極大陸の昭和基地および周辺地域にインフラサウンドセンサを多地点アレイ
配置しデータを長期間にわたって蓄積している。しかし、インフラサウンドの音波源位置推
定を行ったり、防災に活用するための到来方向自動探知の確立には、例えば、桜島噴火など
のイベントが観測された時刻直後での相関処理が必要になり、波形を人間の目視により確
認している現時点ではリアルタイム処理の自動化には至っていない。
1
1.2 研究の目的
本研究では、インフラサウンドの膨大なデータの中から地球物理学的現象の抽出に有用
な部分の選別や到来方向探知を自動化する上で不可欠なイベント検出を行うため、今回は
波形が特徴的な N 型波形イベントにターゲットを絞り、同イベントを自動検出するアルゴ
リズムを構築し、検出精度の検証を目的とする。また今後インフラサウンドの解析を行う際
有用な、イベントを検出した日時やイベントのスペクトグラムを表示しインフラサウンド
を視覚的に確認きる機能を持つ解析ソフトウェアの開発を目指す。
2
第2章
インフラサウンドのイベント検出
2.1 インフラサウンドとは
人間の可聴域は 20 Hz から 20 kHz とされており、可聴域より高い周波数 20 Hz 以上の
超音波をウルトラサウンド(Ultrasound)と呼び、可聴域より低い周波数 20 Hz 以下の可聴
下音をインフラサウンド(Infrasound)と呼ぶ。図 2.1 に示すように電磁波を分類する際、の
可視光より波長の短い紫外線を(Ultra violet)と呼び、波長の長い赤外線を(Infrared)と呼ぶ
ことに似た定義である。
図 2.1 光と音の関係(田平、2005)
地球上で音波が大気中を伝播する際には直流成分に近い周波数領域までの振動が存在し、
鉛直方向へと伝播する音波の場合には重力による影響を受けるため事実上の下限周波数が
存在し、音波遮断周波数と呼ばれる。15 ℃の等温大気を仮定すると、この下限周波数は
0.00321 Hz となる。空気を媒質としたとき音波は、気温や気圧、湿度によって音速が変化
するため、反射や屈折を行うが、インフラサウンドの場合、周波数が音波遮断周波数に近く
なると気体分子が重力や風の影響を受け伝播速度は音速とはかけ離れていく。また音波が
空気中を伝わる際に受ける減衰には、「幾何学的減衰」「非線形減衰」「空気の粘性による減
衰」と呼ばれる 3 種類が存在する。空気の粘性による減衰は音波の周波数に依存し高周波で
あるほど強く現れる。このためインフラサウンドの場合は「空気の粘性による減衰」を受け
づらく理論上 1 kHz の音波が 1 m 伝播する際に受ける減衰は、1 Hz の音波(インフラサウ
ンド)が 1000 km 伝播する際に受ける減衰と等しい。そのため長距離伝播する特性がある。
またインフラサウンドは火山の噴火や津波などの自然現象やロケット発射などの人為的な
爆発的振動の際に発生するためリモートセンシング手段としても注目されている。
2.2 インフラサウンドによる地球物理学的イベントの計測
著者の所属研究室(宇宙地球探査システム研究室)では 2005 年よりインフラサウンドにつ
いて研究が行われており、インフラサウンドの観測や解析手法の開発、インフラサウンドセ
ンサの開発が行われてきた。現在使われているインフラサウンドセンサの主流は図 2.2 に示
3
す絶対圧を測定できる米国 ParoScientific 製の水晶振動子式圧力計 Nano-Baro(6000-16B)
や図 2.2 に示す相対圧を測定する米国 Chaparral Physics Consultant 製のマイクロホン型
インフラサウンドセンサ Model2、Model25 などであるが、いずれも輸入コストを含む価格
が約 80 万円と高額であるため、多地点への配備には低価格なセンサの開発が必要となり、
西山(2007)において図 2.4 に示す膜面ピエゾ素子を用いた圧電型センサが、山田(2009)にお
いて図 2.5 に示す膜面と半導体レーザおよびリニア光検出器(PSD)を用いた非接触型の光学
式センサの開発がそれぞれ行われた。また各センサの性能評価や防災などに貢献できるイ
ンフラサウンドの到来方向探知を視野に入れ、3地点によるロケット打ち上げ時の騒音の
観測が鈴木(2009) によって行われた。さらに小松(2012)において 3 台のセンサによりイン
フラサウンドの到来方向探知を行うアレイ観測を多地点に配置した多地点アレイ観測が行
われ音波源位置の推定が行われ、その際音波の到来方向を自動計算する自動相関機能を持
ったソフトウェアが開発された。
図 2.2ParoScientifi 製センサ
図 2.3 Chaparral Physics Consultant 製センサ
図 2.4 ピエゾ素子を用いて圧電式センサ
2.3 インフラサウンドの観測手法
4
図 2.5 非接触型光学式センサ
本研究室におけるインフラサウンドの観測方法を図 2.6 に示す。インフラサウンドセンサ
はローパスフィルタの役割をする8本のポーラスパイプ(気泡質のパイプ)の接続されたハ
ウジングに入れ屋外に設置される。センサから出力された信号は、白山工業製のデータロガ
ーにより圧縮 win 形式バイナリファイルとして記録する方法か SAYA 製 A/D 変換ボードを
介して観測用パソコンに非圧縮 win 形式バイナリファイルとして記録する方法がある。
図 2.6 インフラサウンド観測の様子
2.4 N 型波形イベントの特徴
火山の噴火や落雷などのように瞬時的かつ爆発的なイベントの際、発せられるインフラ
サウンドには波形が N の形をした特徴を持つ N 型波形と呼ばれる信号があり、ここでは
N 型波形イベントと呼ぶ。図 2.7 に N 型波形イベントの例を示す。爆発的なイベントの
際発せられた音波は大気中を伝播する際、気圧の高い押し波の部分は音速が早く、気圧の
低い引き波の部分は音速が遅くなる。この特徴により波形が N の形となる。
5
図 2.7 N 型波形イベントの例
2.5 N 型波形イベントの検出方法
本研究で検出しようとしている N 型波形イベントは火山の噴火や落雷などの突発的な地
球物理学的規模のイベントにより発生するものであり、検出される波形は周期性がほとん
どなくインパルス的であり、瞬間的に幅広い周波数領域が含まれる特徴がある。そのため
横軸に時間、縦軸に周波数、強度を色で表したスペクトログラムと呼ばれる表示方法によ
って検出信号を表した場合、縦方向の線となって現れる。図 2.8 に桜島の噴火の際に発生
したと考えられるインフラサウンドの N 型波形イベントについて、波形とスペクトログラ
ムを示す。図 2.8 の中央のグラフに示されたスペクトログラムからもイベントが縦方向の
線となって現れていることがわかる。これに対して周期性のある波はスペクトログラム上
で横方向の縞模様となって現れ、ホワイトノイズのような周期性のない信号は不規則な模
様となって現れる。この特徴を利用し N 型波形イベントを電源ノイズや環境ノイズなどと
分類する。
今回はインフラサウンドのデータを短い時間範囲に分割し、各時間ブロックのデータに
対して FFT (Fast Fourier Transform)により周波数領域のデータに変換しその全スペクト
ル成分の平均値を求める。図 2.8 の下のグラフの白色の線で示されているのがスペクトル
成分の平均値である。ある時間ブロック内に N 型波形のイベントを含んでいる場合はこの
値が急激に大きくなる。また N 型波形のイベントではない信号と区別するため、それぞれ
の時間ブロックでのスペクトル成分の平均値を任意のブロック数に対して移動平均フィル
タをかけ、その値でスペクトル成分の平均値を割りその大きさからイベントの有無を判断
6
する。図 2.8 の下のグラフの赤色の線で示されているのがスペクトル成分を平均値の移動
平均した値であり、N 型波形イベントが発生した際には、赤色の線で示された移動平均の
値に対し白色の線で示された該当時間ブロックにおけるスペクトル成分の平均値が大きく
上回っていることがわかる。
図 2.8 N 型波形イベントの波形とスペクトログラム
7
第3章
イベント検出ソフトウェア
3.1 win フォーマットについて
インフラサウンドセンサにより検出された信号はデータロガー、あるいは A/D 変換ボー
トを介して観測用 PC により win フォーマットと呼ばれるバイナリファイルに保存される。
この win フォーマットは、東京大学地震研究所により開発された多チャンネル地震波形デ
ータ用ファイル形式であり主に国内の地震計のデータを収録するために利用されている。
win フォーマットの仕様について以下にまとめる。
・符号付の 32bit 分解能をもつ
・1 秒ブロックごとに可変長の差分方式で記録されサンプリングレート、チャンネル数、サ
ンプリングサイズの変更が可能
・1秒後ごとに時刻が記録される
・65536 までのチャンネル数に対応
win ファイルのファイル構造について図 3.1 のブロック図に示す。
図 3.1 win ファイルの構造
win 形式の 1 秒間ごとに無数の秒ブロックに分かれている。秒ブロックの内部は先頭から、
10 byte の秒ヘッダと測定しているチャンネル数と同じだけのチャンネルブロックからな
8
る。秒ヘッダは先頭の 4 byte がビッグエンディアン方式でその秒ブロックの容量を表し、
残りの 6 byte で時刻情報を表している。この 6 byte は 16 進数で表した場合、10 進数とし
て読める値となり、YY,MM,DD,HH,MM,SS の順に並んでいる。N 型イベント検出ソフト
ウェアのイベントの検出時刻の表示はこのデータを基に算出している。チャンネルブロッ
クの内容は、4 byte チャンネルヘッダとその後に続く可変長の観測データからなる。チャ
ンネルヘッダには、最初の 2 byte にそれぞれのチャンネルに割り振られたチャンネル番号
がビックエンディアン方式により記載され後に続く 4 bit で 1 サンプリングあたりの差分デ
ータの容量が記載される。ここに 0 が記載されている場合は 1 サンプリングあたりの容量
は、4 bit となる。チャンネルヘッダの残りの 12 bit はビックエンディアン方式でサンプリ
ングレートが記載されている。チャンネルブロックのチャンネルヘッダの後に続く観測デ
ータはチャンネルブロックに記載されたサンプリングレートを基にした1秒間の観測デー
タであり最初のサンプルのみ符号付 4 byte のデータで記録されそう後のサンプリルからは、
差分形式によりチャンネルヘッダに記載されたデータ長となる。また 1 サンプルのサイズ
は 4 byte を超えることがなく比圧縮 win ファイルの場合はそのサイズが 4 byte で固定と
なる。観測データは 1 byte を超える場合、すべてビッグエンディアン方式となる。図 3.3 に
バイナリエディタで参照した win ファイルの内容を示す。
図 3.2 win ファイルの内容
9
3.2 ソフトウェアの仕様
インフラサウンドのイベント検出とその統計的研究を視野に入れた自動観測を行ううえ
で、以下のような機能が求められる。
・データロガーや A/D 変換ボードによって記録された圧縮、非圧縮 win ファイルの読み込
み
・長時間にわたるインフラサウンド波形の連続描画
・各時刻における周波数と強度の情報を示すスペクトログラムの描画
・N 型波形イベントの自動検出とデータ抽出および時刻の決定
本研究で開発した N 型イベント検出ソフトウェアは、良好な操作性を持たせるとともに
データを読み込んだ後にもパラメータ変更等を可能とするため GUI (Graphical User
Interface) 使用して IDL (Interactive Data Language)を用いて開発した。
本ソフトウェアは、プログラム起動時にまず解析するべき win ファイルが保存されたデ
ィレクトリを指定する。同ディレクトリ内のファイルが全て読み込まれた後、(1)サンプリ
ングレート、(2)FFT (Fast Fourier Transform)を行う際の窓となるサンプル数、(3)FFT の
窓を重ねる回数(オーバーラップ数)、(4)移動平均フィルタの平均をかける点数(入力値×2+
1 の点数)、(5)イベント検出の基準となるトリガ値の 5 つの値を GUI 上で入力することに
より、N 型波形イベントを検出した時刻、検出したイベント数を表示する、また画面上にイ
ンフラサウンドの波形とスペクトログラム、スペクトルの平均とその移動平均が描画され
る。図 3.3 に本ソフトウェアの GUI 画面を示す。
10
図 3.3 GUI 画面(2012 年 1 月 11 日の 24 時間分尾長事例)
① サンプリングレート
② FFT の窓となるサンプル数
③ FFT の窓を重ねる個数
④ 移動平均フィルタの点数
⑤ イベント検出の基準となるトリガ値
⑥ 検出されたイベント時刻
⑦ イベントを検出した個数
⑧ 検出されたイベント時刻の表示
ソフトウェアの N 型イベント検出までの流れ(フローチャート)を図 3.4 に示す。まず、読
み込んだ波形データは、(2)FFT の窓となる時間範囲に一致するように指定したサンプル数
(図 3.3 の例であれば 1024)や(3)オーバーラップのかけ方のパラメータ設定(図 3.3 の例なら
ば 2)を元に分割する。この分割した各ブロックに順次 FFT をかけることにより、短い時間
範囲におけるスペクトル強度のデータに変換した。スペクトログラムの描画にはこの計算
11
結果を用いた。さらに各時間ブロックにおけるスペクトル強度の平均をとりブロックごと
のスペクトル強度の和に比例した値を求めた。また、電源ノイズなど長期間にわたり混入す
るノイズを省くためにスペクトル強度の平均値を各ブロックに対して移動平均フィルタを
かけ、移動平均フィルタをかけた後の値でスペクトル強度の平均の値を割った値が、最初に
(5)トリガ値を超える条件を、今回作成したアルゴリズムにおける N 型イベント検出の基準
とした。
図 3.4 N 型イベント検
以下に本ソフトウェアの仕様をまとめた。
・東京大学地震研究所 win-format に準拠して作成されたバイナリファイルの読み込み
・ほぼ無制限の長時間インフラサウンド波形の描画
12
・FFT によるスペクトログラムの描画
・N 型波形イベントの検出と検出時刻表示
・リストのクリックによるイベントを検出した時刻でのグラフの拡大
・イベント検出のためのパラメータの変更可能
・描画したグラフの保存
・イベントの検出時間、設定したパラメータのログデータの
13
第4章
アルゴリズムの検証
4.1 検証方法
今回作製した N 型波形イベント検出ソフトウェアの検出精度を小松(2012)により計測さ
れたデータを用いて検証する。ここでは 2011 年 8 月 15 日に高知県香南市で行われた夏
祭りの打ち上げ花火によって発せられたインフラサウンド信号をターゲットに行った観測
時のデータを用いた。図 4.1 に、打ち上げ花火による 2 分間のインフラサウンド波形デー
タを示す。波形の鋭く立ち上がった部分が花火により発生したインフラサウンド信号であ
り、目視にて 11 例のN型波形が明瞭に確認できる。ソフトウェアの検証にはこの波形に
疑似的に生成したノイズを加算したデータを作り、同データに対する自動処理にて検出数
が 11 例となるよう図 3.3 の⑤で示すトリガ値を調整しながら加えるノイズの大きさを増や
していき、ノイズを加算する前の検出時刻と加算後の検出時刻の変化を比較した。この際
加算したノイズの大きさは 11 例のイベント検出時刻が、ノイズを加算していないときの
検出時刻から 11 例すべてが変化してしまうまで、29 通りについて微調整して検証した。
図 4.1 打ち上げ花火によるインフラサウンドの波形
加算したノイズは、図 4.2 に示すように 2009 年 4 月 14 日に南極の昭和基地で観測され
たブリザードと思われる激しい風ノイズのスペクトルに形状が似るよう正規乱数をもとに
パワースペクトルが 1/√fとなるノイズを生成し環境ノイズを模擬した(図 4.3)。
14
図 4.2 昭和基地の環境ノイズのスペクトル
図 4.3 合成したノイズのスペクトル
検証の際、N型波形イベント検出ソフトウェアの各パラメータは、(2)FFT を行う点数を
64(点)とし、(3)オーバーラップは FFT の窓が 32 点ずつ重なるよう 2(分割)とし、(4)スペ
クトル強度の平均値に対する移動平均フィルタの点数は 11(点)となるよう設定した。また
検証に用いた打ち上げ花火によるインフラサウンドを観測したデータは、(1)サンプリング
レートが 200 Hz であった。
4.2 検証結果
加算した疑似ノイズの大きさと検出時刻の関係を表 4.1 に示す。ノイズレベルが 0 のと
きの検出時刻はノイズを加算していないときの検出時刻である。
15
表 4.1 ノイズレベルと検出時刻の関係
noiselevel
0
1000000
2000000
3000000
4000000
5000000
6000000
7000000
8000000
9000000
10000000
11000000
12000000
13000000
14000000
20000000
25000000
30000000
35000000
40000000
45000000
50000000
60000000
70000000
80000000
90000000
100000000
120000000
140000000
170000000
200000000
8.48
8.48
8.48
8.48
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
1.28
1.28
1.28
1.28
1.28
1.28
1.28
1.28
1.28
1.28
1.28
13.28
13.28
13.28
9.92
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.28
13.44
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
9.92
23.20
23.20
23.20
13.28
23.36
23.36
22.24
22.24
22.24
22.24
22.24
22.24
22.40
22.40
22.40
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
22.24
30.24
30.24
30.24
23.20
30.24
30.24
25.12
25.12
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
23.84
検出時刻[second]
39.20 46.08 58.08 62.24
39.20 46.08 58.08 62.24
39.20 46.08 58.08 62.24
30.24 39.20 46.08 58.08
35.36 39.20 46.08 58.08
35.36 39.20 46.08 58.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 35.36 39.20 46.08
30.24 31.84 35.36 39.20
30.24 31.84 35.36 39.20
30.24 31.84 35.20 39.20
25.12 30.24 31.84 35.20
25.12 30.24 31.84 35.20
25.12 30.24 31.84 35.20
25.12 30.24 31.84 35.20
25.12 30.24 31.84 35.04
25.12 25.76 30.24 31.84
25.12 25.76 27.52 30.24
25.12 25.76 27.52 30.24
25.12 25.76 27.52 30.24
25.12 25.76 27.52 30.24
25.12 25.76 27.52 30.24
25.12 25.76 30.24 31.20
25.12 25.76 27.52 31.20
97.28
73.60
73.60
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.24
62.40
35.04
31.84
31.84
31.84
31.84
31.84
31.84
31.84
103.84
103.84
103.84
73.60
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
62.40
35.04
35.04
35.04
35.04
35.04
35.04
35.04
116.16
116.16
116.16
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
103.84
104.00
104.00
73.28
73.28
73.28
73.28
73.28
73.28
73.28
73.28
また、加算したノイズの大きさと、S/N 比、トリガの値、ノイズを加算していない状態
と比べ検出時刻が変化していない個数との関係を表 4.2 に示す。この際の S/N 比は花火の
音波のスペクトル強度の積分と加算したノイズのスペクトル強度の積分値の比を二乗根し
たときの値である。
16
表 4.2 S/N 比と正しい検出時刻の関係
正しい検出時刻数
S/N比
triggerlevel
0 Infinity
1.8000
11 /11
1000000
2.81735
1.4300
10 /11
2000000
1.99217
1.5200
10 /11
3000000
1.62660
1.5300
9 /11
4000000
1.40868
1.4500
7 /11
5000000
1.25996
1.4300
7 /11
6000000
1.15018
1.4084
6 /11
7000000
1.06486
1.4050
6 /11
8000000
0.99608
1.4050
6 /11
9000000
0.93912
1.4200
6 /11
10000000
0.89092
1.4200
6 /11
11000000
0.84946
1.4200
6 /11
12000000
0.81330
1.4400
6 /11
13000000
0.78139
1.4400
6 /11
14000000
0.75297
1.4400
6 /11
20000000
0.62998
1.4200
5 /11
25000000
0.56347
1.4000
5 /11
30000000
0.51438
1.3925
5 /11
35000000
0.47622
1.3910
4 /11
40000000
0.44546
1.3910
3 /11
45000000
0.41999
1.3905
3 /11
50000000
0.39843
1.3900
2 /11
60000000
0.36372
1.3870
1 /11
70000000
0.33674
1.3860
1 /11
80000000
0.31499
1.3830
1 /11
90000000
0.29697
1.3820
1 /11
100000000
0.28174
1.3820
1 /11
120000000
0.25719
1.3820
1 /11
140000000
0.23811
1.3820
1 /11
170000000
0.21608
1.3840
1 /11
200000000
0.19922
1.3820
0 /11
noiselevel
S/N 比と N 型波形イベントを正確に検出できた個数の関係を図 4.4 のグラフに示す。
17
12
正しい検出時刻数
10
8
6
4
2
0
0
0.5
1
1.5
2
2.5
3
S/N比
図 4.4 S/N 比と正しい検出時刻数の関係
図 4.4 から S/N 比が 0.5 を下回ったあたりから N 型波形イベントを正確に検出できた個
数が急激に低下し始め S/N 比 0.2 前後でまったく検出できなくなることがわかる。また、
S/N 比が 1 程度であれば半分の個数は検出できることがわかった。
18
第5章
考察
5.1 N 型波形検出精度の評価
本ソフトウェアによる検出精度を詳細に評価するため 4.2 節に示した花火による 11 例の
N 型波形イベントのデータに対して本ソフトウェアで描画した波形、スペクトログラム各
時間ブロックにおけるスペクトル強度の平均値とその移動平均の図を 3 例表示して説明す
る。
ノイズを加算していないときの N 型イベント検出ソフトウェアで描写した波形とスペク
トログラム、各時間ブロックにおけるスペクトル強度の平均値とその移動平均を図 4.5 に示
す。図中の上に示されたグラフが波形であり、中央に示されたグラフがスペクトログラム、
下に示されたグラフの白い線が各時間ブロックにおけるスペクトル強度の平均、赤い線が
その移動平均である。波形上の赤い縦線が引かれている部分はソフトウェアがイベントを
検出した位置である。
19
図 5.1 ノイズを加算していないときの検出結果
図 5.1 よりイベントがソフトウェアにより正確に検出されていることがわかる。またスペ
クトログラムからも N 型イベントが明確に確認できる。
ノイズレベルを 14000000 まで増やし検出精度がほぼ半分である 6/11 まで落とした時の
N 型波形検出ソフトウェアで描画させた波形をスペクトログラム、各時間ブロックにおけ
るスペクトル強度の平均値とその移動平均を図 5.2 に示す。
図 5.2 ノイズレベル 14000000 のときの検出結果
図 5.2 と図 5.1 を比較すると、比較的振幅の小さな N 型イベントは正確に検出されない
が比較的振幅の大きい N 型イベントは検出できていることが確認できる。また、ノイズの
付加による波形の鈍化によりイベントの検出時刻がずれる例が生じることもわかる(表 4.1
20
参照)。
図 5.3 にノイズレベルを 170000000 まで増やした際の検出結果を示す。
図 5.3 ノイズレベル 170000000 のときの検出結果
このノイズレベルが今回行った検証で 11 例の N 型イベントの中で最後の 1 例(矢印)を正
確に検出できる現問いの値である。図 5.3 の波形からの目視による N 型イベントの確認は
明らかに困難である。
5.2 今後の展望
第5章 1 節による検証結果の図 5.3 に示す波形からソフトウェアがイベントを多く検出
しているは振幅が 0 に近い部分にしているように見える。また同図の下のグラフに示され
21
たスペクトル強度の平均値が、上のグラフに示された波形の絶対値に近い形をしているよ
うに見える。このことから今回の検証で加算したノイズの低周波成分の周期が FFT の窓
となる時間より著しく大きいためオフセット成分のみが強く影響していると考えられる。
これにより移動平均フィルタを書けた値でスペクトル強度の平均値を割ると微小なイベン
トによる値の変化を拾いにくくなり波形の山や谷となっている部分では検出率が低下した
と考える。またオフセット成分の値を無視することや特定の周波数領域に対する重み付け
を行うことにより検出精度を向上できると考えられる。また、N 型波形イベントのスペク
トルは発生したイベント、音波源との距離などにより波形や波長、周波数成分が異なるた
め N 型イベント検出ソフトウェアのパラメータや平均を取るスペクトル強度の重み付けな
どによりイベントの種類の自動分類に応用できると考えられる。
22
第6章
結論
第4章 4.2 節による検証結果よりの本研究で作成したソフトウェア S/N 比が 1 以下の場
合においても N 型波形イベントを半数以上検出できておりインフラサウンドの到来方向探
知を自動化するための判断基準となりえるアルゴリズムと考える。また N 型波形イベント
の目視による確認の難しい S/N 比においても検出できた N 型波形イベントがあることから
膨大な量のデータの中から地球物理学的に有用な現象の抽出を行うことがある程度可能と
考えられる。よって以上のことから本ソフトウェアは、当初の目的を達成出来たと結論づけ
る。
23
謝辞
本研究を行うにあたり、指導教員でありいつも丁寧なご指導、ご鞭撻を下さいました高
知工科大学 システム工学群 山本真行 教授に心から感謝申し上げます。
同じ研究室で共に研究に取り組んだ同輩の斉藤氏、藤津氏、吉永氏、平松氏及び修士課
程 1 年の水本氏、修士課程 2 年の山崎氏、池原氏、河野氏、学部 3 年の皆様ならびに柿並
義宏 助教(現 台湾国立中央大学)に多くのアドバイス頂いたことに感謝いたします。最後
に本研究にかかわるすべての方に心より感謝申し上げます。
24
参考文献
(論文)
・西山好則、新方式インフラサウンドセンサの開発、平成 18 年度高知工科大学卒業研究報
告、2007.
・山田龍樹、インフラサウンドによる微小圧力変動の検出方式検討と食う神経の開発、平成
20 年度高知工科大学卒業研究報告 2009.
・小松孝康、インフラサウンド多地点アレイ観測システムの構築と音波源一の推定、平成 21
年度高知工科大学大学院特別研空報告、2012
・鈴木敏史、ロケット打ち上げにより励起されたインフラサウンドの計測とデータ解析用ソ
フトウェアの開発、平成 20 年度高知工科大学大学院特別研究報告、2009.
(Web サイト)
・東京大学地震研究所、winformat
http://eoc.eri.u-tokyo.ac.jp/cgi-bin/show¥man?winformat、015 年 9 月 13 日参照
・田平誠、インフラサウンドの世界
http://www.senior.aichi-edu.ac.jp/mtahira/IFS/IFS_introduction.htm、2016 年 9 月 10 日
参照
25