忍者ブログ
日々の研究生活のメモ書きなど
数時間前にLIGOから世界初の重力波直接検出が発表されました
検出に関する論文はこれ
Observation of Gravitational Waves from a Binary Black Hole Merger

ダウンロードは有料なので、できない人はLIGO公式ページで落とすこともできます
Observation of Gravitational Waves from a Binary Black Hole Merger





また、詳細だけを知りたいorグラフだけを見たいときはこちらをどうぞ↓
Data release for event GW150914


重力波信号が含まれているデータも公開されていて、「Gravitational-Wave Strain Data」という項目からダウンロードできます
この記事では、その実際のデータをgnuplotでプロットしてみたいと思います

※ 残念ながら、この通りにやっても肉眼では重力波信号はおそらく見えません・・・
よく考えたら当然なんですが・・・
どうやってFIG 1を書いた方はこちらを読んでもらえると分かると思います
https://losc.ligo.org/s/events/GW150914/GW150914_tutorial.html
(今はなんかアクセスできないが、たぶん鯖落ちだと思われ)


Data release for event GW150914
このページにアクセスして、

2016-02-12 8.14.19



ここをクリックすると、こんな感じの画面が表示されるはず

こんなの



3種類のデータが提供されていて、
hdf5 : 使ったことないのでわからん、pythonやMATLAB他で使えるらしい
gwf : フレームファイル形式、読み込むために frame library というパッケージが必要、重力波業界で広く用いられているフォーマット
txt.gz : txtファイルを圧縮したもの
左がHanfordという検出器、右がLiivingstonという検出器のデータ
あと上の表は、データのサンプリングレートが4096Hz
下が16384Hz
おそらく元々は16384Hzだったはずなので、上はダウンサンプルという処理をしたあとのデータだと思われ

では早速ダウンロードしてみる・・・
今回は
・4096Hzサンプリング(上の表)
・32second
のデータを使ってみる


ダウンロードしたtxt.gzを解凍する
gzip -dc L-L1_LOSC_4_V1-1126259446-32.txt.gz > L-L1_LOSC_4_V1-1126259446-32.txt
gzip -dc H-H1_LOSC_4_V1-1126259446-32.txt.gz > H-H1_LOSC_4_V1-1126259446-32.txt
中身を見てみるとわかるが、このままでは重力波信号が含まれている主干渉計信号しかない・・・
つまり、どの値がどの時刻に対応しているかわからない・・・・

幸い「このデータはGPS時刻 1126259446 からのデータで、重力波信号は112625962 に来ました」と書いてあるので、
一番始めの値が時刻はGPS=1126259446に対応してて、そこから1秒間に4096回ずつ、サンプリングされているということ



このデータを使って、時系列のグラフを描いてみる
gnuplotを起動して次のようなコマンドを打つ
set xlabel "GPS[s]"
set ylabel "strain signal"
plot "H-H1_LOSC_4_V1-1126259446-32.txt" w d title "LIGO-Hanford"
こんな感じ(あとで全部一括で行うスクリプトも書いておきます)

H-H1_LOSC_4_V1-1126259446-32L-L1_LOSC_4_V1-1126259446-32



次に、x軸の値を調整する(今はデータの通し番号、上から何行目かが示されている)
set xrange [0:32]
set xlabel "GPS[s] - ${GPS}"
set output "H-H1_LOSC_4_V1-1126259446-32_v1.png"
plot "H-H1_LOSC_4_V1-1126259446-32.txt" u ($0/4096.0:1 w d title "LIGO-Hanford"
こんな感じで $0 をサンプリングレートで割ればよい

H-H1_LOSC_4_V1-1126259446-32_v1L-L1_LOSC_4_V1-1126259446-32_v1



論文によると、今回の重力波の検出器への到来時間差が書いてあった
GW150914 arrived first at L1 and 6.9 +0.5/-0.4 ms later at H1
つまり先にLivingstonに入射してから、約6.9ミリ秒後にH1に入射したと・・・
じゃあこの時間だけずらして、重ねて書いてみる

さらに、重力波が来たとされる時刻 GPS=1126259462 にも印をつけてみる
こんな感じ
time_shifted_t6.9e-3



よくわからないので、ズームしてみる
time_shifted_t6.9e-3_zoom



このデータの中に重力波があるらしいが、重力波信号は雑音信号に埋もれていて見えない・・・

このページ、
Data release for event GW150914 のFIG3を見てもらうとわかるように、いろんなピークが見える
これらは電源やレーザー干渉計の鏡を吊るすことに起因するノイズ(ラインノイズ、他にもたくさんの由来がある)で、
このデータはこれらのノイズがdominantになってしまっている・・・・

重力波信号が入っていないわけではないよ?
ラインノイズというのはある中心周波数をもったsin波、みたいなもなので、つまり周波数領域で見ると1本のピーク見える
今回の重力波信号は、検出器の感度が良い周波数よりもさらに低周波(10Hzよりもっと低い)からISCO(Inner Most Stable Circular Orbit)と呼ばれる周波数程度までの周波数成分をもっている
FIG 1の一番下のスペクトログラム、と呼ばれるグラフを見てもらうとそのことがよくわかる

このグラフは横軸時間 縦軸が周波数なので、
時間が経つに連れて重力波の周波数がどんどん上がっていっているのが見て取れる
その音が人の耳ではどのように聞こえるか? それはご自身で聞いてみてください!(GW150914ページの一番下からダウンロードできる)
こういう風に周波数が時間とともに上がっていくので、chirp(さえずり) signalと呼ばれている

話を戻して、
え、じゃあどのようにしてFIG1は作ったの?という疑問が出てくるがそれはcaptionに書いてある
For visualization, all time series are filtered with a 35–350 Hz band-pass filter to suppress large fluctuations outside the detectors’ most sensitive frequency band, and band-reject filters to remove the strong instrumental spectral lines seen in the Fig. 3 spectra.
つまり、まずデータに35Hz~350Hzのバンドパスフィルターをかける、これによって、重力波があまり含まれていない低周波成分と高周波数成分を落とす
次に、ラインノイズに汚染された周波数帯域を除去するフィルターをかける(notch filter?)
こういう処理をして、やっとFIG 1のように肉眼で重力波信号だけを見て取れるようなグラフになる・・・・

まぁ重力波信号は肉眼で見て探すわけではないので、例えこのように肉眼で見えなくてもご安心くださいw




ちなみに、この重力波は地球にいつ到来したかというと・・・
GW150914 イベントと呼ばれている通り、まず2015年9月14日らしい

さらに時刻はさっきのGPS時刻からわかる
GPS Time Converter
ここの下のボックスに 1126259462 と入力すると・・・・
GPS換算



UTCで、2015年9月14日 9:50:45 らしい
つまり日本時間(JST=UTC+9時間)では 2015年9月14日18:50:45
(ここ間違ってたので修正しました)

twilog で、その時間に自分が何をしていたか調べたみたら・・・
どうやらアニメ「ガッチャマン クラウズ」を見て、岩崎琢氏の作曲センスに感動していたらしい

へぇ〜

この記事の続きを読むに上記のグラフを描くためのスクリプトを添付しておきます
たぶん、txtファイルを解凍したあとに、コピペでグラフが描けるはず







#!/bin/sh

FS=4096 # [Hz]
GPSSTART=1126259446 #[s]
GPS=1126259462 # [s]
TSHIFT=6.9e-3 #[ms]

gnuplot <<EOF
set term png

set ylabel "strain signal"

set output "H-H1_LOSC_4_V1-1126259446-32.png"
plot "H-H1_LOSC_4_V1-1126259446-32.txt" w d title "LIGO-Hanford"
set output "L-L1_LOSC_4_V1-1126259446-32.png"
plot "L-L1_LOSC_4_V1-1126259446-32.txt" w d lc 3 title "LIGO-Livingston"

set xrange [0:32]
set xlabel "GPS[s] - ${GPSSTART}"

set output "H-H1_LOSC_4_V1-1126259446-32_v1.png"
plot "H-H1_LOSC_4_V1-1126259446-32.txt" u (\$0/${FS}):1 w d title "LIGO-Hanford"
set output "L-L1_LOSC_4_V1-1126259446-32_v1.png"
plot "L-L1_LOSC_4_V1-1126259446-32.txt" u (\$0/${FS}):1 w d lc 3 title "LIGO-Livingston"


set arrow nohead from 16, -3e-18 to 16, 8e-19 lw 2 lc 6

set output "time_shifted_t${TSHIFT}.png"
plot [:][-3e-18:8e-19] "H-H1_LOSC_4_V1-1126259446-32.txt" u (\$0/${FS} - ${TSHIFT} ):1 w d title "LIGO-Hanford",\
"L-L1_LOSC_4_V1-1126259446-32.txt" u (\$0/${FS}):1 w d lc 3 title "LIGO-Livingston"

set output "time_shifted_t${TSHIFT}_zoom.png"
plot [15.4:16.6][-3e-18:8e-19] "H-H1_LOSC_4_V1-1126259446-32.txt" u (\$0/${FS} - ${TSHIFT} ):1 w d title "LIGO-Hanford",\
"L-L1_LOSC_4_V1-1126259446-32.txt" u (\$0/${FS}):1 w d lc 3 title "LIGO-Livingston"

EOF



PR
この記事にコメントする
Name
Title
Color
E-Mail
URL
Comment
Password   Vodafone絵文字 i-mode絵文字 Ezweb絵文字
プロフィール
HN:coffee
職業:物理屋(自称)
趣味:映画鑑賞
自己紹介:
#include <stdio.h>
#include "MyProfile.h"

#define TWITTER coffee_pote

#define WISH_LIST
amazonのほしい物リスト
#ifdef RICH_FLAG
// ↑いつも支援いただきありがとうございます m(_ _)m
#endif


int main(void){

printf("\n");
printf("D論・・・? あぁそんな子もいましたね(執筆中)\n");
printf("\n");
printf("猿でもわかるgnuplot を執筆中(こっちの執筆は半年以上何も進んでいない・・・・)\n");
/* 最終更新 2017/07/19 */
return 0;

}
カウンター
ブログ内検索
ツイートするボタン
Flickr

Template "simple02" by Emile*Emilie
忍者ブログ [PR]