从图像中堆叠星 PSF;对齐子像素中心
Stacking star PSFs from an Image; aligning sub-pixel centers
我有一个 (1727,1853) 大小的数组(图像),我在其中识别了星星来模拟点扩散函数。数组的每个索引对应一个图像坐标,但是,每个星的质心由一个子像素坐标给出。我必须执行以下操作
制作每颗星星的二维切片。我已经使用 numpy 的数组切片完成了这项工作。然而,它按索引切片,我有亚像素质心坐标,因此我制作的任何类型的切片都会使星星偏离中心。
在我制作每颗星的二维切片后,我必须将这些阵列堆叠在一起以制作点扩散函数模型。这很简单,只要每个星星的子像素中心对齐即可。
我的问题是对齐这些子像素坐标并将每个 2D 切片堆叠在一起的最有效(和正确)方法是什么?
我希望这是清楚的。任何帮助将非常感激。下面是其中一颗星星的二维切片(不是很好),但是它偏离了中心,因为按索引的 numpy 切片和星星的质心具有子像素坐标。
您可以表达每个 'slice' 中像素相对于星形质心的中心坐标,然后计算加权二维直方图。
首先,一些示例数据:
import numpy as np
from matplotlib import pyplot as plt
# pixel coordinates (integer)
x, y = np.mgrid[:100, :100]
# centroids (float)
cx, cy = np.random.rand(2, 9) * 100
# a Gaussian kernel to represent the PSF
def gausskern(x, y, cx, cy, sigma):
return np.exp(-((x - cx) ** 2 + (y - cy) ** 2) / (2 * sigma ** 2))
# (nstars, ny, nx)
stars = gausskern(x[None, ...], y[None, ...],
cx[:, None, None], cy[:, None, None], 10)
# add some noise for extra realism
stars += np.random.randn(*stars.shape) * 0.5
fig, ax = plt.subplots(3, 3, figsize=(5, 5))
for ii in xrange(9):
ax.flat[ii].imshow(stars[ii], cmap=plt.cm.hot)
ax.flat[ii].set_axis_off()
fig.tight_layout()
加权二维直方图:
# (nstars, ny, nx) pixel coordinates relative to each centroid
dx = cx[:, None, None] - x[None, ...]
dy = cy[:, None, None] - y[None, ...]
# 2D weighted histogram
bins = np.linspace(-50, 50, 100)
h, xe, ye = np.histogram2d(dx.ravel(), dy.ravel(), bins=bins,
weights=stars.ravel())
fig, ax = plt.subplots(1, 1, subplot_kw={'aspect':'equal'})
ax.hold(True)
ax.pcolormesh(xe, ye, h, cmap=plt.cm.hot)
ax.axhline(0, ls='--', lw=2, c='w')
ax.axvline(0, ls='--', lw=2, c='w')
ax.margins(x=0, y=0)
我有一个 (1727,1853) 大小的数组(图像),我在其中识别了星星来模拟点扩散函数。数组的每个索引对应一个图像坐标,但是,每个星的质心由一个子像素坐标给出。我必须执行以下操作
制作每颗星星的二维切片。我已经使用 numpy 的数组切片完成了这项工作。然而,它按索引切片,我有亚像素质心坐标,因此我制作的任何类型的切片都会使星星偏离中心。
在我制作每颗星的二维切片后,我必须将这些阵列堆叠在一起以制作点扩散函数模型。这很简单,只要每个星星的子像素中心对齐即可。
我的问题是对齐这些子像素坐标并将每个 2D 切片堆叠在一起的最有效(和正确)方法是什么?
我希望这是清楚的。任何帮助将非常感激。下面是其中一颗星星的二维切片(不是很好),但是它偏离了中心,因为按索引的 numpy 切片和星星的质心具有子像素坐标。
您可以表达每个 'slice' 中像素相对于星形质心的中心坐标,然后计算加权二维直方图。
首先,一些示例数据:
import numpy as np
from matplotlib import pyplot as plt
# pixel coordinates (integer)
x, y = np.mgrid[:100, :100]
# centroids (float)
cx, cy = np.random.rand(2, 9) * 100
# a Gaussian kernel to represent the PSF
def gausskern(x, y, cx, cy, sigma):
return np.exp(-((x - cx) ** 2 + (y - cy) ** 2) / (2 * sigma ** 2))
# (nstars, ny, nx)
stars = gausskern(x[None, ...], y[None, ...],
cx[:, None, None], cy[:, None, None], 10)
# add some noise for extra realism
stars += np.random.randn(*stars.shape) * 0.5
fig, ax = plt.subplots(3, 3, figsize=(5, 5))
for ii in xrange(9):
ax.flat[ii].imshow(stars[ii], cmap=plt.cm.hot)
ax.flat[ii].set_axis_off()
fig.tight_layout()
加权二维直方图:
# (nstars, ny, nx) pixel coordinates relative to each centroid
dx = cx[:, None, None] - x[None, ...]
dy = cy[:, None, None] - y[None, ...]
# 2D weighted histogram
bins = np.linspace(-50, 50, 100)
h, xe, ye = np.histogram2d(dx.ravel(), dy.ravel(), bins=bins,
weights=stars.ravel())
fig, ax = plt.subplots(1, 1, subplot_kw={'aspect':'equal'})
ax.hold(True)
ax.pcolormesh(xe, ye, h, cmap=plt.cm.hot)
ax.axhline(0, ls='--', lw=2, c='w')
ax.axvline(0, ls='--', lw=2, c='w')
ax.margins(x=0, y=0)