Я хотел бы построить два профиля через точку наивысшей интенсивности в двумерном массиве с нюхами, который является изображением капли (то есть линии, проходящей через большую полуось, и другой линии, проходящей через полуосновную ось). , Капля вращается на угол theta против часовой стрелки от стандартной оси X и является асимметричной.
Это массив 600x600 с максимальной интенсивностью 1 (только на один пиксель), который расположен прямо в центре в точке (300, 300). Угол поворота от оси x (который затем определяет положение большой полуоси при повороте на этот угол) составляет theta = 89.54 degrees. Я не хочу использовать scipy.ndimage.rotate, потому что он использует сплайн-интерполяцию, и я не хочу изменять какие-либо из моих значений пикселей. Но я полагаю, что метод интерполяции ближайшего соседа будет в порядке.
Я попытался сгенерировать линии, соответствующие большой и вспомогательной осям по всему изображению, но результат был совсем не правильным (пик был намного меньше 1), поэтому, возможно, я сделал что-то не так. Код для этого ниже:
import numpy as np
import matplotlib.pyplot as plt
from scipy import ndimage
def profiles_at_angle(image, axis, theta):
theta = np.deg2rad(theta)
if axis == 'major':
x_0, y_0 = 0, 300-300*np.tan(theta)
x_1, y_1 = 599, 300+300*np.tan(theta)
elif axis=='minor':
x_0, y_0 = 300-300*np.tan(theta), 599
x_1, y_1 = 300+300*np.tan(theta), -599
num = 600
x, y = np.linspace(x_0, x_1, num), np.linspace(y_0, y_1, num)
z = ndimage.map_coordinates(image, np.vstack((x,y)))
fig, axes = plt.subplots(nrows=2)
axes[0].imshow(image, cmap='gray')
axes[0].axis('image')
axes[1].plot(z)
plt.xlim(250,350)
plt.show()
profiles_at_angle(image, 'major', theta)
Я сделал что-то явно не так в моем коде выше? Или как еще я могу это сделать? Спасибо.
Редактировать: вот несколько примеров изображений. Извините за плохое качество; мой браузер зависал каждый раз, когда я пытался загрузить их в любое место, поэтому мне приходилось делать фотографии с экрана.
Рисунок 1: Это результат моего кода выше, который явно ошибочен, так как пик должен быть на 1. Я не уверен, что я сделал неправильно, хотя.
Рисунок 2: Я сделал этот график ниже, просто проведя профили по стандартным осям x и y, игнорируя любое вращение (это выглядит хорошо только по совпадению, потому что реальный угол поворота настолько близок к 90 градусам, что я смог просто поменять метки и получить это). Я хочу, чтобы мой результат выглядел примерно так, но с учетом поправочного угла поворота.
Редактировать: Было бы полезно запустить тесты для этого метода, используя данные, очень похожие на мои (это 2D гауссовский с почти такими же параметрами):
image = np.random.random((600,600))
def generate(data_set):
xvec = np.arange(0, np.shape(data_set)[1], 1)
yvec = np.arange(0, np.shape(data_set)[0], 1)
X, Y = np.meshgrid(xvec, yvec)
return X, Y
def gaussian_func(xy, x0, y0, sigma_x, sigma_y, amp, theta, offset):
x, y = xy
a = (np.cos(theta))**2/(2*sigma_x**2) + (np.sin(theta))**2/(2*sigma_y**2)
b = -np.sin(2*theta)/(4*sigma_x**2) + np.sin(2*theta)/(4*sigma_y**2)
c = (np.sin(theta))**2/(2*sigma_x**2) + (np.cos(theta))**2/(2*sigma_y**2)
inner = a * (x-x0)**2
inner += 2*b*(x-x0)*(y-y0)
inner += c * (y-y0)**2
return (offset + amp * np.exp(-inner)).ravel()
xx, yy = generate(image)
image = gaussian_func((xx.ravel(), yy.ravel()), 300, 300, 5, 4, 1, 1.56, 0)
image = np.reshape(image, (600, 600))