画像に関する Tips¶

観測領域全面で 1 枚の fits 画像を作りたい場合¶

HSC の視野が大きいため、HSC pipeline で生成される coadd 画像は patch 毎に出力されます。しかし、場合によっては 1 枚の大きな fits 画像を生成して使う場合があります。例えば、次の python スクリプトを使うと全ての観測領域を含む large.fits を生成してくれます。

# python スクリプトの実行例
# 1枚にまとめたい coadd 画像を指定します。
# hscPipe 3.x
python make-large-image.py $home/hsc/rerun/dith_16h_test/deepCoadd/HSC-I/0/*.fits

# hscPipe 4.x
python make-large-image.py $home/hsc/rerun/dith_16h_test/deepCoadd/HSC-I/0/*/calexp*.fits

スクリプトの中身は以下の通りです。

#!/usr/bin/env python
# encoding: utf-8

# AFW (Application FrameWork) のうち画像に関するライブラリ
import lsst.afw.image as afwImage
# AFW のうちジオメトリに関するライブラリ
import lsst.afw.geom  as afwGeom

import sys

def main():

        if len(sys.argv) <= 1:
                print "一つ以上のパッチを指定してください"
                return 1

        largeImage = makeLargeImageFromPatches(sys.argv[1:])

        print "保存中 (しばらくお待ちください)..."
        largeImage.writeFits("large.fits")


def makeLargeImageFromPatches(filenames):

        # 画像のメトリックの情報だけ取得する
        bboxes = []
        for filename in filenames:
                print "バウンディングボックスを取得:", filename, "..."
                bboxes.append(
                        # 親 (= トラクト) におけるこのパッチのバウンディング・ボックスを取得
                        afwImage.ExposureF(filename).getBBox(afwImage.PARENT)
                )

        # 全ての BBox を含むような BBox を求める
        #     getMin/Max 系のメソッドの戻り値はインクルーシヴ
        #     (領域 = [min, max]) cf. getBegin/End (下で使用)
        x0 = min(bbox.getMinX() for bbox in bboxes)
        y0 = min(bbox.getMinY() for bbox in bboxes)
        x1 = max(bbox.getMaxX() for bbox in bboxes)
        y1 = max(bbox.getMaxY() for bbox in bboxes)

        width  = x1 - x0 + 1
        height = y1 - y0 + 1

        # 大きな画像を用意する
        print "画像を用意: ({}, {})".format(width, height)
        exposure = afwImage.ExposureF(width, height)

        # トラクト内における画像の位置を設定
        exposure.setXY0(afwGeom.Point2I(x0, y0))

        # Exposure から画像部分だけ抜き出す
        largeImage = exposure.getMaskedImage()

        # largeImage にパッチを貼り付けていく
        for filename in filenames:
                print "パッチを取得:", filename, "..."

        # 再びパッチをロード
                patch = afwImage.ExposureF(filename)
                bbox = patch.getBBox(afwImage.PARENT)

                # このパッチの、largeImage の中における位置
                #     getBegin/End 系のメソッドの戻り値は one past the end
                #     (領域 = [begin,end) ) cf. getMin/Max (上で使用)
                xp0 = bbox.getBeginX() - x0
                yp0 = bbox.getBeginY() - y0
                xp1 = bbox.getEndX  () - x0
                yp1 = bbox.getEndY  () - y0

                # パッチを張りつける ([]の中の記法については numpy のドキュメントなどを参考に)
                largeImage[xp0:xp1, yp0:yp1] = patch.getMaskedImage()[:]

                # 全てのパッチが同じ WCS その他の情報をもっているものだと盲信して
                # どれか一個を大きな画像にセットする
                if not exposure.hasWcs() and patch.hasWcs():
                        exposure.setCalib   (patch.getCalib   ())
                        exposure.setFilter  (patch.getFilter  ())
                        exposure.setMetadata(patch.getMetadata())
                        exposure.setWcs     (patch.getWcs     ())

        return exposure


if __name__ == "__main__":
        main()

例えば、このスクリプトで生成された large.fits を ds9 で見ると図1 のようになります。

図1 large.fits の画像