English

回転している中性子星: パルサーのタイミング解析

1. はじめに

1.1 学習課題

すざく衛星搭載の硬X線観測装置(HXD)が観測した代表的な二つのパルサーである 「かにパルサー」と「PSR1509-58」の タイミング解析を通じて、波の位相について学習してみましょう。

1.2 予備知識

1.3 利用する環境

ツール
  • fv
データ
  • SUZAKU/FITS files (教材からリンクされています)

1.4 関連天体

パルサー・中性子星 かにパルサー
パルサー・中性子星 PSR1509-58

1.5 課題作成

作成者 更新日
海老沢 研 2010-02-03

2. 実習 (1) 「かにパルサー」のタイミング解析

かに星雲は1054年に起きた超新星爆発の残骸で、その中心に残ったのが回転している中性子星である 「パルサー」です。


(http://chandra.harvard.edu/photo/2009/crab/crab.jpgより。X線: NASA/CXC/SAO/F.Seward; 可視光: NASA/ESA/ASU/J.Hester & A.Loll; 赤外線: NASA/JPL-Caltech/Univ. Minn./R.Gehrz)

上の図は、赤外線、可視光、X線で撮ったかに星雲の写真を重ねたものです。中心にブルーで見えているのがX線データで、真ん中で明るく光っている点がパルサーです。中性子星は地球のように磁極を持っていて(ただし、磁場の強さは地球の 約一兆倍です)、二つの磁極が見え隠れすることによって、パルサーは灯台のように 点滅しています。


パルサーが自転により周期的に発光する原理(クリックするとアニメーションが見られます。http://www.nasa.gov/mov/283513main_pulsar_512x288.movより。)

日本の「すざく」衛星が実際に2005年9月15日に観測した「かにパルサー」のデータを解析してみましょう。 (ただし、データのうち、この学習に必要ない部分は除いてあります。) 精密な解析より、この観測をしたときのかにパルサーの回転周期は0.03357701秒 であることが分かっています。

crab.fitsをファイルに保存してください (たいていのブラウザでは、右クリックして「別名でリンク先を保存」、または「対象をファイルに保存」を選べば保存できます)。このとき、ファイルの保存先にちょっと注意が必要です。 以下でこのファイルを読むために使う「fv」というプログラムはもともとアメリカ製なので、 日本語名を扱うのが苦手です。crab.fitsが含まれているフォルダーに日本語が含まれないようにしてください。 たとえば、Windowsだったらc://の下に置けば大丈夫です。

crab.fitsは、天文学の世界で標準的に 使われているFITSフォーマットで書かれています。FITSファイルを扱うための"fv"という プログラムで、このファイルを見てみましょう。"fv"のメインメニューから"Open File..." を選び、crab.fitsというファイルを指定します。

ファイルは3つの部分(Extension)に分かれていますが、 "EVENTS" extensionにX線データが入っています。"Events"という 文字の右側の"All"をクリックすると、"EVENTS" extensionの中身を見ることができます。

TIMEという列に、X線光子を検出した時刻(観測の始めが原点)が書かれています。

"Binary Table of crab.fits[1]..."と書いてあるウィンドウで、

  Tools --> Histogram...
を選んでください。  "fv:Histogram"という新しいウィンドウが開いたら、そこでX軸のcolumn Name に"TIME"を指定します(Y軸は空白のままにしておきます)。 Min=0, Max=0.3357701 (10周期), binsize=0.001 と入れて、"Make"をクリックしてみてください。下図のように、X軸は観測開始からの時間、 縦軸は1ビンあたりのX線光子の数という図(ライトカーブ)が表示されます。 (ここで、時間軸を等間隔に区切ったときのひとつひとつの単位を"ビン"と呼びます。 これは英語で、bin="ふたつきの大箱"のことです。)

これで、パルサーが10回、点滅する様子がわかったでしょうか?

とてもこれだけではわからないでしょう。それは何故かというと、 検出器に入ってくる光の数が少ないからです。 binsize=0.001というのは、時間ビンの幅が0.001秒ということで、 グラフの縦軸を見ると、0.001秒のあいだ検出された光子は、0,1,2個の いずれかだということがわかります。こんなに光子数が少なくては、 パルスは見えません。 どうしたらパルスがよく見えるようになるでしょうか?

この観測の開始から終了までの経過時間は、約22,000秒です。その間に、 パルサーは65万回転以上しています。パルサーが一回転する間の 強度曲線を65万回転分重ね合わせれば、65万倍明るいパルサーを一周期 観測したのと同じ効果が得られます。さっそく、それを見てみましょう。

観測の先頭から数えて、 各光子が到達したときのパルス数(以下、Nとします)を数えます。それには各光子の到達時刻、 TIMEを、周期0.03357701で割ってやれば良いです。Nが実数であることに注意してください。  

"Binary Table of crab.fits[1]..."と書いてあるウィンドウに戻り、

  Tools --> Calculator...
を選んでください。  "Calculator"という新しいウィンドウが開いたら、上部の二つのダイアログボックスで、"="の 左にN、右にTIME/0.03357701と入れて、N=TIME/0.03357701で定義される新たなNと言う名前の列を 作成します。右側の"Calculate"をクリックしてください。 "Binary Table of crab.fits[1]..."ウィンドウに"N"という列が追加されたのがわかると思います。 ここでNが実数であることに注意してください。 Nの小数部分はもちろん0から1までの値を取りますが、これが各光子の到達時刻における 位相(phase)です。 各光子の位相は、光子がパルスの中のどの位置で検出されたかを示しています。 位相を0から1までの実数で表すこともありますし、0°から360°、あるいは0ラジアンから2πラジアンまでの 角度で表すこともあります。

次に、Phaseを求めて新たな列にしましょう。

 "Calculator"ウィンドウに戻り、上部の二つのダイアログボックスで Phase=N-floor(N) と言う式を作り、"Calculate"をクリックします。 floor(N)は、Nを超えない最大の整数を表します。よって、N-floor(N)は、Nの小数点以下の部分です。 "Binary Table of crab.fits[1]..."ウィンドウに"Phase"という 列が追加されたこと、NとPhaseの値を比較して、PhaseにはNの小数点以下の値が入っていることを確認してください。

次に、パルスの位相を合わせて、65万回の回転を重ね合わせてみます。  

 

"fv:Histogram"ウィンドウに戻ります(すでに閉じてしまっていたら、 ふたたび、"Binary Table of crab.fits[1]..."と書いてあるウィンドウで、 "Tools --> Histogram..." を選んでください)。  X軸にPhaseを指定してください。 Min=0, Max=1.0, binsize=0.025 と入れて、"Make"をクリックすると、横軸はPhase(位相)、縦軸は1位相ビンあたり のX線光子という図(周期で折り畳んだライトカーブ)が表示されます。これが かにパルサーが一回転する間に私たちが観測するX線の強度変化を示しています。 0.03357701秒の周期を0.025位相ビンごとに分けたから、 一位相ビンの長さは0.00084秒です。縦軸を見ると、 一位相ビンに光子が18000から21000個入っていることがわかります。 中性子星が一回転する間に二つの磁極が見え隠れするので、 二つのピークが存在することがわかります。

上手のような「ヒストグラム」ではなく、折れ線グラフが表示される場合もありますが、データとしては全く同じです。 ヒストグラム表示にするには、グラフ表示ウィンドウのEditメニューから、Preferencesを選択してください。 POW Preferencesという新しいウィンドウが開きます。 Linesというタブを選び、ConnectでNormalの代わりにHistogramを指定します(下図参照)。その後で下部のSaveを クリックします。すると次回からは折れ線グラフではなく、ヒストグラム表示になるはずです。

実習 (2) 「PSR1509-58」のタイミング解析

次に、PSR1509-58という名前のパルサーの解析をしてみましょう。 PSRというのはパルサー(Pulsar)を意味していて、1509と-58というのは それぞれ空の上での経度、緯度を表しています。 すざく衛星が2005年8月23日に観測したPSR1509-58のデータ、 PSR1509-58.fitsをファイルに保存してください

このパルサーの周期は約0.151364秒です。また、観測の開始から終了までの時間は 12万6千秒でした。上記の「かにパルサー」の例のように、データをパルス周期で重ね合わせる事でパルス波形が はっきりと見えてくるはずですが、ここでは、どれだけの精度でパルスの 周期が決定可能かどうかを考えてみましょう。

まず、周期をP、観測開始から終了までの時間をTとすると、その中にパルスはT/P個あることが わかると思います。これをNとしましょう。正しい周期を用いて、N個の パルスを重ね合わせると正しいパルス波形が得られます。 まちがって、Nから少しずれた数のパルスを重ね合わせてしまうと、位相がずれてパルスの波形が壊れてしまいます。 周期がΔPだけ間違っていたとすると、一周期の間に、位相は

ΔP/P
の割合だけずれてしまいます。さらにN回のパルスでは
N ΔP/P
だけ位相がずれます。これが大きすぎると、たくさんのパルスを重ねあわせても正しいパルス波形は求まりません。 この許される位相のずれとして、どの値をとるのが良いか、はっきりと決まっている訳ではありませんが、 半周期分ずれたら波形の崩れがわかることを考える、1/2を取るのが適当でしょう。 したがって、
N ΔP/P = 1/2
と置き、N=T/Pの関係から、
ΔP =    (1/2) P2/T
が得られます(微分を使った説明はこちら)。 今の例では、Pは約0.151364秒、Tは12万6千秒だから、期待される精度ΔPは~9x10-8秒となります(約10-7秒ですね)。 この観測から決定されるパルサーの周期を0.151364?秒と書いたとき、?の桁は決まるが、それ以上 正確には周期は決定できないと言うことになります。この"?"に0から9までの数字をいれてデータを 重ね合わせてみて、正しい周期を求めてみましょう。重ね合わせたパルス波形がもっともきれいにでるものが正しい答えです。以下のようなパルス波形が得られたでしょうか?(正しい周期の答え


パルス周期の決定精度の微分をつかった考え方
すでに微分を習っている人にとっては、以下のように考えた方がわかりやすいかもしれません。

NP = T
から、両辺の微小な変化量は、Tが不変であることを使って、
ΔN P + N ΔP = 0
となることがわかります。これを変形して、
ΔP = - ΔN P/N =  - ΔN P2 /T.
この式が、正しいパルス数からのずれΔNと、正しい周期からのずれΔPの関係を与えます。 ΔNの値として1/2を取ると、
ΔP =   - (1/2) P2/T.
が得られます。上で得られた式と符号は逆ですが、ズレの大きさはの絶対値だけが問題なので、ここでは符号は気にしなくて良いです。
PSR1509-58の正しい周期
正しい周期は0.1513648秒です。?の部分に0から9までの数字を順番に入れて、並べて比較すると違いがよくわかるでしょう。

Last Modified: 2010-2-4