2017-05-31 88 views
0

我无法使用起始点和停止点的坐标切片图像。到目前为止,我有以下代码在两个aribtray点上切片图像

hdulist = fits.open(filename) 
hdr= hdulist[0].header 

import numpy as np 
import scipy 
import matplotlib.pyplot as plt 
from matplotlib.pyplot import figure, show 
from astropy.io import fits 
from scipy import interpolate 

data=hdulist[0].data 

#Make a line with "num" points 
D, B = input('Enter the coordinates of the starting point:').split(',') 
E, C = input("Enter the coordinates of the stopping point: ").split(',') 
x0= float(D) 
x1= float(E) 
y0= float(B) 
y1= float(C) 

x = np.arange(data.shape[1]) 
y = np.arange(data.shape[0]) 
#length = int((np.hypot(x1-x0, y1-y0))) (can be used instead of num_points) 
num_points = 1000 
xvalues = np.linspace(x0, x1, num_points) 
yvalues = np.linspace(y0, y1, num_points) 
f = scipy.interpolate.interp2d(x, y, data) #default is linear 


# Extract the values along the line 
profile = f(xvalues, yvalues) #this gives me a 2D array, I think it needs to be 1D 
#c = profile.flatten() 
print(profile.shape) 

'轮廓'不是线性的,而是立方体。有没有办法让我的轮廓线性化,这样我就可以在起点和终点之间的点上切分图像?我只需要制作'轮廓'1D而不是2D。

我要绘制这样的:

import numpy as np 
from numpy import random 
from matplotlib.pyplot import figure, show 

vels = np.linspace(0, 530, len(profile)) 
fig = figure() 
frame = fig.add_subplot(1,1,1) 
frame.plot(vels, profile) 
frame.set_ylabel('y-axis') 
frame.set_xlabel('x-axis') 
frame.grid(True) 
show() 
print(vels.shape) 
print(profile.shape) 
print(len(profile)) 

我的代码不工作,因为我得到的阴谋没有显示线路的片,但立方体的一个切片。

回答

0

interp2D的文档看来,似乎网格是从插值构建的。就直觉而言,在我看来,你需要该网格的diagonal。建立快速的实验与您的代码的改编:

import numpy as np 
import scipy 
import matplotlib.pyplot as plt 
from matplotlib.pyplot import figure, show 
from scipy import interpolate 

X = np.arange(-5, 5, 0.25) 
Y = np.arange(-5, 5, 0.25) 
X, Y = np.meshgrid(X, Y) 
R = np.sqrt(X**2 + Y**2) 
data = np.sin(R) 

x0 = 2 
x1 = 23 
y0 = 1 
y1 = 36 

rx0 = 4 
rx1 = 2 
ry0 = 7 
ry1 = 32 

x = np.arange(data.shape[1]) 
y = np.arange(data.shape[0]) 
num_points = 1000 
xvalues = np.linspace(x0, x1, num_points) 
yvalues = np.linspace(y0, y1, num_points) 
f = scipy.interpolate.interp2d(x, y, data) #default is linear 


# Extract the values along the line 
profile = f(xvalues, yvalues) 

xvalues2 = np.linspace(rx0, rx1, num_points) 
yvalues2 = np.linspace(ry0, ry1, num_points) 
profile2 = f(xvalues2, yvalues2) 

plt.subplot(121) 
plt.imshow(data.T, origin="lower", interpolation="nearest") 
plt.scatter([x0, x1], [y0, y1]) 
plt.plot([x0, x1], [y0, y1]) 
plt.scatter([rx0, rx1], [ry0, ry1], c="r") 
plt.plot([rx0, rx1], [ry0, ry1], c="r") 
# plt.show() 

diag = np.diag(profile) 
diag2 = np.diag(profile2) 
plt.subplot(122) 
plt.plot(np.arange(diag.shape[0]), diag) 
plt.plot(np.arange(diag2.shape[0]), diag2, c="r") 
plt.show() 

这将返回以下:

Slicing a 2D image

注:我没考虑到的坐标到2D图(这就是为什么这两个线看起来相同的大小)。

+0

非常感谢,我的代码现在按照我想要的方式工作! – Thomas