1

等高線図を作成し、画像に等高線を重ねる必要があります。aplpy ライブラリを使用して、天体画像の輪郭を重ねました。APlpy Web サイト ( https://github.com/aplpy/aplpy-examples/tree/master/data ) で 2MASS データをダウンロードし、次のコードを記述しました。

import numpy as np
import aplpy
import atpy
from pyavm import AVM
import scipy
import matplotlib.pyplot as plt
import scipy.ndimage
from matplotlib import cm
import montage_wrapper
import matplotlib.colors as colors
from scipy import stats

def get_contour_verts(cn):
    contours = []
    # for each contour line
    for cc in cn.collections:
        paths = []
        # for each separate section of the contour line
        for pp in cc.get_paths():
            xy = []
            # for each segment of that section
            for vv in pp.iter_segments():
                xy.append(vv[0])
            paths.append(np.vstack(xy))
        contours.append(paths)
    return contours

# Convert all images to common projection
aplpy.make_rgb_cube(['2MASS_h.fits', '2MASS_j.fits', '2MASS_k.fits'], '2MASS_rgb.fits')

# Customize it
aplpy.make_rgb_image('2MASS_rgb.fits','2MASS_rgb_contour.png',embed_avm_tags=True)

# Launch APLpy figure of 2D cube
img = aplpy.FITSFigure('2MASS_rgb_contour.png')
img.show_rgb()

# Modify the tick labels for precision and format
img.tick_labels.set_xformat('hhmm')
img.tick_labels.set_yformat('ddmm')
# Move the tick labels
img.tick_labels.set_xposition('top')
img.tick_labels.set_yposition('right')
img.tick_labels.set_font(size='small', family='sans-serif', style='normal')

data = scipy.loadtxt('sources.txt')
m1=data[:,0]
m2=data[:,1]
xmin = m1.min()
xmax = m1.max()
ymin = m2.min()
ymax = m2.max()
#Gridding the data 

X, Y = np.mgrid[xmin:xmax:100j, ymin:ymax:100j]
positions = np.vstack([X.ravel(), Y.ravel()])
values = np.vstack([m1, m2])
kernel = stats.gaussian_kde(values)
Z = np.reshape(kernel(positions).T, X.shape)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.imshow(np.rot90(Z), cmap=plt.cm.gist_earth_r, extent=[xmin, xmax, ymin, ymax])
ax.plot(m1, m2, 'k.', markersize=2)
ax.set_xlim([xmin, xmax])
ax.set_ylim([ymin, ymax])

CS = plt.contour(X, Y, Z)
Contour_arrays=get_contour_verts(CS)

#Adding contour plots
for contours_at_level in Contour_arrays:
    radec = [img.pixel2world(*verts.T) for verts in contours_at_level]
    new_radec=[]
    for coor in radec:
        #new_radec.append(np.column_stack((coor[0], coor[1])))
        new_radec.append(np.vstack((coor[0], coor[1])))
    print new_radec
    img.show_lines(new_radec,color='red', zorder=100)

img.save('tutorial.png')

今でも全然通用しないようです。

4

1 に答える 1

2

直接「輪郭を FITS ファイルとして保存」することはできませんが、他にも試すことができる方法があります。

matplotlib._cntrここで説明されているようにツールを使用できます: Python: matplotlib.pyplot.contour() から等高線を見つけて Figure 座標の端点を取得し、WCS を使用してピクセル座標とワールド座標の間で変換します。 s には、それぞれ 2 つの配列を受け入れるaplpy.FITSFigure便利な関数F.world2pixelおよびがあります。F.pixel2world

F.pixel2world(arange(5),arange(5))

したがって、ウィンドウに表示されているグリッドと同じFITSFigureグリッドで作業している場合は、ポイントをワールド座標に変換し、次のようにプロットできますshow_lines

ra,dec = F.pixel2world(xpix,ypix)
F.show_lines([[ra,dec]])

または、より現実的な輪郭のケースでは、リンクされた記事からコードをコピーします。

import numpy as np
import aplpy

def get_contour_verts(cn):
    contours = []
    # for each contour line
    for cc in cn.collections:
        paths = []
        # for each separate section of the contour line
        for pp in cc.get_paths():
            xy = []
            # for each segment of that section
            for vv in pp.iter_segments():
                xy.append(vv[0])
            paths.append(np.vstack(xy))
        contours.append(paths)

    return contours

# This line is copied from the question's code, and assumes that has been run previously
CS = plt.contour(Xgrid, Ygrid, Hsmooth,levels=loglvl,extent=extent,norm=LogNorm())

contours = get_contour_verts(CS) # use the linked code


gc = aplpy.FITSFigure('2MASS.fits',figsize=(10,9))

# each level comes with a different set of vertices, so you have to loop over them
for contours_at_level in Contour_arrays:
    clines = [cl.T for cl in contours_at_level]
    img.show_lines(clines,color='red', zorder=100)

それでも最も簡単な方法は、上記の等高線を生成する元の FITS ファイルがある場合は、そのファイルを直接使用することです。または、FITS ファイルがない場合は、独自のヘッダーを作成して pyfits で作成できます。

于 2014-02-07T23:34:55.523 に答える