分段仿射变换+扭曲输出看起来很奇怪
Piecewise Affine Transform+warp output looks strange
我有一张图像，我正尝试使用 skimage.PiecewiseAffineTransform 和 skimage.warp 进行变形。我有一组控制点 (true) 映射到一组新的控制点 (ideal)，但扭曲没有返回我期望的结果。
在这个例子中，我有一个简单的波长梯度，我试图 'straighten out' 到列中。 （您可能会问我为什么要查找轮廓和插值，但那是因为我实际上是将这段代码应用到更复杂的用例中。我只是想重现这个导致相同奇怪输出的简单示例的所有代码。）
为什么我的输出图像只是将输入图像扭曲成正方形并插入？我正在使用 Python 2.7.12 和 matplotlib 1.5.1。这是代码。
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform
true = np.array([range(i,i+10) for i in range(20)])
ideal = np.array([range(10)]*20)
# Find contours of ideal and true images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
for lam in ideal[0]:
try:
# Get the isowavelength contour in the true and ideal images
tc = measure.find_contours(true, lam)[0]
ic = measure.find_contours(ideal, lam)[0]
nc = np.ones(ic.shape)
# Use the y coordinates of the ideal contour
nc[:, 0] = ic[:, 0]
# Interpolate true contour onto ideal contour y axis so there are the same number of points
nc[:, 1] = np.interp(ic[:, 0], tc[:, 0], tc[:, 1])
# Add the control points to the appropriate list
true_control_pts.append(nc.tolist())
ideal_control_pts.append(ic.tolist())
except (IndexError,AttributeError):
pass
true_control_pts = np.array(true_control_pts)
ideal_control_pts = np.array(ideal_control_pts)
length = len(true_control_pts.flatten())/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]).astype(np.int32)
plt.figure()
plt.imshow(image, origin='lower', interpolation='none')
plt.title('Input image')
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none')
plt.title('Should be parallel lines')
输出如下：
如有任何帮助，我们将不胜感激！
我会把我所有的都给你。我不认为我能够用给出的测试数据确定这一点。在这么小的图像中将 45 度角映射到一条直线意味着需要发生很多运动，而不是大量的数据作为基础。我确实发现了一些在下面的代码中修复的特定错误，标记为 /* */（因为这不是您通常在 Python 文件中看到的东西，该标记应该突出 :)）。
请在您的真实数据上试用此代码，如果有效请告诉我！对于这个输入数据集，有一些非零输出。
主要问题是：
interp数据需要排序- 许多控制点都是无意义的，因为此数据集中没有足够的数据（例如，将大跨度的列映射到单个点）
- 图像中的值范围已关闭（因此您的原始输出实际上几乎都是 0 - 色块的值在 1e-9 左右）。
我添加的最重要的东西是 3D 图，显示了 "true" 控制点如何映射到 "ideal" 控制点。这为您提供了一个调试工具来显示您的控制点映射是否符合您的预期。正是这个情节让我遇到了 interp 问题。
顺便说一句，请使用 "before" 和 "after" 之类的名称，而不是 ideal 和 true :)。试图记住哪个是哪个至少让我绊倒了一次。
第一次尝试
import pdb
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform, img_as_float
from mpl_toolkits.mplot3d import Axes3D # /**/
#/**/
# From by
# https://whosebug.com/users/1355221/dansalmo
def flatten_list(L):
for item in L:
try:
for i in flatten_list(item): yield i
except TypeError:
yield item
#end flatten_list
true_input = np.array([range(i,i+10) for i in range(20)]) # /** != True **/
ideal = np.array([range(10)]*20)
#pdb.set_trace()
# Find contours of ideal and true_input images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
OLD_true=[] # /**/ for debugging
OLD_ideal=[]
for lam in [x+0.5 for x in ideal[0]]: # I tried the 0.5 offset just to see,
try: # but it didn't make much difference
# Get the isowavelength contour in the true_input and ideal images
tc = measure.find_contours(true_input, lam)[0]
ic = measure.find_contours(ideal, lam)[0]
nc = np.zeros(ic.shape) # /** don't need ones() **/
# Use the y /** X? **/ coordinates of the ideal contour
nc[:, 0] = ic[:, 0]
# Interpolate true contour onto ideal contour y axis so there are the same number of points
# /** Have to sort first - https://docs.scipy.org/doc/numpy/reference/generated/numpy.interp.html#numpy-interp **/
tc_sorted = tc[tc[:,0].argsort()]
# /** Thanks to by
# https://whosebug.com/users/208339/steve-tjoa **/
nc[:, 1] = np.interp(ic[:, 0], tc_sorted[:, 0], tc_sorted[:, 1],
left=np.nan, right=np.nan)
# /** nan: If the interpolation is out of range, we're not getting
# useful data. Therefore, flag it with a nan. **/
# /** Filter out the NaNs **/
# Thanks to by
# https://whosebug.com/users/449449/eumiro
#pdb.set_trace()
indices = ~np.isnan(nc).any(axis=1)
nc_nonan = nc[indices]
ic_nonan = ic[indices]
# Add the control points to the appropriate list.
# /** Flattening here since otherwise I wound up with dtype=object
# in the numpy arrays. **/
true_control_pts.append(nc_nonan.flatten().tolist())
ideal_control_pts.append(ic_nonan.flatten().tolist())
OLD_true.append(nc) # /** for debug **/
OLD_ideal.append(ic)
except (IndexError,AttributeError):
pass
#pdb.set_trace()
# /** Make vectors of all the control points. **/
true_flat = list(flatten_list(true_control_pts))
ideal_flat = list(flatten_list(ideal_control_pts))
true_control_pts = np.array(true_flat)
ideal_control_pts = np.array(ideal_flat)
# Make the vectors 2d
length = len(true_control_pts)/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
#pdb.set_trace()
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]) / 30.0 # /**.astype(np.int32)**/
# /** You don't want int32 images! See
# http://scikit-image.org/docs/dev/user_guide/data_types.html .
# Manually rescale the image to [0.0,1.0] by dividing by 30. **/
image_float = img_as_float(image) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title('Input image')
fig.show() # /** I needed this on my test system **/
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# /** since we started with float, and this is float, too, the two are
# comparable. **/
pdb.set_trace()
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none') # /**note: float**/
plt.title('Should be parallel lines')
fig.show()
# /** Show the control points.
# The z=0 plane will be the "true" control points (before), and the
# z=1 plane will be the "ideal" control points (after). **/
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
fig.show()
for rowidx in range(length):
ax.plot([true_control_pts[rowidx,0], ideal_control_pts[rowidx,0]],
[true_control_pts[rowidx,1], ideal_control_pts[rowidx,1]],
[0,1])
input() # /** because I was running from the command line **/
第二次尝试
越来越近：
下面是看起来更有希望的控制点映射视图：
您可以看到它试图稍微旋转图像，这正是我希望从此数据集中看到的结果。
代码
import pdb
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform, img_as_float
from mpl_toolkits.mplot3d import Axes3D # /**/
#/**/
# From by
# https://whosebug.com/users/1355221/dansalmo
def flatten_list(L):
for item in L:
try:
for i in flatten_list(item): yield i
except TypeError:
yield item
#end flatten_list
#/**/
# Modified from by
# https://whosebug.com/users/2588210/christian-k
def equispace(data, npts):
x,y = data.T
xd =np.diff(x)
yd = np.diff(y)
dist = np.sqrt(xd**2+yd**2)
u = np.cumsum(dist)
u = np.hstack([[0],u])
t = np.linspace(0,u.max(),npts)
xn = np.interp(t, u, x)
yn = np.interp(t, u, y)
return np.column_stack((xn, yn))
true_input = np.array([range(i,i+10) for i in range(20)]) # /** != True **/
ideal = np.array([range(10)]*20)
#pdb.set_trace()
# Find contours of ideal and true_input images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
OLD_true=[] # /**/ for debugging
OLD_ideal=[]
for lam in [x+0.5 for x in ideal[0]]: # I tried the 0.5 offset just to see,
try: # but it didn't make much difference
# Get the isowavelength contour in the true_input and ideal images
tc = measure.find_contours(true_input, lam)[0]
# /** So this might not have very many numbers in it. **/
ic = measure.find_contours(ideal, lam)[0]
# /** CAUTION: this is assuming the contours are going the same
# direction. If not, you'll need to make it so. **/
#nc = np.zeros(ic.shape) # /** don't need ones() **/
# /** We just want to find points on _tc_ to match _ic_. That's
# interpolation _within_ a curve. **/
#pdb.set_trace()
nc_by_t = equispace(tc,ic.shape[0])
ic_by_t = equispace(ic,ic.shape[0])
### /** Not this **/
## Use the y /** X? **/ coordinates of the ideal contour
#nc[:, 0] = ic[:, 0]
#
## Interpolate true contour onto ideal contour y axis so there are the same number of points
#
## /** Have to sort first - https://docs.scipy.org/doc/numpy/reference/generated/numpy.interp.html#numpy-interp **/
#tc_sorted = tc[tc[:,0].argsort()]
# # /** Thanks to by
# # https://whosebug.com/users/208339/steve-tjoa **/
#
#nc[:, 1] = np.interp(ic[:, 0], tc_sorted[:, 0], tc_sorted[:, 1],
# left=np.nan, right=np.nan)
# # /** nan: If the interpolation is out of range, we're not getting
# # useful data. Therefore, flag it with a nan. **/
# /** Filter out the NaNs **/
# Thanks to by
# https://whosebug.com/users/449449/eumiro
#pdb.set_trace()
#indices = ~np.isnan(nc).any(axis=1)
#nc_nonan = nc[indices]
#ic_nonan = ic[indices]
#
# Add the control points to the appropriate list.
## /** Flattening here since otherwise I wound up with dtype=object
## in the numpy arrays. **/
#true_control_pts.append(nc_nonan.flatten().tolist())
#ideal_control_pts.append(ic_nonan.flatten().tolist())
#OLD_true.append(nc) # /** for debug **/
#OLD_ideal.append(ic)
true_control_pts.append(nc_by_t)
ideal_control_pts.append(ic_by_t)
except (IndexError,AttributeError):
pass
pdb.set_trace()
# /** Make vectors of all the control points. **/
true_flat = list(flatten_list(true_control_pts))
ideal_flat = list(flatten_list(ideal_control_pts))
true_control_pts = np.array(true_flat)
ideal_control_pts = np.array(ideal_flat)
# Make the vectors 2d
length = len(true_control_pts)/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
#pdb.set_trace()
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]) / 30.0 # /**.astype(np.int32)**/
# /** You don't want int32 images! See
# http://scikit-image.org/docs/dev/user_guide/data_types.html .
# Manually rescale the image to [0.0,1.0] by dividing by 30. **/
image_float = img_as_float(image) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title('Input image')
fig.show() # /** I needed this on my test system **/
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# /** since we started with float, and this is float, too, the two are
# comparable. **/
pdb.set_trace()
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none') # /**note: float**/
plt.title('Should be parallel lines')
fig.show()
# /** Show the control points.
# The z=0 plane will be the "true" control points (before), and the
# z=1 plane will be the "ideal" control points (after). **/
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
fig.show()
for rowidx in range(length):
ax.plot([true_control_pts[rowidx,0], ideal_control_pts[rowidx,0]],
[true_control_pts[rowidx,1], ideal_control_pts[rowidx,1]],
[0,1])
input() # /** because I was running from the command line **/
说明
以不同的方式考虑插值：X 和 Y 坐标实际上并不是您想要映射的，对吗？我想你想要映射的是沿着轮廓的距离，这样对角线轮廓就会被拉伸成垂直轮廓。这就是这些行的作用：
nc_by_t = equispace(tc,ic.shape[0])
ic_by_t = equispace(ic,ic.shape[0])
我们沿每个等高线创建 ic.shape[0] 个等距点，然后将这些点相互映射。 equispace 修改自 this answer。因此，这会拉伸较短的轮廓以适应较长的轮廓，无论哪种方式，都由 ic 中的轮廓长度定义的点数。实际上，您可以使用此技术使用任意数量的点；我刚刚测试了 100 点，得到了基本相同的结果。
再次在您的实际数据上进行测试。想想你的插值参考究竟是什么。它真的是 X 或 Y 坐标吗？是沿着轮廓的距离吗？等高线的百分比（如上）?
好的，多一个主意
对于这个特定的测试用例，更多的数据会有帮助吗？是的！
我使用较大的图像来确定控制点，然后在该区域的中心映射较小的图像。
代码
（此时一团糟 - 见 =========== 标记）
import pdb
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform, img_as_float
from mpl_toolkits.mplot3d import Axes3D # /**/
#/**/
def test_figure(data,title):
image_float = img_as_float(data) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title(title)
fig.show()
#/**/
# From by
# https://whosebug.com/users/1355221/dansalmo
def flatten_list(L):
for item in L:
try:
for i in flatten_list(item): yield i
except TypeError:
yield item
#end flatten_list
#/**/
# Modified from by
# https://whosebug.com/users/2588210/christian-k
def equispace(data, npts):
x,y = data.T
xd =np.diff(x)
yd = np.diff(y)
dist = np.sqrt(xd**2+yd**2)
u = np.cumsum(dist)
u = np.hstack([[0],u])
t = np.linspace(0,u.max(),npts)
xn = np.interp(t, u, x)
yn = np.interp(t, u, y)
return np.column_stack((xn, yn))
# ====================== BIGGER
true_input = np.array([range(i,i+20) for i in range(30)]) # /** != True **/
ideal = np.array([range(20)]*30)
# ======================
test_figure(true_input / 50.0, 'true_input')
test_figure(ideal / 20.0, 'ideal')
#pdb.set_trace()
# Find contours of ideal and true_input images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
OLD_true=[] # /**/ for debugging
OLD_ideal=[]
for lam in [x+0.5 for x in ideal[0]]: # I tried the 0.5 offset just to see,
try: # but it didn't make much difference
# Get the isowavelength contour in the true_input and ideal images
tc = measure.find_contours(true_input, lam)[0]
# /** So this might not have very many numbers in it. **/
ic = measure.find_contours(ideal, lam)[0]
# /** CAUTION: this is assuming the contours are going the same
# direction. If not, you'll need to make it so. **/
#nc = np.zeros(ic.shape) # /** don't need ones() **/
# /** We just want to find points on _tc_ to match _ic_. That's
# interpolation _within_ a curve. **/
#pdb.set_trace()
nc_by_t = equispace(tc,ic.shape[0])
ic_by_t = equispace(ic,ic.shape[0])
### /** Not this **/
## Use the y /** X? **/ coordinates of the ideal contour
#nc[:, 0] = ic[:, 0]
#
## Interpolate true contour onto ideal contour y axis so there are the same number of points
#
## /** Have to sort first - https://docs.scipy.org/doc/numpy/reference/generated/numpy.interp.html#numpy-interp **/
#tc_sorted = tc[tc[:,0].argsort()]
# # /** Thanks to by
# # https://whosebug.com/users/208339/steve-tjoa **/
#
#nc[:, 1] = np.interp(ic[:, 0], tc_sorted[:, 0], tc_sorted[:, 1],
# left=np.nan, right=np.nan)
# # /** nan: If the interpolation is out of range, we're not getting
# # useful data. Therefore, flag it with a nan. **/
# /** Filter out the NaNs **/
# Thanks to by
# https://whosebug.com/users/449449/eumiro
#pdb.set_trace()
#indices = ~np.isnan(nc).any(axis=1)
#nc_nonan = nc[indices]
#ic_nonan = ic[indices]
#
# Add the control points to the appropriate list.
## /** Flattening here since otherwise I wound up with dtype=object
## in the numpy arrays. **/
#true_control_pts.append(nc_nonan.flatten().tolist())
#ideal_control_pts.append(ic_nonan.flatten().tolist())
#OLD_true.append(nc) # /** for debug **/
#OLD_ideal.append(ic)
true_control_pts.append(nc_by_t)
ideal_control_pts.append(ic_by_t)
except (IndexError,AttributeError):
pass
#pdb.set_trace()
# /** Make vectors of all the control points. **/
true_flat = list(flatten_list(true_control_pts))
ideal_flat = list(flatten_list(ideal_control_pts))
true_control_pts = np.array(true_flat)
ideal_control_pts = np.array(ideal_flat)
# Make the vectors 2d
length = len(true_control_pts)/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
#pdb.set_trace()
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]) / 30.0 # /**.astype(np.int32)**/
# /** You don't want int32 images! See
# http://scikit-image.org/docs/dev/user_guide/data_types.html .
# Manually rescale the image to [0.0,1.0] by dividing by 30. **/
# ======================
# /** Pad from 10x20 to 20x30 just for grins **/
#pdb.set_trace()
image = np.concatenate( (np.zeros((20,5)),image,np.zeros((20,5))), 1)
# now it's 20x20
image = np.concatenate( (np.zeros((5,20)),image,np.zeros((5,20))), 0)
# ======================
#Plot it
image_float = img_as_float(image) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title('Input image')
fig.show() # /** I needed this on my test system **/
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# /** since we started with float, and this is float, too, the two are
# comparable. **/
pdb.set_trace()
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none') # /**note: float**/
plt.title('Should be parallel lines')
fig.show()
# /** Show the control points.
# The z=0 plane will be the "true" control points (before), and the
# z=1 plane will be the "ideal" control points (after). **/
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
fig.show()
for rowidx in range(length):
ax.plot([true_control_pts[rowidx,0], ideal_control_pts[rowidx,0]],
[true_control_pts[rowidx,1], ideal_control_pts[rowidx,1]],
[0,1])
input() # /** because I was running from the command line **/
我有一张图像，我正尝试使用 skimage.PiecewiseAffineTransform 和 skimage.warp 进行变形。我有一组控制点 (true) 映射到一组新的控制点 (ideal)，但扭曲没有返回我期望的结果。
在这个例子中，我有一个简单的波长梯度，我试图 'straighten out' 到列中。 （您可能会问我为什么要查找轮廓和插值，但那是因为我实际上是将这段代码应用到更复杂的用例中。我只是想重现这个导致相同奇怪输出的简单示例的所有代码。）
为什么我的输出图像只是将输入图像扭曲成正方形并插入？我正在使用 Python 2.7.12 和 matplotlib 1.5.1。这是代码。
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform
true = np.array([range(i,i+10) for i in range(20)])
ideal = np.array([range(10)]*20)
# Find contours of ideal and true images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
for lam in ideal[0]:
try:
# Get the isowavelength contour in the true and ideal images
tc = measure.find_contours(true, lam)[0]
ic = measure.find_contours(ideal, lam)[0]
nc = np.ones(ic.shape)
# Use the y coordinates of the ideal contour
nc[:, 0] = ic[:, 0]
# Interpolate true contour onto ideal contour y axis so there are the same number of points
nc[:, 1] = np.interp(ic[:, 0], tc[:, 0], tc[:, 1])
# Add the control points to the appropriate list
true_control_pts.append(nc.tolist())
ideal_control_pts.append(ic.tolist())
except (IndexError,AttributeError):
pass
true_control_pts = np.array(true_control_pts)
ideal_control_pts = np.array(ideal_control_pts)
length = len(true_control_pts.flatten())/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]).astype(np.int32)
plt.figure()
plt.imshow(image, origin='lower', interpolation='none')
plt.title('Input image')
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none')
plt.title('Should be parallel lines')
输出如下：
如有任何帮助，我们将不胜感激！
我会把我所有的都给你。我不认为我能够用给出的测试数据确定这一点。在这么小的图像中将 45 度角映射到一条直线意味着需要发生很多运动，而不是大量的数据作为基础。我确实发现了一些在下面的代码中修复的特定错误，标记为 /* */（因为这不是您通常在 Python 文件中看到的东西，该标记应该突出 :)）。
请在您的真实数据上试用此代码，如果有效请告诉我！对于这个输入数据集，有一些非零输出。
主要问题是：
interp数据需要排序- 许多控制点都是无意义的，因为此数据集中没有足够的数据（例如，将大跨度的列映射到单个点）
- 图像中的值范围已关闭（因此您的原始输出实际上几乎都是 0 - 色块的值在 1e-9 左右）。
我添加的最重要的东西是 3D 图，显示了 "true" 控制点如何映射到 "ideal" 控制点。这为您提供了一个调试工具来显示您的控制点映射是否符合您的预期。正是这个情节让我遇到了 interp 问题。
顺便说一句，请使用 "before" 和 "after" 之类的名称，而不是 ideal 和 true :)。试图记住哪个是哪个至少让我绊倒了一次。
第一次尝试
import pdb
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform, img_as_float
from mpl_toolkits.mplot3d import Axes3D # /**/
#/**/
# From by
# https://whosebug.com/users/1355221/dansalmo
def flatten_list(L):
for item in L:
try:
for i in flatten_list(item): yield i
except TypeError:
yield item
#end flatten_list
true_input = np.array([range(i,i+10) for i in range(20)]) # /** != True **/
ideal = np.array([range(10)]*20)
#pdb.set_trace()
# Find contours of ideal and true_input images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
OLD_true=[] # /**/ for debugging
OLD_ideal=[]
for lam in [x+0.5 for x in ideal[0]]: # I tried the 0.5 offset just to see,
try: # but it didn't make much difference
# Get the isowavelength contour in the true_input and ideal images
tc = measure.find_contours(true_input, lam)[0]
ic = measure.find_contours(ideal, lam)[0]
nc = np.zeros(ic.shape) # /** don't need ones() **/
# Use the y /** X? **/ coordinates of the ideal contour
nc[:, 0] = ic[:, 0]
# Interpolate true contour onto ideal contour y axis so there are the same number of points
# /** Have to sort first - https://docs.scipy.org/doc/numpy/reference/generated/numpy.interp.html#numpy-interp **/
tc_sorted = tc[tc[:,0].argsort()]
# /** Thanks to by
# https://whosebug.com/users/208339/steve-tjoa **/
nc[:, 1] = np.interp(ic[:, 0], tc_sorted[:, 0], tc_sorted[:, 1],
left=np.nan, right=np.nan)
# /** nan: If the interpolation is out of range, we're not getting
# useful data. Therefore, flag it with a nan. **/
# /** Filter out the NaNs **/
# Thanks to by
# https://whosebug.com/users/449449/eumiro
#pdb.set_trace()
indices = ~np.isnan(nc).any(axis=1)
nc_nonan = nc[indices]
ic_nonan = ic[indices]
# Add the control points to the appropriate list.
# /** Flattening here since otherwise I wound up with dtype=object
# in the numpy arrays. **/
true_control_pts.append(nc_nonan.flatten().tolist())
ideal_control_pts.append(ic_nonan.flatten().tolist())
OLD_true.append(nc) # /** for debug **/
OLD_ideal.append(ic)
except (IndexError,AttributeError):
pass
#pdb.set_trace()
# /** Make vectors of all the control points. **/
true_flat = list(flatten_list(true_control_pts))
ideal_flat = list(flatten_list(ideal_control_pts))
true_control_pts = np.array(true_flat)
ideal_control_pts = np.array(ideal_flat)
# Make the vectors 2d
length = len(true_control_pts)/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
#pdb.set_trace()
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]) / 30.0 # /**.astype(np.int32)**/
# /** You don't want int32 images! See
# http://scikit-image.org/docs/dev/user_guide/data_types.html .
# Manually rescale the image to [0.0,1.0] by dividing by 30. **/
image_float = img_as_float(image) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title('Input image')
fig.show() # /** I needed this on my test system **/
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# /** since we started with float, and this is float, too, the two are
# comparable. **/
pdb.set_trace()
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none') # /**note: float**/
plt.title('Should be parallel lines')
fig.show()
# /** Show the control points.
# The z=0 plane will be the "true" control points (before), and the
# z=1 plane will be the "ideal" control points (after). **/
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
fig.show()
for rowidx in range(length):
ax.plot([true_control_pts[rowidx,0], ideal_control_pts[rowidx,0]],
[true_control_pts[rowidx,1], ideal_control_pts[rowidx,1]],
[0,1])
input() # /** because I was running from the command line **/
第二次尝试
越来越近：
下面是看起来更有希望的控制点映射视图：
您可以看到它试图稍微旋转图像，这正是我希望从此数据集中看到的结果。
代码
import pdb
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform, img_as_float
from mpl_toolkits.mplot3d import Axes3D # /**/
#/**/
# From by
# https://whosebug.com/users/1355221/dansalmo
def flatten_list(L):
for item in L:
try:
for i in flatten_list(item): yield i
except TypeError:
yield item
#end flatten_list
#/**/
# Modified from by
# https://whosebug.com/users/2588210/christian-k
def equispace(data, npts):
x,y = data.T
xd =np.diff(x)
yd = np.diff(y)
dist = np.sqrt(xd**2+yd**2)
u = np.cumsum(dist)
u = np.hstack([[0],u])
t = np.linspace(0,u.max(),npts)
xn = np.interp(t, u, x)
yn = np.interp(t, u, y)
return np.column_stack((xn, yn))
true_input = np.array([range(i,i+10) for i in range(20)]) # /** != True **/
ideal = np.array([range(10)]*20)
#pdb.set_trace()
# Find contours of ideal and true_input images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
OLD_true=[] # /**/ for debugging
OLD_ideal=[]
for lam in [x+0.5 for x in ideal[0]]: # I tried the 0.5 offset just to see,
try: # but it didn't make much difference
# Get the isowavelength contour in the true_input and ideal images
tc = measure.find_contours(true_input, lam)[0]
# /** So this might not have very many numbers in it. **/
ic = measure.find_contours(ideal, lam)[0]
# /** CAUTION: this is assuming the contours are going the same
# direction. If not, you'll need to make it so. **/
#nc = np.zeros(ic.shape) # /** don't need ones() **/
# /** We just want to find points on _tc_ to match _ic_. That's
# interpolation _within_ a curve. **/
#pdb.set_trace()
nc_by_t = equispace(tc,ic.shape[0])
ic_by_t = equispace(ic,ic.shape[0])
### /** Not this **/
## Use the y /** X? **/ coordinates of the ideal contour
#nc[:, 0] = ic[:, 0]
#
## Interpolate true contour onto ideal contour y axis so there are the same number of points
#
## /** Have to sort first - https://docs.scipy.org/doc/numpy/reference/generated/numpy.interp.html#numpy-interp **/
#tc_sorted = tc[tc[:,0].argsort()]
# # /** Thanks to by
# # https://whosebug.com/users/208339/steve-tjoa **/
#
#nc[:, 1] = np.interp(ic[:, 0], tc_sorted[:, 0], tc_sorted[:, 1],
# left=np.nan, right=np.nan)
# # /** nan: If the interpolation is out of range, we're not getting
# # useful data. Therefore, flag it with a nan. **/
# /** Filter out the NaNs **/
# Thanks to by
# https://whosebug.com/users/449449/eumiro
#pdb.set_trace()
#indices = ~np.isnan(nc).any(axis=1)
#nc_nonan = nc[indices]
#ic_nonan = ic[indices]
#
# Add the control points to the appropriate list.
## /** Flattening here since otherwise I wound up with dtype=object
## in the numpy arrays. **/
#true_control_pts.append(nc_nonan.flatten().tolist())
#ideal_control_pts.append(ic_nonan.flatten().tolist())
#OLD_true.append(nc) # /** for debug **/
#OLD_ideal.append(ic)
true_control_pts.append(nc_by_t)
ideal_control_pts.append(ic_by_t)
except (IndexError,AttributeError):
pass
pdb.set_trace()
# /** Make vectors of all the control points. **/
true_flat = list(flatten_list(true_control_pts))
ideal_flat = list(flatten_list(ideal_control_pts))
true_control_pts = np.array(true_flat)
ideal_control_pts = np.array(ideal_flat)
# Make the vectors 2d
length = len(true_control_pts)/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
#pdb.set_trace()
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]) / 30.0 # /**.astype(np.int32)**/
# /** You don't want int32 images! See
# http://scikit-image.org/docs/dev/user_guide/data_types.html .
# Manually rescale the image to [0.0,1.0] by dividing by 30. **/
image_float = img_as_float(image) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title('Input image')
fig.show() # /** I needed this on my test system **/
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# /** since we started with float, and this is float, too, the two are
# comparable. **/
pdb.set_trace()
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none') # /**note: float**/
plt.title('Should be parallel lines')
fig.show()
# /** Show the control points.
# The z=0 plane will be the "true" control points (before), and the
# z=1 plane will be the "ideal" control points (after). **/
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
fig.show()
for rowidx in range(length):
ax.plot([true_control_pts[rowidx,0], ideal_control_pts[rowidx,0]],
[true_control_pts[rowidx,1], ideal_control_pts[rowidx,1]],
[0,1])
input() # /** because I was running from the command line **/
说明
以不同的方式考虑插值：X 和 Y 坐标实际上并不是您想要映射的，对吗？我想你想要映射的是沿着轮廓的距离，这样对角线轮廓就会被拉伸成垂直轮廓。这就是这些行的作用：
nc_by_t = equispace(tc,ic.shape[0])
ic_by_t = equispace(ic,ic.shape[0])
我们沿每个等高线创建 ic.shape[0] 个等距点，然后将这些点相互映射。 equispace 修改自 this answer。因此，这会拉伸较短的轮廓以适应较长的轮廓，无论哪种方式，都由 ic 中的轮廓长度定义的点数。实际上，您可以使用此技术使用任意数量的点；我刚刚测试了 100 点，得到了基本相同的结果。
再次在您的实际数据上进行测试。想想你的插值参考究竟是什么。它真的是 X 或 Y 坐标吗？是沿着轮廓的距离吗？等高线的百分比（如上）?
好的，多一个主意
对于这个特定的测试用例，更多的数据会有帮助吗？是的！
我使用较大的图像来确定控制点，然后在该区域的中心映射较小的图像。
代码
（此时一团糟 - 见 =========== 标记）
import pdb
import matplotlib.pyplot as plt
import numpy as np
from skimage import measure, transform, img_as_float
from mpl_toolkits.mplot3d import Axes3D # /**/
#/**/
def test_figure(data,title):
image_float = img_as_float(data) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title(title)
fig.show()
#/**/
# From by
# https://whosebug.com/users/1355221/dansalmo
def flatten_list(L):
for item in L:
try:
for i in flatten_list(item): yield i
except TypeError:
yield item
#end flatten_list
#/**/
# Modified from by
# https://whosebug.com/users/2588210/christian-k
def equispace(data, npts):
x,y = data.T
xd =np.diff(x)
yd = np.diff(y)
dist = np.sqrt(xd**2+yd**2)
u = np.cumsum(dist)
u = np.hstack([[0],u])
t = np.linspace(0,u.max(),npts)
xn = np.interp(t, u, x)
yn = np.interp(t, u, y)
return np.column_stack((xn, yn))
# ====================== BIGGER
true_input = np.array([range(i,i+20) for i in range(30)]) # /** != True **/
ideal = np.array([range(20)]*30)
# ======================
test_figure(true_input / 50.0, 'true_input')
test_figure(ideal / 20.0, 'ideal')
#pdb.set_trace()
# Find contours of ideal and true_input images and create list of control points for warp
true_control_pts = []
ideal_control_pts = []
OLD_true=[] # /**/ for debugging
OLD_ideal=[]
for lam in [x+0.5 for x in ideal[0]]: # I tried the 0.5 offset just to see,
try: # but it didn't make much difference
# Get the isowavelength contour in the true_input and ideal images
tc = measure.find_contours(true_input, lam)[0]
# /** So this might not have very many numbers in it. **/
ic = measure.find_contours(ideal, lam)[0]
# /** CAUTION: this is assuming the contours are going the same
# direction. If not, you'll need to make it so. **/
#nc = np.zeros(ic.shape) # /** don't need ones() **/
# /** We just want to find points on _tc_ to match _ic_. That's
# interpolation _within_ a curve. **/
#pdb.set_trace()
nc_by_t = equispace(tc,ic.shape[0])
ic_by_t = equispace(ic,ic.shape[0])
### /** Not this **/
## Use the y /** X? **/ coordinates of the ideal contour
#nc[:, 0] = ic[:, 0]
#
## Interpolate true contour onto ideal contour y axis so there are the same number of points
#
## /** Have to sort first - https://docs.scipy.org/doc/numpy/reference/generated/numpy.interp.html#numpy-interp **/
#tc_sorted = tc[tc[:,0].argsort()]
# # /** Thanks to by
# # https://whosebug.com/users/208339/steve-tjoa **/
#
#nc[:, 1] = np.interp(ic[:, 0], tc_sorted[:, 0], tc_sorted[:, 1],
# left=np.nan, right=np.nan)
# # /** nan: If the interpolation is out of range, we're not getting
# # useful data. Therefore, flag it with a nan. **/
# /** Filter out the NaNs **/
# Thanks to by
# https://whosebug.com/users/449449/eumiro
#pdb.set_trace()
#indices = ~np.isnan(nc).any(axis=1)
#nc_nonan = nc[indices]
#ic_nonan = ic[indices]
#
# Add the control points to the appropriate list.
## /** Flattening here since otherwise I wound up with dtype=object
## in the numpy arrays. **/
#true_control_pts.append(nc_nonan.flatten().tolist())
#ideal_control_pts.append(ic_nonan.flatten().tolist())
#OLD_true.append(nc) # /** for debug **/
#OLD_ideal.append(ic)
true_control_pts.append(nc_by_t)
ideal_control_pts.append(ic_by_t)
except (IndexError,AttributeError):
pass
#pdb.set_trace()
# /** Make vectors of all the control points. **/
true_flat = list(flatten_list(true_control_pts))
ideal_flat = list(flatten_list(ideal_control_pts))
true_control_pts = np.array(true_flat)
ideal_control_pts = np.array(ideal_flat)
# Make the vectors 2d
length = len(true_control_pts)/2
true_control_pts = true_control_pts.reshape(length,2)
ideal_control_pts = ideal_control_pts.reshape(length,2)
#pdb.set_trace()
# Plot the original image
image = np.array([range(i,i+10) for i in range(20)]) / 30.0 # /**.astype(np.int32)**/
# /** You don't want int32 images! See
# http://scikit-image.org/docs/dev/user_guide/data_types.html .
# Manually rescale the image to [0.0,1.0] by dividing by 30. **/
# ======================
# /** Pad from 10x20 to 20x30 just for grins **/
#pdb.set_trace()
image = np.concatenate( (np.zeros((20,5)),image,np.zeros((20,5))), 1)
# now it's 20x20
image = np.concatenate( (np.zeros((5,20)),image,np.zeros((5,20))), 0)
# ======================
#Plot it
image_float = img_as_float(image) #/** make sure skimage is happy */
fig = plt.figure()
plt.imshow(image_float, origin='lower', interpolation='none')
plt.title('Input image')
fig.show() # /** I needed this on my test system **/
# Warp the actual image given the transformation between the true and ideal wavelength maps
tform = transform.PiecewiseAffineTransform()
tform.estimate(true_control_pts, ideal_control_pts)
out = transform.warp(image, tform)
# /** since we started with float, and this is float, too, the two are
# comparable. **/
pdb.set_trace()
# Plot the warped image!
fig, ax = plt.subplots()
ax.imshow(out, origin='lower', interpolation='none') # /**note: float**/
plt.title('Should be parallel lines')
fig.show()
# /** Show the control points.
# The z=0 plane will be the "true" control points (before), and the
# z=1 plane will be the "ideal" control points (after). **/
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
fig.show()
for rowidx in range(length):
ax.plot([true_control_pts[rowidx,0], ideal_control_pts[rowidx,0]],
[true_control_pts[rowidx,1], ideal_control_pts[rowidx,1]],
[0,1])
input() # /** because I was running from the command line **/