すばる中秋の名月の学校
宿題までのウォーミングアップ
市川幸平
2015
年
8
月
1
はじめに
本書の目的は、「すばる中秋の名月の学校(以下、名月の学校)」の参加者が、模擬研究を行うに あたって必要になるであろう必要最低限のUnix知識を前もって取得しておくことにある。天文学 の研究では大量のデータから必要な情報を持ってきたり、大量の計算を行ったりすることは日常茶 飯事である。このような時に1つずつ人間がファイルを弄ったり計算をしていては終わりが来な い*1し、人間は残念ながらよくミスをする。特定のことを大量にやる、という作業は人間は苦手な のだ。逆に、コンピュータはそういうのが得意で、そのコンピュータに指令を送るのに便利なのが Unixのコマンドである。最終的にコマンドの塊をスクリプトやプログラムとしてまとめておけば、 その部分にミスがないかどうかを確認すればよく、ミスが起きていないかの確認作業も場所が絞ら れてグッと楽になる。このように、Unixコマンドは研究の効率を大幅にあげてくれる要素がふん だんに盛り込まれていて、研究の大きな友となってくれるわけだ。その友に出会える第一歩が本書 の目指すところである。すでにUnixのコマンド操作やPerl, Python*2などに慣れ親しんでいる人であれば、ほとんどの
ページをすっ飛ばして、最後の宿題の部分だけ取り組んでもらっても構わない。ただ、はじめて Unixを触る人であれば、ここで紹介しているコマンドを一緒に動かしてみよう。このページを最 後まで読んで、一緒にコマンドを実行した人であれば、最終ページにある宿題はそこまで難しくは ならないように設定してある。
2
用意するもの
この文書を読み進めるにあたって、まず準備してほしいものを羅列してみた。 • PC • ターミナル • gnuplot *1 100万天体に対して同じ計算をすることにして、人間であれば一回の計算に10秒かかるとしよう。それだけで 106× 10 = 107秒かかる。1日のうち人間は1/3くらいしか働けないし疲れてくるので、3× 107∼ 107.5秒くら いはかかるだろう。1日は24× 3600 = (30 − 6)(30 + 6) × 100 = 86400秒しかないし、1年は107.5秒くらいし かない。 *2Perl, Pythonともにプログラミング言語の名前。両言語とも、天文学で良く使われるプログラミング言語。最近は Pythonの方が人気な印象。• すばる春の学校「Unixのど基礎」を読む*3 今回、名月の学校に参加する皆さんであれば、自分のPCは所有しているだろうという仮定して も問題ないはず。なので、一つ目はすでにクリアである。また、これから行う作業は、2つ目に挙 げているターミナル上で行うものがほとんどなので、ターミナルは準備が必須である。とは言うも のの、Linux, Macの場合は、すでにターミナルはインストールされているので、ターミナルを開 いておこう*4。Windowsの場合は、少し準備が必要。よく行われているのは、Cygwinをインス トールして対処する方法。これによって、ターミナルと同等の作業を容易に行うことができるよう になる。「Cygwin インストール」 などで検索してみれば、方法はいくつでも見つかるはず。その 場合は、後述するgnuplotもCygwinと同時にインストールしておくと楽である。もしどうしても
インストールがうまくいかない場合は“unix emulator online” などと検索すると、ターミナルっ
ぽいものをウェブ上で使用することができるサイトがいくつか見つかるはず。
gnuplotは今回使用するプロッタ*5である。matplotlib*6が好きだったり、SuperMongo*7 が好 きな人は、もちろんそれらを使用しても構わない。が、特に今までこだわりが無かった人は、まず はgnuplotから始めてみよう。満足がいかない点が出てきたあたりで、他のプロッタも使用してみ て比較する。そうやって自分の好きなプロッタに今後出会っていけばよい。ちなみに、gnuplotの インストールもそんなに難しくない。各自のOSに合わせてインストールを行うこと*8。 最後の「Unixのど基礎」で解説されたことは、今後はすでに既知のものとする。特にUnixを初 めて使う人は、「Unixのど基礎」を読んでいなければ、「3章からいきなり難易度が上がったな!」 と感じることだろう。必ず一度は目を通し、自分で手を動かしておくこと。ターミナルがあなたの 目の前で今開かれているなら、時間はそんなにかからないだろう。「Unixの基礎」で書かれている ディレクトリの概念は特に大事だし、ファイルを開いたりするエディタとしてemacsも紹介され ている*9。これらの言葉や使い方を知っているだけで、今後の解説の理解度がグッと変わってくる はずだ。
3
ターミナルを開いて遊んでみる
3.1
ファイルのダウンロード
では、まずはターミナルを開いて少し遊んでみる。ターミナルを開いたら、おそらくは自分の ホームディレクトリにいるはず。自分が今どのディレクトリにいるのかの確認は、pwdコマン *3 http://www.naoj.org/Observing/DataReduction/mtk/spring12/references/unix_class.pdf Unixのど基礎に何度か出てくる別冊は http://www.naoj.org/Observing/DataReduction/mtk/spring12/references/linuxbasic10s.pdf に あ る。 *4もしターミナルがどこにあるかわからない場合はウェブ等で検索してみよう *5図を描画するソフト/ライブラリのこと。 *6天文業界で良く使われるプロッタその1。Pythonパッケージの一つであり、Pythonユーザであれば違和感なく習 得が可能である。 *7天文業界で良く使われるプロッタその2。 *8これも、わからないければ自分で調べてみたり、近くにいる詳しそうな友人や先輩に質問してみよう。 *9世の中には数多くのエディタが存在するが、少なくとも一つは使えるようにしておこう。そこで不満が出始めたら別 のエディタと比較検討して乗り換えていこう。emacs, vi, geditなどは周りにもユーザーがいることだろう。大事な のは、彼らが言うことを鵜呑みにするのではなく、自分にマッチしているかどうかを自分で判断するという姿勢であ る。そうでないと、気がつかないうちに不毛なエディタ戦争に巻き込まれてしまっているかもしれない。ド*10で可能だ。僕の場合、 $ pwd /home/ichikawa /home/ichikawaにいる。このディレクトリ名は各自異なるはずだが、自分の好きなディレクトリ 下で作業をしてほしい。僕はこの/home/ichikawaディレクトリで作業をすることにする。さて、 自分のいるディレクトリ*11が確認できたところで、 $ wget http://www.kusastro.kyoto-u.ac.jp/~ichikawa/Autumn2015.tar.gz と打ってみよう。データのダウンロードが始まるはず。ダウンロードが済んだのを、lsコマン ド*12で確認しよう。 $ ls Autumn2015.tar.gz lsコマンドは特定のファイルやディレクトリを探すのにも便利である。 $ ls Autumn2015.tar.gz Autumn2015.tar.gz とファイル名が返ってくれば、カレントディレクトリに“Autumn2015.tar.gz” があるという意味 になる。 ところで、さきほどのコマンドのファイル名、まさか全部手でタイプはしていないだろうか。 $ ls A まで打った時点で、一度タブキーを押してみよう。すると、 $ ls Autumn2015.tar.gz と勝手にファイル名を補完してくれる。Unixコマンドは、タブキーを押すと、ユニーク性が保た れる部分までは勝手にファイル名などの補完を行ってくれる。人間がタイプする手間を省いてくれ る優しさあふれる設計をしてくれているのだ。これは、ディレクトリの移動をする際や、長いファ イル名などを打ちたいときなどに大活躍する。 さて、話を元に戻そう。ダウンロードしたファイルは“.tar.gz”という拡張子を持っていて、こ のままではファイルの中身を見ることができない。解凍をするにはtarコマンドを使う。 $ tar zxvf Autumn2015.tar.gz とすれば、解凍がはじまる*13。lsしてみると、新しいディレクトリが表示されるはず。 $ ls Autumn2015/ Autumn2015.tar.gz
*10print working directoryの略。
*11自分が今いるディレクトリをカレントディレクトリという。
*12listの略。
*13きちんとtar zxvf Aまでタイプところでタブキーを押しただろうか。もしそうでないなら、君は少なくとも数秒は
もしwgetを打ってもエラーが返ってくるようだったら、その場合はwgetがインストールされ ていないかもしれない。そのチェックにはwhichコマンドを使ってみよう。 $ which wget /usr/bin/wget などとパス*14が返ってくれば、wgetはきちんとインストールされている。wgetコマンドはきち んと動くはずだ。もう一度URLなどが間違っていないかチェックしよう。一方で、whichコマン ドを打って $ which wget $ などと何も返事が返って来ない場合は、wgetコマンドはインストールされていない*15。wgetコ マンドを自分でインストールするか、直接URLを打ち込んでダウンロードを行い、カレントディ レクトリまでファイルを移動しよう。そのためにはmvコマンドが便利だ。あなたが指定している ダウンロードディレクトリが”YourDownloadDirectory”という名前だとすると、そのディレクト リからカレントディレクトリ“.”までmvしたいので、 $ mv YourDownloadDirectory/Autumn2015.tar.gz . $ tar zxvf Autumn2015.tar.gz 最期に“.” が書かれているのを見逃さないようにしよう。tar で解凍して現れるディレクトリ “Autumn2015”内に、今回のウォーミングアップで使用するデータをすべて入れておいた。もし 以下で説明するダウンロードがうまくいかなかった場合は、このファイルの中に入っているものを 適宜cpして使用しよう。
3.2
今回扱うデータのダウンロード
さて、今回のウォーミングアップでは、実際の観測サーベイのカタログデータを利用してみる。使用するカタログは、Swift /BAT 70ヶ月カタログ(以後、BATカタログと言う) で、これは硬X
線*16(E = 14− 195 keV) で全天を70ヶ月間サーベイして受かった*17天体をまとめている。
まずは、BAT カタログのウェブサイト*18 に行き、カタログを見てみよう。図 1のようなカ
タログを紹介するウェブサイトにたどり着くはずだ。今回はテキストデータを扱うので、“BAT
70-month Survey Products” と書かれているところの“Catalog: ASCII (∼300 KB)” の部分を クリックしてみよう。すると、よくわからないテキスト文字列の羅列が表示されるだろう。よく
見てみると、これは、区切り文字*19 が | (バー) で、左から順に天体の番号、天体名 (SWIFT
*14/usr/bin/wgetはwgetがどこに存在するかを表しているファイルの住所のようなもの。このようなものをパスと 言う。詳細を知りたければ検索してみよう。
*15あるいは、wgetコマンドが存在する場所にパスが通っていない。
*16硬X線は“hard X-ray”の日本語訳。天文業界では、X線のうち、高エネルギー側を“hard (硬)” X線、低エネ ルギー側を“soft (軟)” X線と呼ぶ。soft/hardの境界となるエネルギーは必ずしも決まっていないが、感覚的には E = 2 keVを境にsoftからhardと言う人が増える。
*17“detect”の和訳として、「受かる」という言葉を天文学者はよく使う。たとえば、「NGC 1068でHα輝線が受かっ た」などのように使う。天文学者以外ではあまり使わないかもしれない。
*18http://swift.gsfc.nasa.gov/results/bs70mon/ *19セパレータと言ったりもする。
図1 Swift/BAT 70ヶ月カタログに載っている天体の種類と位置。わからない単語は自分で 調べてみよう。http://swift.gsfc.nasa.gov/results/bs70mon/より取得。 JXXXX.X-XXXXで表示されているもの)、赤経(RA)、赤緯(DEC)、· · · と並んでいるのがわか る。つまり、一列ごとに各天体ごとの情報がずらりと並んでいるカタログであることがわかってく る。これが今回、我々がウォーミングアップに使用するカタログデータである。 では、次にこれを自分のPCにダウンロードしてみよう。ダウンロードするには、Macであれ ば“Cmd+s”であったり、コピペする方法でもよいのだけれども*20、今回はターミナルでコマン ドを使ってダウンロードを試してみる。すでに試したwgetコマンドを使えば*21、ターミナル上で $wget http://swift.gsfc.nasa.gov/results/bs70mon/inc/data/BAT_70m_catalog_20nov2012.txt
とすると、今自分がいるディレクトリに“BAT 70m catalog 20nov2012.txt”というファイルが存
在しているはず。ファイル確認の方法はlsコマンドだった。
$ls
Autumn2015/ Autumn2015.tar.gz BAT_70m_catalog_20nov2012.dat
どうせなので、BATカタログデータを新しいディレクトリを作って、その下におくことにしよう。
新しいディレクトリの作成はmkdirコマンド*22を使う。
$ mkdir BAT70mon $ ls
Autumn2015/ Autumn2015.tar.gz BAT70mon/ BAT_70m_catalog_20nov2012.dat
*20あまりおすすめしないけれども。
*21もしwgetがうまくいかなかった場合は、Autumn2015/の下に同様のファイルを置いてあるので、自分でcpして くること。
ファイルを移動しよう。ファイルの移動はmvコマンド*23が使える。ファイル名が長いので、タ ブキーを使って、補完機能を使っていこう。
$ mv BAT_70m_catalog_20nov2012.dat BAT70mon/ $ ls
Autumn2015/ Autumn2015.tar.gz BAT70mon/
lsしてみると、BAT カタログデータが表示されていない。これで、BAT カタログデータは、 “BAT70mon”ディレクトリの下に移動したことが確認できた。さて、我々もBAT70monディレ クトリの下に行きたい。ディレクトリの移動は“cd”コマンド*24で可能だ。 $ cd BAT70mon/ $ pwd /home/ichikawa/BAT70mon $ ls BAT_70m_catalog_20nov2012.dat BAT70monの下にきちんとBATカタログデータが入っていることがわかる。
3.3
カタログを少し見てみる
さて、BATカタログの中身を見てみよう。データの中身を見る場合は“less”コマンドが使える。 $ less BAT_70m_catalog_20nov2012.dat| BAT_NAME| RA| DEC| SNR| COUNTERPART_NAME| OTHER_NAME| CTPT_RA|CTPT_DEC| FLUX| FLUX_LO| FLUX_HI|CONTA|GAMM|GAMM|GAMM| CHI_SQ_R|REDSHIF| LUM|ASSOC_STREN| CL2| TYPE
| | deg| deg| | | | deg| deg|10^-12er|10^-12er|10^-12er| | | | | | |log[e| | |
1|SWIFT J0001.0-0708 | 0.261| -7.123| 6.10|2MASX J00004876-0709117 | | 0.2032| -7.1532| 13.03| 9.05| 17.56| 0.00|2.17|1.75|2.69| 0.60| 0.0000| NULL| 0| 2|Galaxy
2|SWIFT J0001.6-7701 | 0.326|-77.001| 5.41|Fairall 1203 | | 0.4419|-76.9540| 10.10| 6.61| 14.08| 0.00|2.02|1.54|2.60| 0.70| 0.0584|43.92| 0| 4|Sy1
3|SWIFT J0002.5+0323 | 0.664| 3.332| 5.10|NGC 7811 | | 0.6103| 3.3519| 11.69| 7.36| 16.64| 0.00|1.82|1.31|2.39| 0.50| 0.0255|43.24| 0| 4|Sy1.5 4|SWIFT J0003.3+2737 | 0.862| 27.676| 5.03|2MASX J00032742+2739173 ... ... どうやら、一行目には各列にどのようなデータが入っているかの説明が書かれているようだ。二行 目は、degなどと書かれているのをみると、どうやら各列の単位が書かれている。あと、6列目に は天体名が入っている。でもこいつら、列が多すぎていかんせん見づらい。たとえば、6列目だけ を抜き出したいときにはどうしたらいいんだろうか。その場合にはawkというプログラミング言 語が非常に便利だ。 $ which awk /usr/bin/awk などときちんとパスが返ってくれば、君のPCでもawkは使える。さっそく、awkを使って6列 目だけを抜き出そう。 *23moveの略 *24change directoryの略。
$ awk -F"|" ’{print $6}’ BAT_70m_catalog_20nov2012.txt ... 2MASX J23470479+5142179 PKS 2349-01 IGR J23524+5842 2MASX J23525143-1704370 H 2356-309 PKS 2356-61 $ このコマンドには少し解説が必要だろう。今、右端にファイル名が書かれている。このファイルの 中にあるデータの6列目を抜き出したい。それがawk ’{print $6}’の部分だ。ただ、これだけでは 今回はうまくいかない。awkは基本の区切り文字は空白“ ”を採用しているので、区切り文字は| であることを伝える必要がある。それがオプション部分の-F”|”だ。まとめると、上のコマンドは
「BAT 70m catalog 20nov2012.txt というファイルの6列目を、区切り文字|で(一行ずつ)表示 せよ。」となる。
次に、今回のカタログの座標を切り出してみよう。このカタログには座標は3列目(R.A.*25)、4
列目 (Dec*26) と8列目 (ctpt RA)、9列目 (ctpt DEC)に入っている。元論文*27を読む限り、8
列目、9列目の座標を使った方が良さそうだ。なので、今後はカタログから8列目と9列目を抜き
出そう。またもやawkを使って、
$ awk -F"|" ’{print $8,$9}’ BAT_70m_catalog_20nov2012.txt ... 358.1120 58.7000 358.2143 -17.0770 359.7830 -30.6280 359.7682 -60.9165 とすればよい。さきほどの$6の部分が、$8, $9となっているのがわかる。これで、「8列目と9列 目を表示する」という意味になるのだ。非常に直感的でわかりやすい。 さて、どうせなので、後でいくらでも確認できるように、この座標を新たなファイルに保存して しまおう。これは簡単で、>という記号を使う。ファイル名は“ra vs dec.txt”としよう。
$ awk -F"|" ’{print $8,$9}’ BAT_70m_catalog_20nov2012.txt > ra_vs_dec.txt $ ls
BAT_70m_catalog_20nov2012.dat ra_vs_dec.txt $ less ra_vs_dec.txt
CTPT_RA CTPT_DEC
*25R.A. は Right Ascension の 略 。日 本 語 で は 赤 経 と い う 。後 に 出 て く る 赤 緯 と 共 に 用 い る こ と で 、天 球 上 の 二 次 元 面 の 座 標 を 表 現 す る こ と が で き る 。天 文 学 で よ く 用 い ら れ る 座 標 系 の 一 つ で あ る 。詳 細 は
http://www.kusastro.kyoto-u.ac.jp/~iwamuro/LECTURE/OBS/coord.htmlなどを参考にしてみよう。
*26Declinationの略。赤緯という。
deg deg 0.2032 -7.1532 0.4419 -76.9540 0.6103 3.3519 0.8643 27.6548 1.0082 70.3217 1.5813 20.2029 ... と新たに“ra vs dec.txt”というファイルが作られているのがわかる。このように>を用いて (標 準) 出力先をファイルに書き換えることを「リダイレクト」と言う。この機能は非常によく使うの で、必ず覚えておこう。 さて、さきほどのファイル (ra vs dec.txt) を見てみると最初の三行ほどは座標以外のものが 入ってしまっている。これらを除きたいときは、たとえば以下のようにすればよい。
$ awk ’NR>3{print}’ ra_vs_dec.txt > ra_vs_dec_nohead.txt
NRというのは行数を表す。awk は一行ごとに処理をするため、一行目の処理のときにはNR は1、二行目の処理のときには NR は2,m行目の処理のときには NR は整数 m が入ってい る。’NR>3{print}’の部分は、「NRが3より大きいときにのみ、printを行う。」という意味にな る。awkの場合、printの後に何もかかなければ、すべての列を表示する、という意味になる。ま とめると、このコマンドは「ra vs dec.txtの4列目以降をすべて表示せよ」という意味であること がわかる*28。 さて、本当に余分な部分は切り取られただろうか。その確認方法は、lessコマンドで可能であっ た。実は、別の方法もある。要は最初の数行が切り取られているかを確認すればいいわけだから、 headコマンドも使うことができる。headコマンドは、ファイルの先頭から10行を表示するコマ ンドだ。 $ head ra_vs_dec_nohead.txt 0.2032 -7.1532 0.4419 -76.9540 0.6103 3.3519 0.8643 27.6548 1.0082 70.3217 1.5813 20.2029 2.2982 -0.6152 2.6292 10.9749 4.2853 81.5856 5.2814 -19.1682 見事、最初の3行分は切り取られているのがわかる。ファイルの上からn行を取り出したい時 は、オプション−nをつければよく、 *28実はこのコマンドは、{print}部分を省いても、同様の操作が行える。
$ head -n ra_vs_dec_nohead.txt などとすればよい。
3.4
プロットしてみる
さて、座標を手に入れたのだから、それをプロットしてみたくなるものだろう。ここでとうとう gnuplotの登場だ。gnuplotと打ってみよう。 $ gnuplot G N U P L O TVersion 4.6 patchlevel 1 last modified 2012-09-26 Build System: Darwin x86_64
Copyright (C) 1986-1993, 1998, 2004, 2007-2012 Thomas Williams, Colin Kelley and many others
gnuplot home: http://www.gnuplot.info faq, bugs, etc: type "help FAQ"
immediate help: type "help" (plot window: hit ’h’)
Terminal type set to ’x11’ gnuplot>
みたいなものが出てくれば、gnuplotにログイン完了だ。
さて、まずは簡単なものを試してみよう。
gnuplot> plot sin(x)
と打てば、図2が表示されるだろう。 -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 -10 -5 0 5 10 sin(x) 図2 gnuplotで作成したy = sin(x)の図。
このように、gnuplotは最初に用意されている関数をもとにプロットをすることができる。
gnuplotは我々が用意したファイルをプロットすることもできる。次は以下のように打ってみ
よう。
gnuplot> plot "ra_vs_dec_nohead.txt" u 1:2
すると、図3のような図が現れる。このコマンドは、 ”ra vs dec nohead.txt”の1列目を横軸、
2列目を縦軸にして、データをプロットせよ、という意味をもつ。実はu 1:2の“u”は“using”の
略で、ユニーク性が保たれていれば、gnuplotはコマンドの省略が可能である。gnuplot では様々
な条件を“set ∼”と書くが、“se ∼”でも同様の挙動をしてくれる。gnuplotは文字を打つ手間を
省いてくれる優しさも持ち合わせているのだ。 -100 -80 -60 -40 -20 0 20 40 60 80 100 0 50 100 150 200 250 300 350 400 "ra_vs_dec_nohead.txt"
図3 天球上のBATカタログ天体の分布。横軸はR.A.、縦軸はDec
あらためて図3をみてみよう。基本的にはまんべんなく天体は分布しているのがわかるが、座標 (180,−60)くらいから(300, 40)くらいまでにかけて、一筋の天体の塊が見えている。これはなん だろうか。答えは後でわかるはず。 図の見栄えをもう少しマシにしよう。まず、右上の“ra vs dec.txt”の部分は図が見えづらくな るので省きたい。横軸、縦軸の範囲も座標が分布する0≤ x ≤ 360と−90 ≤ y ≤ 90 にしてみた い。その場合は、続いて以下を打ってみよう。 gnuplot> se xrange[0:360] gnuplot> se yrange[-90:90]
gnuplot> plot "ra_vs_dec_nohead.txt" u 1:2 notitle
図4のようになるはずだ。この他にも、色を変更したり、プロット点の種類を変更したり、大き さを変更したりも可能である。また、軸の文字が見づらい場合も好きなように変更できる。興味が ある人は自分で調べてみよう。 さて、最後に図を保存してみよう。自分の見やすい形式 (.pngとか) にしてくれてもよい。こ の文書の執筆時では、天文学の業界では保存形式は.epsあるいは.pdfが主流である。なので今回 は.eps形式の保存方法を紹介しよう。他の形式も保存方法はほとんど変わらないので、ぜひ試して もらいたい。
-80 -60 -40 -20 0 20 40 60 80 0 50 100 150 200 250 300 350
図4 天球上のBATカタログ天体の分布改訂版。横軸はR.A.、縦軸はDec
Terminal type set to ’epscairo’
Options are ’ transparent fontscale 0.5 size 5.00in, 3.00in ’ gnuplot> se output "radec.eps"
gnuplot> replot gnuplot> exit
最初に保存したいeps形式を指定 (se term eps) して、その次に保存するファイル名を指定 (se
output ”radec.eps”)している。その後、もう一度プロットしてやる (replot)。最後に、gnuplot
をログアウト(exit)。これで、図4と同じものが”radec.eps”として保存されているはずだ。
$ ls
BAT_70m_catalog_20nov2012.dat ra_vs_dec.txt ra_vs_dec_nohead.txt radec.eps
ちなみに、eps ファイルはMacユーザーや一部の Linuxディストリビューションであれば、
openコマンドで見ることが可能である。 $ open radec.eps openコマンドに相当するコマンドはディストリビューションごとにいくつかあるようなので、 Linuxユーザーは各自チェックされたし。
3.5
どれくらい遠い天体がいるのだろう
BATカタログには天体の距離情報として、18列目に赤方偏移z*29の情報が入っている。また、 19列目にはX線光度*30をlogスケールで表した量が入っている。この2つの量を取得してプロッ トをすれば、このカタログには、どの距離にどれくらいの明るさのものが存在するかが見て取れる*29赤方偏移, 英語でredshiftは慣例的にz を用いる。定義は、本来の静止波長λiniに対して、観測波長がλobs
のとき、λobs = (1 + z)λini で定義される。つまり、波長の伸び具合をz で表している。z の値が大きいほ
ど我々から遠い場所にあることを示していて、z = 0.1でだいたい10 億光年ほど離れている。世の中には z を入れるとどれくらい離れているかを計算してくれるウェブサイトはたくさんあり、たとえばcosmo calc http://www.astro.ucla.edu/~wright/CosmoCalc.htmlなどはよく見かけるウェブサイトである。
だろう。もう皆さん慣れてきたと思うけれども、こういう時はawkが便利だ。
$ awk -F"|" ’{print $18,$19}’ BAT_70m_catalog_20nov2012.dat ... 0.1640 44.95 -1.0000 NULL 0.1651 45.16 0.0963 44.52 となる。どうやら、いくつかの天体にはzの情報に−1.0000が入っているようだ。その時には2 列目は”NULL”が入っている。このような天体は、信頼できるz情報が得られないほど近傍の天体 の場合が多い。今回は、そういう天体は何天体いるのかを一応調べておこう。そのような条件の天 体数をカウントするだけのシンプルな作業だ。もちろん手で数えてはいけない。人間は必ずミスを する。なので、Unixコマンドに頼ろう。そういう時は、wcコマンドが活躍してくれる。 $ wc -l BAT_70m_catalog_20nov2012.dat 1213 BAT_70m_catalog_20nov2012.dat のように、ファイルの行数を数えてくれる。これを使うと、
$ awk -F"|" ’$19~/NULL/{print $18,$19}’ BAT_70m_catalog_20nov2012.dat | wc -l 432 と書ける。少々難しいかもしれないが、これは、|の左半分で「19列目にNULLという単語が入っ ている行のみ、18、19列目を表示する」ということを行い、その表示をそのまま|の右側に引き 継ぎ*31、wc -lコマンドで表示された行数を数えている。そして、返された数が432である。つま り、1213天体*32中、432天体には信頼できるz情報が入っていない。 さて、当初の目的は「このカタログにはどの程度のzにどの程度の光度の天体が存在しているの か?」であった。上記の432天体は超近傍にいることはわかっているので、今回はこれを省いて ファイルに出力して、プロットをしてやろう。
$ awk -F"|" ’NR>3 && $19!~/NULL/{print $18,$19}’ BAT_70m_catalog_20nov2012.dat > z_vs_lx.txt
awkの条件にいろいろ書かれているが、これは「4列目以降 (NR> 3)かつ、19列目にNULLと いう単語が入っていない行 ($19!∼/NULL/) のみ、18、19列目を表示する。」というものである。 これで、一列目にz、二列目にX線光度log LXが入っている“z vs lx.txt”を作ることができた。 さて、どうせなので、最遠方天体のzを調べてみよう。こういうときは、sortコマンドが便利だ。 $ sort z_vs_lx.txt ... 3.1040 48.36 3.2000 47.98 3.2680 48.42 3.3660 48.01 *31このように出力を次のコマンドに引き継ぐ役割を持つ|を、「パイプ」という。 *32正確には、最初の3行には天体情報が入っていないので、1210天体
3.6680 48.14 と、1列目に入っているzを小さい方から順に並べてくれている。どうやら最大はz = 3.668のよ うだ。同じように、X線光度の最小値、最大値も調べてみよう。X線光度は2列目のデータなの で、sortコマンドに2列目でsortするよう-kオプションを加えてやれば良い。 $ sort -k2 z_vs_lx.txt ... 3.2680 48.42 2.4000 48.45 2.3450 48.55 2.9793 48.57 2.5000 48.63
X線光度の最大値はlog LX/erg s−1 = 48.63のようだ。太陽光度はおよそ4× 1033 erg s−1なの
で、この天体は太陽よりおよそ1015倍は明るいことになる。 次に、X線光度の最小値を探そう。sortとheadを合わせて $ sort -k2 z_vs_lx.txt | head -5 0.0007 39.81 0.0011 40.79 0.0029 41.02 0.0015 41.03 0.0033 41.47 としてもよいし、sortのオプション-r (逆順表示) を用いて、 $ sort -k2 -r z_vs_lx.txt ... 0.0033 41.47 0.0015 41.03 0.0029 41.02 0.0011 40.79 0.0007 39.81 としてもよい。どうやらX線光度の最小値はlog LX/erg s−1 = 39.81だ。それでも太陽光度のだ いたい106倍である。 さて、最後に“z vs lx.txt”に入っているデータをgnuplotでプロットしてみよう。 $ gnuplot gnuplot> plot "z_vs_lx.txt" u 1:2 図5を見てみると、ほとんどの天体はz < 0.5に分布していて、まばらにz = 4.0近くまで天体 が分布しているのがわかる。桁を超える幅広いパラメータ範囲を表示するにはlog表示が便利なの で、x軸をlog表示してみよう。gnuplotにはそのような機能も備わっている。
39 40 41 42 43 44 45 46 47 48 49 0 0.5 1 1.5 2 2.5 3 3.5 4 "z_vs_lx.txt" 図5 BATカタログ天体のz (x軸)とX線光度(y軸)の関係 gnuplot> se logscale x gnuplot> plot "z_vs_lx.txt" u 1:2 としてやると、図6のようなプロットが見られるはずだ。 39 40 41 42 43 44 45 46 47 48 49 0.0001 0.001 0.01 0.1 1 10 "z_vs_lx.txt" 図6 BATカタログ天体のz (x軸)とX線光度(y軸)の関係 図6を見てあげると、z = 0.01− 0.3にかけて多くの天体が分布しているのがわかる。また、皆 さんも気づいたかもしれないが、図の右下には天体がまったく存在していない。これは遠方宇宙 (z∼ 1)には、近傍宇宙 (たとえばz∼ 0.01) で見られるような光度 (log LX= 42− 44) の天体は まったくいないという意味なのだろうか。必ずしもそれは真とは限らない。BATサーベイも含め て、その装置ごとに限界等級や限界fluxが決まっているため、より遠方のものほど、もともと明る いものしか受かって来ないのだ。なので、右下の領域(z大かつlog LX小) に存在するような天体 はBATサーベイではそもそも受からず、議論をすることはできない。一方、左上の領域 (z小か つlog LX大) も天体はほとんどいないが、これは純粋にそのような天体はレアであることがわか る。ここに書いてある説明はやや難しいかもしれないが、興味がある人は、横軸をz,縦軸をX線
fluxとしてプロットをしてみるとよい。fluxと光度 (luminosity) の違いや、限界fluxの理解の助
3.6
筋模様の天体を暴く
さて、図3や図4などに見られた筋模様はどのような天体が分布して作っているのだろう。あり
がたいことに、BATカタログは天体を種族ごとにわけてidしており、その情報はカタログの21列
目に入っている。種族ごとにR.A., Dec情報を取得して、同様にplotしてみよう。たとえば、い
わゆる活動銀河核 (Active Galactic Nuclei; AGN) と呼ばれる種族は、カタログではid =4–8に
分類されている。これを取得してみよう。
$ awk -F"|" ’$21>=4 && $21<=8 {print $8,$9}’ BAT_70m_catalog_20nov2012.dat > ra_vs_dec_AGN.txt
として、gnuplotでplotしてみる(自分でやってみよう)。すると、図7のようになるはずだ。 -100 -80 -60 -40 -20 0 20 40 60 80 100 0 50 100 150 200 250 300 350 400 "ra_vs_dec_AGN.txt"
図7 BATカタログ内のAGNの分布。横軸はR.A.、縦軸はDec。
図7を見てやる限り、AGNは一様に分布しているように見えるが、どうやら筋模様に該当する
部分にはあまり存在していない。AGNは筋模様の起源ではなさそうだ。AGNは銀河系の外に存
在する天体なので、筋模様は銀河系内の天体が作り出しているのかもしれない。同様にして、いく
つかの種族の座標を取り出してファイルに保存し、gnuplotで表示してみるとよい。よいawkの訓
練になるはずだ。
銀河系内の天体としてはX線連星 (X-ray binary) やパルサー、激変星(cataclysmic variable;
CV) などが考えられる。これら系内天体候補をすべて含むためには、id = 1と9≤ id ≤ 15の種 族の座標を取り出せば良い。取り出す方法は今の読者ならうまくできるはずだ。わからない場合は 自分で調べてみよう。調べる力も大事だ*33。座標を表示してみると、図8のようになる。 図8を見てやれば、見事、系内天体が筋模様を作っていたことがわかる。実は筋模様に独立し て、図の左下あたりに2つほど点の分布 ((80,−70)と(20,−75)あたり) があることがわかるが、 この分布はどこ由来なのかを調べてみるのも勉強になるだろう。
3.7
まとめ
以上、一つのカタログを元に、天体の様々な分布を見てやることができた。Unixコマンドやawkの大事さが伝わっただろうか。ここからはさらにUnixやawkの知識をどんどんつけていく
-100 -80 -60 -40 -20 0 20 40 60 80 0 50 100 150 200 250 300 350 400
図8 BATカタログ内の系内天体の分布。横軸はR.A.、縦軸はDec。
のもよいし、少し背伸びをしたければ、Pythonというプログラミング言語を使えば、皆さんがや
りたいことの幅をさらに広げてくれるだろう。
実際の研究では今回よりもさらに条件を絞ったり、より定量的な議論を行ったり、他のカタログ とマッチングを取ったりしていろいろなことを行うわけだが、そのあたりは実際に模擬研究で体験 できるはずだ。
4
宿題
1. 今回紹介したUnixのコマンド、およびUnixの大事な概念であるリダイレクト、パイプを
自分で復習せよ。
2. BATカタログデータ“BAT 70m catalog 20nov2012.dat” には最初三行だけ天体情報が
入っていない。Unixコマンドwcを用いてファイルの行数を数え、tailコマンドを用い
て*34、最初の三行を削った新たなファイル“BAT 70m catalog nohead.dat”を作成せよ。
作成後に、wcコマンドを用いて、ファイルの行数が三行少なくなっていることを確認せよ。
3. 世の中にはコメントアウトという概念があり、コメントアウトしたところは飛ばしてデータ
を読み込むオプションが様々なプログラミング言語で用意されている。カタログデータの最
初3行のみ、行の先頭に#を追加し、ファイル名“BAT 70m catalog cmthead.dat”として
保存せよ。
4. BATカタログデータからz > 1を満たす天体のみを抜き出し、“BAT 70m catalog zgt1.dat”
として保存せよ。また、天体数も数えて答えよ。
5. “BAT 70m catalog zgt1.dat”にはどのような種族が入っているだろうか。2種類に絞られ
るはずだが、その種族の名前(BATカタログウェブサイトでいうところのsource type)を
答えよ。また、それぞれの種族の天体数を数えよ。
6. “BAT 70m catalog zgt1.dat”に入っている2種族を分類し*35、それぞれ“BAT 70m catalog zgt1 ty7.dat” 、“BAT 70m catalog zgt1 ty8.dat” として保存せよ。
7. “BAT 70m catalog zgt1 ty7.dat”および“BAT 70m catalog zgt1 ty8.dat”に入っている
天体の種族の座標分布(横軸をR. A.,縦軸をDec)をそれぞれの種族の分布がわかるように 一つの図に表示し、保存せよ。余裕があれば、各自gnuplotのオプションを調べ、横軸や縦 軸のタイトルを加えたりフォントサイズを大きくしたりなどして、図を見やすくする工夫を してみよ。 *34tail: ファイルの末尾を表示するコマンド。 *35ヒント: たとえば、awkで条件を付けて表示してやればよい。
•
参考文献・ウェブサイト
世の中にはUnixやawkの解説書はたくさんあるが、まずはgoogleで検索して、解答にたどり着
く力を得られるようにしよう。あなたがMacユーザーで、Unixに初めて触れて、ターミナルでさ
らにどんなことができるのかを知りたいのであれば、入門Unix for Mac OS X*36という本が
ある。
また、awkをもっと詳しく学びたいという人もいるかもしれない。その場合にはプログラミング
言語AWK*37 という素晴らしい本がある。この本の後半はこのご時世、もはやawkでやる必要
があるのか、という気もするが、前半は便利なワンライナー*38がたくさん載っている。
また、本書で紹介した BAT カタログ内の天体を細かく調べてみたい人がいるかもしれな
い。世の中には天体の検索方法はたくさんある。まずは、NASA Extragalactic Database
(NED)*39で調べてみるとよいだろう。天体の座標やzだけでなく、天体がどのような種族 (銀河
なのかBlazarなのか等) なのか、過去にどの波長でその天体が観測されていて、その値はいくら
で· · · という情報まで載っている。また、SIMBAD Astronomical Database*40も同様に有
用だ。Simbadでは天体検索をすると、検索結果の右側の画面*41で、各波長の画像が非常に簡単に
見ることができる。しかし肝に命じてもらいたいのは、「NEDやSimbadに情報がない」からと
言って、それを根拠に結論を出すのは危険である。NEDやSimbadはあくまで過去の論文やサー
ベイの結果を集約して表示しているのであって、過去の観測のすべてを網羅しているのではないこ
とをきちんと理解しておこう。(ここから書くことは読み飛ばしてもらってかまわないが) 論文
でAGNの可否をNEDなどに頼っている場合が見られるが、「AGNという情報がない」からと
言って「この天体はAGNではない」という結論に必ずしもなるわけではない。気をつけよう。
•
著者紹介
市川 幸平
(Ichikawa, Kohei)
https://motigomeman.wordpress.com/国立天文台ハワイ観測所プロジェクト研究員。京都大学大学院 宇宙物理学教室で博士号(理学) を 取得。超巨大ブラックホールの成長現場として活動銀河核に興味を持ち、活動銀河核(AGN)がど のようにして誕生してはブラックホールに物を落とし、最期を迎えるのかを研究している。 *36http://www.oreilly.co.jp/books/9784873116440/ *37https://www.usp-lab.com/book.awk.html *38意味がわからなければ調べるのである。 *39http://ned.ipac.caltech.edu/ *40http://simbad.u-strasbg.fr/simbad/ *41AladinLiteという
5
解答
1. 略
2. $ wc BAT_70m_catalog_20nov2012.dat
1213 20988 301789 BAT_70m_catalog_20nov2012.dat
$ tail -1210 BAT_70m_catalog_20nov2012.dat > BAT_70m_catalog_nohead.txt $ wc BAT_70m_catalog_nohead.txt
1210 20956 301290 BAT_70m_catalog_nohead.txt
3. $ awk ’NR<=3{print "#"$0}NR>3’ BAT_70m_catalog_20nov2012.dat > BAT_70m_catalog_cmthead.txt
4. $ awk -F"|" ’$18>1.0’ BAT_70m_catalog_nohead.txt > BAT_70m_catalog_zgt1.txt $ wc -l BAT_70m_catalog_zgt1.txt
30 BAT_70m_catalog_zgt1.txt
より、は30天体存在する。
5. $ awk -F"|" ’{print $21}’ BAT_70m_catalog_zgt1.txt | sort | uniq 7
8
より、source typeは7: Blazar/BL Lacと8: QSOである。また、それぞれの天体数は
$ awk -F"|" ’$21==7{print $21}’ BAT_70m_catalog_zgt1.txt | wc -l 11
$ awk -F"|" ’$21==8{print $21}’ BAT_70m_catalog_zgt1.txt | wc -l 19
より、種族7は11天体、種族8は19天体である。
6. $ awk -F"|" ’$21==7 {print}’ BAT_70m_catalog_zgt1.txt > BAT_70m_catalog_zgt1_ty7.txt $ awk -F"|" ’$21==8 {print}’ BAT_70m_catalog_zgt1.txt > BAT_70m_catalog_zgt1_ty8.txt 7. $ gnuplot
gnuplot> se datafile sepa "|"
gnuplot> plot "BAT_70m_catalog_zgt1_ty7.txt" u 8:9, "BAT_70m_catalog_zgt1_ty8.txt" u 8:9 gnuplot> se term eps
Terminal type set to ’epscairo’
Options are ’ transparent fontscale 0.5 size 5.00in, 3.00in ’ gnuplot> se output "radec_zgt1.eps"
gnuplot> replot gnuplot> exit
-60 -40 -20 0 20 40 60 80 100 0 50 100 150 200 250 300 350 400 "BAT_70m_catalog_zgt1_ty7.txt" u 8:9 "BAT_70m_catalog_zgt1_ty8.txt" u 8:9