HSC 解析講習会用チュートリアル

2019年度 HSC 画像解析講習会用のチュートリアルです。実習では Pipeline を用いた解析と、最終的なデータで疑似カラー画像作成と色等級図作成を行います。解析に使うコマンドの中には今回の講習会に特化したコマンドオプションもありますので、注意してください。

HSC pipeline マニュアル (ver. 6用です) と HSC pipeline バイナリ配布所 にも詳細がありますので、参照下さい。 また、上級者には こちら も参考になると思います。Pipeline についてかなり細かく書かれていますが、古いバージョンの情報もありますので、参照する際にはご注意下さい。

今回使用するデータは 3 バンド (i, z, y)で、中心 6 CCD を 5 ショット分です(Hubble Deep Field North, HDFN領域で取得されたデータです)。バイアス、ダーク、ドームフラットについても同様です。

また、HSC pipeline のバージョンは最新の 7.9.1 になります。

多波長解析で実習を行う皆様： 今回の講習会では8コアを使うプロセスがあります。稼働状況をご確認の上、解析を進めてください。

まずはじめに、"unset LD_LIBRARY_PATH"を行ってください。新しいターミナルを開いた場合、そのターミナルについても同様に"unset LD_LIBRARY_PATH"を行ってください。 
 unset LD_LIBRARY_PATH

#"echo $LD_LIBRARY_PATH"をして、何も表示されないことを確認しておくと良いと思います。
 echo $LD_LIBRARY_PATH


# HSC pipeline セットアップコマンドの登録
source /work/hscpipe/7.9.1/bashrc

 # HSC pipeline を呼び出すためのセットアップコマンド (ログインの度に必ず行ってください)
setup-hscpipe

echo $LD_LIBRARY_PATH

Rawdata の準備

Pipeline とアストロメトリ用カタログファイルはインストール済みですので、まず、rawdata をダウンロードしましょう。

 # 個人の解析ディレクトリを作成します。
 cd /work
 mkdir [User Name]
 # 自分の好きなように作って下さい。
 
 # Rawdata をダウンロードし展開
 cd /work/[User Name]/
 curl -L http://hsc.mtk.nao.ac.jp/HSC_Training_download/rawdata_tutorial2018.tar.gz | tar -xzf -
 
 # Rawdata の確認
 cd /work/[User Name]/rawdata/
 ls | wc -l
 # 結果が 240 であれば問題ありません。

 # どのような画像が撮れているのか確認してみましょう。
 ds9 xxxxx.fits
 # どれでもいいので、自分の好きなファイルを選んでください。

ここで、生データと pipeline 内部でのデータ名について復習しておきましょう (データ名について)。解析コマンドの中でデータを指定する際にも必要になります。

続いて、解析用ディレクトリを準備し、レジストリを作成します。

# 解析用ディレクトリ（リポジトリの親ディレクトリ）の作成
 mkdir -p /work/[User Name]/HSC

 # _mapper ファイルの作成 (これがないと処理が進みません)
 echo 'lsst.obs.hsc.HscMapper' > /work/[User Name]/HSC/_mapper

/work/[User Name]/HSC の下を確認すると、_mapper という名前のファイルが作られています。これで解析用ディレクトリの準備は完了です。

次に、rawdata を解析用ディレクトリに配置し、登録します。

# 生データのあるディレクトリに移動
cd /work/[User Name]/rawdata

# 解析ディレクトリ下にレジストリを作成
hscIngestImages.py /work/[User Name]/HSC  HSC*.fits --mode=link --create
# --mode は link, move などがあります。link にすると生データへのリンクを作成します。

#レジストリ登録時にできたディレクトリを確認
ls /work/[User Name]/HSC 

/work/[User Name]/HSC の下には、BIAS、DARK、DOMEFLAT、HDFN という名前のフォルダと、registry.sqlite3 というファイルが作られています。

では、レジストリに登録された画像の情報を確認しましょう。 registry.sqlite3 に登録された情報は SQL というデータベース言語を用いて検索することができます (SQL の詳細についてはここでは触れません)。

cd /work/[User Name]/HSC
sqlite3 /work/[User Name]/HSC/registry.sqlite3

sqlite> .header on
sqlite> SELECT visit, filter, field, count(visit)
        FROM raw
        GROUP BY visit, filter, field;

# 以下のように、カラムの右端 (= CCD 数) が全て 6 になっていればOK.
visit|filter|field|count(visit)
212|HSC-Y|DOMEFLAT|6
214|HSC-Y|DOMEFLAT|6
216|HSC-Y|DOMEFLAT|6
218|HSC-Y|DOMEFLAT|6
220|HSC-Y|DOMEFLAT|6
710|HSC-Y|HDFN|6
712|HSC-Y|HDFN|6
714|HSC-Y|HDFN|6
716|HSC-Y|HDFN|6
718|HSC-Y|HDFN|6
950|HSC-I|HDFN|6
952|HSC-I|HDFN|6
954|HSC-I|HDFN|6
...
# ここで見えているカラム名で解析に使うデータを指定します。
# それぞれの filter, field の visit 何番から何番までかをメモしておきましょう。

# 全てのカラムが見たい場合は
sqlite> select * from raw;

# 終了
sqlite> .q

以上で rawdata の準備は終了です。

Brighter-Fatter カーネルの設定

hscPipe 7 では、１次処理の際に行われる brighter-fatter 効果 (明るい星のような高カウントのピクセルは、有効面積が小さくなり、隣のピクセルに電荷がいきやすい。その結果、天体が広がって見える効果) を取り除くカーネルを自身の解析ディレクトリの下にコピーする必要があります。 以下に方法を示します。

# BFKERNEL ディレクトリ (名前は必ず BFKERNEL にして下さい) を /work/[User Name]/HSC/CALIB の下に作成。
cd /work/[User Name]/
mkdir -p /work/[User Name]/HSC/CALIB/BFKERNEL
        
# brighter_fatter_kernel.pkl は $OBS_SUBARU_DIR/hsc/ の下にあるので、そこへ向けてリンクを張ります
ln -s /work/hscpipe/7.9.1/lsst_home/stack/miniconda3-4.5.12-1172c30/Linux64/obs_subaru/7.7.2-hsc+4/hsc/brighter_fatter_kernel.pkl   /work/[User Name]/HSC/CALIB/BFKERNEL/

# リンクができているか、確認
ls /work/[User Name]/HSC/CALIB/BFKERNEL

    

アストロメトリ用カタログファイルへのリンク作成

hscPipe の CCD 解析では、１次処理の後にアストロメトリ較正、等級原点較正、天体の検出と測定を行います。 その際に外部のカタログファイルを参照します。

Note

hscPipe4 以前ではログインの度に setup astrometry_net_data コマンドを実行する必要がありましたが、hscPipe7 ではこのリンクを作成するだけで良いです。

        # 解析ディレクトリの下に ref_cats ディレクトリを作成 (名前は必ず ref_cats にして下さい)。
        mkdir -p /work/[User Name]/HSC/ref_cats
        
        # ref_cats に、リファレンスカタログへのリンクを張る。
        ln -s     /work/ps1_pv3_3pi_20170110     /work/[User Name]/HSC/ref_cats/
    
        # リンクを確認。
        ls  /work/[User Name]/HSC/ref_cats/

HSC filterのtransmission curveデータの設定

hscPipe 7になってからは、解析前にHSC filterのtransmission curveデータの設定を行う必要があります。
以下を実行することで、設定できます。

        #以下を実行すると、transmissionというディレクトリができ、その中に必要なファイルができます。
        installTransmissionCurves.py    /work/[User Name]/HSC
        
        # transmissionディレクトリを確認。
     　　 ls /work/[User Name]/HSC

Yバンド迷光

hscPipe 7になってからは、Yバンドの画像を処理するために迷光パターンのファイルが必要です (https://community.lsst.org/t/y-band-stray-light-correction-for-hsc/2517)。 以下からファイルを取得し、/work/[User Name]/HSC/CALIB/STRAY_LIGHTにおく必要があります。

        #迷光パターンのファイルをダウンロードして解凍。
        curl -L https://hscdata.mtk.nao.ac.jp/hsc_bin_dist/STRAY_LIGHT.tar.xz | tar -xJf - 

        #リンクを張る。
        ln -s   /work/[User Name]/STRAY_LIGHT     /work/[User Name]/HSC/CALIB/

        #リンクの確認。
        ls /work/[User Name]/HSC/CALIB/     
           

較正用データの作成

バイアスデータの作成

バイアスデータは 0 秒積分のデータ（バイアス生データ）から CCD の電荷を空読みした overscan を引くことで生成されます。 overscan 領域は各 CCD 各 CH の左端もしくは右端に 16 ピクセル分ついています。 hscPipe ではまず各 [visit, CCD] 毎に、overscan 領域の短軸方向に平均し、長軸方向に Spline 補間をかけたものをピクセルのカウント値から差引ます。 この各 [visit, CCD] のデータから外れ値を除き、全 visit で平均化したものが最終的なバイアスデータとなります。 バイアスデータ生成は constructBias.py で実行されます。

# バイアスデータ作成
            # 解析ディレクトリは /work/[User Name]/HSC, 1次処理用データディレクトリを /work/[User Name]/HSC/CALIB,
                rerun 名を calib_tutorial とする。
            constructBias.py /work/[User Name]/HSC \
            --calib /work/[User Name]/HSC/CALIB \
            --rerun calib_tutorial --id visit=1360..1368:2 \
            --batch-type=smp --cores=8
            
            # constructBias.py [解析ディレクトリ] --calib [1次処理用データディレクトリ] --rerun [rerun 名] --id [visit]
            # オプション
            #       --calib: 正規表現 [A-Za-z0-9_+-] にマッチするものでなければいけません。カンマやピリオドは含まれないようにして下さい。
            #       --rerun: rerun 名。
            #       --id: データID、ここでは Bias 生データの visit 番号を指定。
                visit=1360..1368:2 は 1360 から２つおきに 1368 まで指定。
            #              一つずつ visit を指定したい場合は visit=1360^1362 のように ^ で連結する。
            #              他にも registry.sqlite3 の中のデーブルのコラム名を使って指定することも可能 (例えば field=BIAS)。
            #       --batch-type [slurm, pbs, smp]: ジョブを投げるバッチ処理の方法を指定
            #       --cores: コア数の指定 (slurm, smp のみ)
            #       完了したらds9で画像を確認してみましょう。 ds9 /work/[User Name]/HSC/rerun/calib_tutorial/BIAS/*/*/BIAS-[ccd].fits
        

そして、出来上がったバイアスデータをレジストリに登録します。これを忘れると、以降の解析がすべて失敗します。

# レジストリに登録
            ingestCalibs.py \
            /work/[User Name]/HSC \
            --calib /work/[User Name]/HSC/CALIB \
            /work/[User Name]/HSC/rerun/calib_tutorial/BIAS/*/*/BIAS*fits --mode=link \
            --validity=30 \
            --create
            
            # --validity: キャリブレーションデータが適用される日数。
            # キャリブレーションデータを取得した日を挟んだ XX 日間に取得されたデータにはこのキャリブレーションデータは有効、という意味。
            # 観測が複数ランに渡っている場合、それぞれの期間で正しく有効期間を設定しないとキャリブレーションデータが使えません。
            # ちゃんと登録できたか、確認してみましょう。sqlite3 /work/[User Name]/HSC/CALIB/calibRegistry.sqlite3 
        

これでバイアスデータの作成は終了です。

ダークデータの作成

ダークデータはある積分時間の間に CCD にたまる暗電流を較正するためのデータで、ダーク生データから Bias データと overscan を引くことで生成されます。 hscPipe ではまず各 [visit, CCD] の ダーク毎に overscan の平均をピクセルのカウント値から差し引いたデータを作ります。 この各 [visit, CCD] のデータの平均から バイアスデータを差し引いたものが最終的なダークデータとなります。 ダークデータ生成は constructDark.py で実行されます。constructDark.py はバッチ処理が可能です。

# ダークデータ作成
            constructDark.py /work/[User Name]/HSC \
            --calib /work/[User Name]/HSC/CALIB \
            --rerun calib_tutorial  --id visit=1466..1474:2 \
            --batch-type=smp --cores=8
            
            # constructDark.py [解析ディレクトリ] --calib [1次処理用データディレクトリ] --rerun [rerun 名] --id [visit, field, ...]
            # オプションは constructBias.py と同じです。
            # 完了したらds9で画像を確認してみましょう。 ds9 /work/[User Name]/HSC/rerun/calib_tutorial/BIAS/*/*/BIAS-[ccd].fits
        

処理が終了したらデータを確認します。 生成された ダークデータは /work/[User Name]/HSC/rerun/calib_tutorial/DARK/[dateObs]/NONE/DARK-[dateObs]-[ccd].fits です。

ダークデータの確認後、登録が必要な場合はバイアスデータと同様にできます。

# ダークデータをレジストリに登録
            # オプションはバイアスデータの時と同じ
            ingestCalibs.py /work/[User Name]/HSC \
            --calib /work/[User Name]/HSC/CALIB \
            /work/[User Name]/HSC/rerun/calib_tutorial/DARK/*/*/DARK-*.fits --mode=link --validity 30
            # ちゃんと登録できたか、確認してみましょう。sqlite3 /work/[User Name]/HSC/CALIB/calibRegistry.sqlite3 
        

フラットデータの作成

ドームフラットからフラットデータを作成します。これはそれぞれのバンドで実行します。

constructFlat.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun calib_tutorial \
        --id field=DOMEFLAT filter=[filter] \ # [filter]にはfilter名をいれてください（HSC-Yの場合は、filter=HSC-Y）。filter ごとに計 3 回走らせます。
        --batch-type=smp --cores=8

# rerun とは HSC pipeline で使われる解析プロセスの概念です。後ほど詳しく説明します。
# --id は入力データを検索する条件を入れます。今回は field と filter を指定していますが、visit 番号を指定することもできます (visit 980 から 994 まで二つ飛び；visit=980..994:2)。
# 基本的には visit 指定を推奨していますが (pipeline マニュアルでは全て visit 指定です)、講習会では直観的に理解しやすい field 指定にしています。
# field などで指定する場合は、 FITS ヘッダが間違っていたときに、どのデータが入力されたのかわかりにくくなるのでご注意下さい。
# --batch-type=smp --cores=8 は「ローカルコンピュータ上で処理を実行、8 スレッドで」という意味です。
# このあたりのパラメータについては後ほど解説します。
# 完了したらds9で画像を確認してみましょう。 ds9 /work/[User Name]/HSC/rerun/calib_tutorial/FLAT/*/*/FLAT-[ccd].fits

/work/[User Name]/HSC/rerun/calib_tutorial/FLAT/[dateObs]/[filter]/[calibVersion] の下にフラットデータ (FLAT-[ccd].fits) ができます。ds9 で好きな画像ファイルを選んで、できた画像を確認してみましょう。 filterごとに計３回、走らせましょう。

次にフラットデータをレジストリに登録します。

# レジストリに登録
  ingestCalibs.py /work/[User Name]/HSC \
  --calib /work/[User Name]/HSC/CALIB/ \
  /work/[User Name]/HSC/rerun/calib_tutorial/FLAT/*/*/FLAT-*.fits --mode=link --validity 30

# --validity: キャリブレーションデータが適用される日数。
# キャリブレーションデータを取得した日を挟んだ XX 日間に取得されたデータにはこのキャリブレーションデータは有効、という意味。
# 観測が複数ランに渡っている場合、それぞれの期間で正しく有効期間を設定しないとキャリブレーションデータが使えません。
# ちゃんと登録できたか、確認してみましょう。sqlite3 /work/[User Name]/HSC/CALIB/calibRegistry.sqlite3 

これでフラットデータの作成は終了です。

フリンジデータの作成

フリンジデータは天体画像のスタックから生成されます。今回のデータのうち、y バンドについてはフリンジデータから画像の干渉縞を除去する必要があります。 Pipeline ではフリンジデータ作成はデフォルトでは y バンドと NB0921 のみ対応しています。

# フリンジデータ(y-band のみ)
constructFringe.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun calib_tutorial \
        --id field=HDFN filter=HSC-Y \
        --batch-type=smp --cores=8

# 天体画像を使うため、field=HDFN を指定。
# field 指定より visit 指定のほうが好きな場合は visit=710..718:2 でも OK。
# 完了したらds9で画像を確認してみましょう。 ds9 /work/[User Name]/HSC/rerun/calib_tutorial/FRINGE/*/*/FRINGE-[ccd].fits

完了すると、 /work/[User Name]/HSC/rerun/calib_tutorial/FRINGE/[dateObs]/[filter]/[calibVersion] の下にフリンジデータ (FRINGE-[ccd].fits) ができます。フリンジデータができていることを確認したら、フラットデータと同様にレジストリに登録しましょう。コマンドはフラットデータの時と同じです。どんな画像ができたのか気になる人は ds9 で確認してみましょう。

# レジストリに登録
 ingestCalibs.py /work/[User Name]/HSC \
 --calib /work/[User Name]/HSC/CALIB/ \
 /work/[User Name]/HSC/rerun/calib_tutorial/FRINGE/*/*/FRINGE-*.fits --mode=link --validity 30

以上でフリンジデータの作成は終了です。

スカイデータの作成

hscPipe 7以前では、特に大きな天体周りのskyを引き過ぎてしまう問題がありました。hscPipe 7からは、Sky データを作成し、それを利用することで、sky引きの方法が改善され、大きな天体周りのsky引きもうまくできるようになりました。 Sky データは constructSky.py を実行することで作成されます。実行方法は以下です。HSC-I, HSC-Zについても同様に行ってください。

# Sky データ作成
        constructSky.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun calib_tutorial --id field=HDFN filter=HSC-Y \
        --batch-type=smp --cores=8     
    

Sky データのファイル名は /work/[User Name]/HSC/CALIB/HSC/rerun/calib_tutorial/[dataObs]/[filter]/SKY-[dateObs]-[filter]-[ccd].fits です。

Sky データも忘れずに登録しましょう。

# Sky データをレジストリに登録
        ingestCalibs.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        /work/[User Name]/HSC/rerun/calib_tutorial/SKY/*/*/SKY-*.fits --mode=link --validity 30
    

キャリブレーションデータレジストリは /work/[User Name]/HSC/CALIB/calibRegistry.sqlite3 です。これも生データ同様、SQL で確認できます。

sqlite3 /work/[User Name]/HSC/CALIB/calibRegistry.sqlite3

# ヘッダ情報読み込み
sqlite> .header on
# テーブルを表示
sqlite> .tables
bias          dark          flat          fringe           sky
# bias, dark, flat, fringe, sky テーブルに登録したデータが入っている
sqlite> select * from flat;
id|filter|ccd|calibDate|validStart|validEnd|
1|HSC-Y|41|2014-03-252014-02-23|2014-04-24
2|HSC-Y|42|2014-03-252014-02-23|2014-04-24
3|HSC-Y|49|2014-03-252014-02-23|2014-04-24
4|HSC-Y|50|2014-03-252014-02-23|2014-04-24
5|HSC-Y|57|2014-03-252014-02-23|2014-04-24
6|HSC-Y|58|2014-03-252014-02-23|2014-04-24
7|HSC-Y|41|2014-03-252014-02-26|2014-04-27
8|HSC-Y|42|2014-03-282014-02-26|2014-04-27
9|HSC-Y|49|2014-03-282014-02-26|2014-04-27
...

# ccd が全て揃っているか、validStart/validEnd は --validity で指定した日数になっているかは見ておきましょう。

CCD ごとの解析

CCD 毎の解析を走らせます。まず、天体データの一次処理 (バイアス引き、ダーク引き、フリンジ除去、フラット割り、その他装置固有の処理) が行われます。そして、sky を引き、天体検出を行います。検出された天体はカタログファイルとマッチングされ、等級原点と座標決めが行われます。最終的に一次処理済み画像と sky 画像などの画像データと、カタログデータが生成されます。

singleFrameDriver.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id field=HDFN \
        --batch-type=smp --cores=8

# --rerun=tutorial を指定しています。ここからは本解析に入るので、先程の calibration とは rerun 名を変えます。
# --id は visit 指定を推奨します。

singleFrameDriver.py は約 30 分程度かかりますので、その間にコマンドのパラメータについて説明します。

こちらのページ -> パラメータの説明

singleFrameDriver.py が終了したら、まず画像データを確認していきましょう。

画像データの出力ディレクトリは /work/[User Name]/HSC/rerun/[rerun]/[pointing]/[filter]/corr、サムネイル画像の出力ディレクトリは /work/[User name]/HSC/rerun/[rerun]/[pointing]/[filter]/thumbs です。ここで、[rerun] は –rerun=tutorial で指定した rerun 名、[pointing] は観測時のポインティングで、y バンドでは 00815、i, z バンドでは 00816 になります。

つまり、y バンドでは画像データは /work/[User Name]/HSC/rerun/tutorial/00815/HSC-Y/corr の下に生成されることになります。

• 1 次処理済天体画像（sky も引き済）　　　：CORR-[visit]-[ccd].fits
• sky 引きで用いられた sky 画像　　　　　　：BKGD-[visit]-[ccd].fits
• overscan 引きの結果サムネイル　　　　　：oss-[visit]-[ccd].png
• Flat データで補正した天体画像サムネイル ：flattened-[visit]-[ccd].png

座学でも説明があった通り、CORR 画像には Science イメージ以外にも Mask イメージと Variance イメージがついています。 それぞれ、CORR-[visit]-[ccd].fits[1], CORR-[visit]-[ccd].fits[2], CORR-[visit]-[ccd].fits[3] です。 ds9 で Mask イメージを開いてみましょう。

# ds9 でマスクイメージを開く
ds9 /work/[User Name]/HSC/rerun/tutorial/00816/HSC-Z/corr/CORR-0001042-041.fits[2]

Mask イメージの値は header を見て確認しましょう。Maskイメージの詳細を知りたい方はパイプラインマニュアルをご覧ください。

こちらで詳細が確認できます。

次に、等級原点と座標決めのために用いられたカタログデータを紹介します。生成された全てのカタログデータは /work/[User Name]/HSC/rerun/[rerun]/[pointing]/[filter]/output の下にあります。

• 等級原点、座標決めのために天体画像から検出された比較的明るい天体のカタログ ： ICSRC-[visit]-[ccd].fits
• 等級原点、座標決めのために天体画像から検出された天体のカタログ　： SRC-[visit]-[ccd].fits
• 外部カタログ天体とクロスマッチさせたリスト　： SRCMATCH-[visit]-[ccd].fits
• 検出天体の測定結果と、それとマッチした外部カタログ天体の情報がまとまったカタログ　： SRCMATCHFULL-[visit]-[ccd].fits

これらカタログデータは .fits の形式ですが画像データではありませんので、ds9 では中身を見ることはできません。 簡単に中身を確認するには

fits viewer やpythonが便利です。fits viewerで中身を確認した場合 —> こちら

以上で CCD ごとの解析は終了です。

skyCorrection

hscPipe 7では、天体データの重ね合わせ（coaddDriver.py）の前に、skyCorrection.pyを行う必要があります。これを行うことで、大きな天体周りのsky引きもうまく行うことが可能となります。実行方法は以下です。

# 解析ディレクトリは ~/HSC, rerun 名は tutorial とします。
    skyCorrection.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial --id  field=HDFN --batch-type=smp --cores=8
    

/work/[User Name]/HSC/rerun/tutorial/00815/HSC-Y/corr の下にskyCorr-[visit]-[ccd].fitsが生成されることになります。

位置合せとフラックススケール決め (モザイク)

tract 定義

全ての天体データを使い、 tract を定義します。天球面をぶつ切りにした領域を tract と呼び、一つの tract は一つの平面に投影されます。 1 視野くらいのデータの場合は、その画像を全て含むような tract を定義しなおします。 ここで使用する makeDiscreteSkyMap.py ではヘッダーの座標情報をもとにデータの全観測領域に収まる最大の正方形を設定し、この領域を 1 tract と定義します。 tract は観測された領域全体から決められるため、必ず全ての天体データを指定し makeDiscreteSkyMap.py を実行してください。

makeDiscreteSkyMap.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id field=HDFN

tract 0 has corners ... (RA, Dec deg) and 4 x 4 patches のように出力されると成功です。ここで、patch とは tract をさらに細かく分けたものです。これで tract 0 番が定義されました。 この tract/patch の座標情報は /work/[User Name]/HSC/rerun/[rerun]/deepCoadd/skyMap.pickle に書かれます。 skyMap.pickle は同一 rerun 内では makeDiscreteSkyMap.py のたびに上書きされるのでご注意ください。

モザイク

モザイクでは、各 [visit, ccd] の等級原点と座標決めの情報に複数の visit の情報を加え、より精密に各 [visit, ccd] の位置合わせとフラックススケールを決定します。

export LD_LIBRARY_PATH="$LD_LIBRARY_PATH:/work/hscpipe/7.9.1/lib64"
# このパスを外すと共同利用室のマシン、多波長解析システムではmosaic.pyが動きません。
mosaic.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id field=HDFN filter=HSC-Y tract=0  --diagnostics --diagDir /work/[User Name]/HSC/rerun/mosaic_diag # filter を書き換えて計 3 回走らせます。

# mosaic.py はシングルプロセスで動きますので、--batch-type, --cores の指定をしてはいけません。
# --diagnostics : 評価結果を保存するときに使用 --diagDir : 評価結果を保存するディレクトリを指定  

モザイクではフラックススケールファイルと WCS 座標ファイルが生成され、画像については特に生成されません。i、zバンドについても filter=HSC-I, HSC-Z として実行します。

では、生成されたファイルを確認していきましょう。/work/[User Name]/HSC/rerun/[rerun]/jointcal-results/[tract] の下に fcr-[visit]-[ccd].fits と wcs-[visit]-[ccd].fits というファイルができているはずです。これらがそれぞれフラックススケールファイルと WCS 座標ファイルです。[tract] は tract 番号で、今回は 0000 になります。

モザイクを行う別のコマンドとして、jointcalがあります (mosaic.py (jointcal.py)を行なった場合はjointcal.py (mosaic.py)を行わないでください)。 hscpipeはmosaic.pyからjointcal.pyに移行しているところですが、現在のところ、jointcalでは大きな枚数のスタックを行う際に問題が生じることが報告されており、mosaic.pyを用いることを推奨しております。

jointcal.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id field=HDFN filter=HSC-Y tract=0
        # filter を書き換えて計 3 回走らせます。
        # jointcal.py はシングルプロセスで動きますので、--batch-type, --cores の指定をしてはいけません。
        

天体データの重ね合わせ (スタック)

モザイクによって生成されたデータを用いて、天球面座標を平面座標に投影したデータを生成し（warp）、観測された全ショットの積分を行います（coadd）。

coaddDriver.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id filter=HSC-Y tract=0 \ # filter はここと
        --selectId field=HDFN filter=HSC-Y \ # ここにあるので注意！！
        --batch-type=smp --cores=8

filter を HSC-I, HSC-Z として同様に走らせて下さい。上記のコマンドでは、filter 情報は二箇所に入りますので両方書き換えて下さい。 Coadd 画像は /work/[User Name]/HSC/rerun/[rerun]/deepCoadd-results/[filter]/[tract]/[patch] の下の calexp-[filter]-[tract]-[patch].fits です。また、warp 画像 (warp-[filter]-[tract]-[patch]-[visit].fits) は/work/[User Name]/HSC/rerun/[rerun]/deepCoadd/[filter]/[tract]/[patch]に生成されます。 warp 画像が足し合わせる前の画像ですので、もし coadd 画像におかしな点が見られた場合は、まず warp 画像を確認してみましょう。 また、coaddDriver.py 以降のプロセッシングではデータは上書きされず、スキップされるので、再度解析を行う場合は古い calexp 画像や warp 画像は移動させてください。

さらに、/work/[User Name]/HSC/rerun/[rerun]/deepCoadd-results/[filter]/[tract]/[patch] には検出された天体のカタログ (det-[filter]-[tract]-[patch].fits) と sky 引きのパターン (bkgd-[filter]-[tract]-[patch].fits) ができます。この det- カタログが次のステップのマルチバンド解析で使われることになります。 coaddDriver.py後の画像を一度に表示させたい場合は、i バンドの場合であれば、

cd /work/[User Name]/HSC/rerun/tutorial/deepCoadd-results/HSC-I/0
ds9   -mosaic   wcs    */calexp-*.fits

で確認することが可能です (ただし、全 104 CCD でこれをやるとかなりのメモリを消費しますのでご注意ください)。

今回は中心の 6 CCD のみを使用しましたが、全ての CCD を使用した場合はこのようになります。

–> こちら

以上でcoaddDriver.pyが完了です。

11/8の実習・講義内容

まず最初に、coaddDriver.pyの残りがある場合は多波長解析システムからコピーしてから、マルチバンド解析(multiBandDriver.py)を走らせてもらいます。

(時間がかかるので、ひとまず走らせてもらいます。順調ならmultiBandDriver.pyは昼過ぎに終わります)

コマンドを叩いたら、別のターミナルを開いて、fits viewerを使ったカタログ閲覧、butlerを使った解析(画像表示)、ds9で擬似カラー画像作成を行います。

休憩を挟んで、11時頃からマルチバンド解析についての座学を行います。

お昼休みのあと13時頃から大規模解析システムの紹介を行います。

最後に色等級図作成を行ってもらって講習会は終了となる予定です。

# coaddDriver.pyの結果のコピーの仕方。昨日コマンドを打ち終わらなかったフィルターについて行ってください。大体5分くらいかかります。
cd /work/[User name]
cp -r /lfs01/kubomr/HSC/rerun/tutorial/deepCoadd/HSC-I HSC/rerun/tutorial/deepCoadd/HSC-I
cp -r /lfs01/kubomr/HSC/rerun/tutorial/deepCoadd-results/HSC-I HSC/rerun/tutorial/deepCoadd-results/HSC-I

マルチバンド解析

Coaddによって生成されたそれぞれのフィルタの天体カタログをマージして、マルチバンドで測光しなおし、新たにカタログを生成します。

multiBandDriver.pyには大体2時間40分ほどかかります。

multiBandDriver.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id tract=0 filter=HSC-I^HSC-Z^HSC-Y \ # 使うフィルタを全選択。
        --batch-type=smp --cores=8

まず最初に各 filter で検出された天体カタログを集めて天体の位置情報のマージカタログを作成します（mergeDet）。この時、デブレンド（2 つの天体が近接していて重なって見える時にその 2 つの天体を分離するプロセス）も行われます。 mergeDet カタログをもとに各 filter の coadd 画像に戻って天体の測光を行います（meas）。 その後、各 filter の meas カタログで得られた天体の位置情報を再度集め直して、全ての filter で共通の座標情報を持つ天体カタログを作成します（ref）。 この ref カタログの位置情報をもとに天体の測光を filter 毎に改めて行い、forced_src- カタログが生成されます。

生成されるファイルを以下にまとめます。

/work/[User Name]/HSC/rerun/[rerun]/deepCoadd-results/merged/[tract]/[patch] ;

• 天体の位置情報のマージカタログ : mergeDet-[tract]-[patch].fits
• meas カタログをもとに生成される天体の位置情報カタログ : ref-[tract]-[patch].fits

/work/[User Name]/HSC/rerun/[rerun]/deepCoadd-results/[filter]/[tract]/[patch] ;

• mergeDet をもとに測光した天体カタログ : meas-[filter]-[tract]-[patch].fits
• 外部カタログとマッチした天体カタログ : srcMatch-[filter]-[tract]-[patch].fits
• srcMatch に座標情報などを追加し、カラムに展開したもの : srcMatchFull-[filter]-[tract]-[patch].fits
• 最終天体カタログ : forced_src-[filter]-[tract]-[patch].fits

各 CCD 画像, カタログを較正しなおす

最後にモザイクによって決められた等級原点と WCS を、各 CCD 画像、カタログに反映させましょう。 WCS とゼロ点の較正された画像/カタログが不要な場合は行う必要はありませんが、講習会なので一通りやっていきます。 現在のhscpipeではcalibrateExposure.py, calibrateCatalog.pyは完全には動作しません。カタログは生成されますが、半分が空になってしまうようです。

# 画像を較正しなおす
calibrateExposure.py /work/[User Name]/HSC\
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id field=HDFN tract=0 \
        -j 8
# -j : スレッド数を指定。-j 8 とすると、8 スレッドで、という意味になる。


# カタログを較正しなおす
calibrateCatalog.py /work/[User Name]/HSC \
        --calib /work/[User Name]/HSC/CALIB \
        --rerun tutorial \
        --id field=HDFN tract=0 \
        -j 8

calibrateExposure.py が終了すると、/work/[User Name]/HSC/rerun/[rerun]/[pointing]/[filter]/corr/[tract] 以下に CALEXP-[visit]-[ccd].fits というファイルが生成されます。 また、calibrateCatalog.py で較正されたカタログは、/work/[User Name]/HSC/rerun/[rerun]/[pointing]/[filter]/output/[tract] 以下の CALSRC-[visit]-[ccd].fits です。

以上で全ての解析が終了です。

疑似カラー画像作成

それでは、解析で出来上がったデータを元に、疑似カラー画像を作成してみましょう。

ds9のバージョンによって動作は異なります。

1. 必要なスクリプトをコピーしてきます。

   cd /work/[User Name]/
   wget http://hsc.mtk.nao.ac.jp/HSC_Training_download/script/2019/imageStitcher1.py
   # スクリプトの中身は各自で確認下さい。
   
   # 日本語が文字化けしている場合は、
   # メニューから「Terminal/Preferences/Encoding/Unicord(UTF-8)」を選び、export LC_ALL=ja_JP.UTF-8 をしてください。
   

2. パッチに分かれていた画像を1つの大きな画像にします。

   python imageStitcher1.py -o HSC-Y.fits \
   HSC/rerun/tutorial/deepCoadd-results/HSC-Y/0/*/calexp-*.fits
   
   python imageStitcher1.py -o HSC-Z.fits \
   HSC/rerun/tutorial/deepCoadd-results/HSC-Z/0/*/calexp-*.fits
   
   python imageStitcher1.py -o HSC-I.fits \
   HSC/rerun/tutorial/deepCoadd-results/HSC-I/0/*/calexp-*.fits
   
   # または、
   
   for f in Y Z I ; do python imageStitcher1.py -o HSC-$f.fits \
   HSC/rerun/tutorial/deepCoadd-results/HSC-$f/0/*/calexp-*.fits ; done
   

3. ds9のカラーモードで表示します。

   ds9 -rgb -red HSC-Y.fits -green HSC-Z.fits -blue HSC-I.fits
   

4. 色を調整します。

• 各フィルタごとに次の操作をして色を調整します

  1. メニューから Frame/RGB... を選ぶ。

     • 「R」「G」「B」からチャンネルを選ぶ。

  2. メニューから Scale/Scale Parameters...(スケール/スケールの設定) を選ぶ。

     • スケールの最小値、最大値を調整する。

  3. メニューから Color/Colormap Parameters... を選ぶ。

     • Contrast, Biasを調整する。

• （おすすめの調整の仕方）

  • フィルタごとに

    1. メニューから Scale/Log を選ぶ。

    2. メニューから Scale/Scale Parameters を選ぶ。

       • 最小値を0、最大値を100くらいにする。

5. png形式で保存します。

   1. メニューから File/Save Image/png を選ぶ。
   2. HDFN.png の名前で保存する。

6. 保存した画像を確認します。

   eog HDFN.png
   
   # または
   
   firefox HDFN.png
   

色等級図作成

最後に、カタログ情報から色等級図を描いてみましょう。

ここでは、butlerと呼ばれている、様々なタイプ（画像、カタログ）のデータを検索したりロードしたりするためのPipeline toolを使って色等級図の作成を行います。 multiband解析の結果を待つ間にcoaddまでの結果を使ってbutlerを使ってみましょう-> butlerで画像上にカタログを表示。 butlerについてさらに知りたい方は こちらのページ を参照してください。

今回 2. で使用する catalog_all.py では、bulter でカタログ内の欲しいデータをどう取ってくるかが書かれています。これを応用することで自分の欲しいデータを正しく取得することができるようになりますので、 是非一度読んでみて下さい。

1. 必要なスクリプトをコピーします。

   cd /work/[User Name]
   wget http://hsc.mtk.nao.ac.jp/HSC_Training_download/script/2019/catalog_all_pipe7.py
   wget http://hsc.mtk.nao.ac.jp/HSC_Training_download/script/2019/catalog.sh
   wget http://hsc.mtk.nao.ac.jp/HSC_Training_download/script/2019/fig_izy.py
   

2. 各バンドで作成したカタログから、butler を用いて天体の ID と PSF flux を取得し、一つのファイルにまとめます。

   python catalog_all_pipe7.py /work/[User Name]/HSC/rerun/tutorial 0 HSC-I > HSC-I.txt
   python catalog_all_pipe7.py /work/[User Name]/HSC/rerun/tutorial 0 HSC-Z > HSC-Z.txt
   python catalog_all_pipe7.py /work/[User Name]/HSC/rerun/tutorial 0 HSC-Y > HSC-Y.txt
   sh catalog.sh
   

3. 一つにまとめたファイルを使って、色等級図を表示します。

   python fig_izy.py
   

4. 表示された図を保存します。

   1. save the figureのボタンを押す。
   2. 保存する名前を決めてsaveボタンを押す。

5. 保存した図を確認します。

   eog hoge.png
   
   # または
   
   firefox hoge.png
   

以上で実習の内容はすべて終了です。お疲れ様でした。