2012-01-25 31 views
5

Tôi đang cố gắng vừa với một Gaussian 2D thành một hình ảnh. Tiếng ồn rất thấp, do đó, nỗ lực của tôi là xoay hình ảnh sao cho hai trục chính không đồng thay đổi, tìm ra mức tối đa và chỉ tính toán độ lệch chuẩn trong cả hai chiều. Vũ khí của sự lựa chọn là python.Tính toán các thành phần riêng biệt của hình ảnh trong python

2d more-or-less gaussian distribution

Tuy nhiên tôi đã bị mắc kẹt tại tìm vector riêng của hình ảnh - numpy.linalg.py giả định các điểm dữ liệu rời rạc. Tôi nghĩ về việc chụp ảnh này là một phân bố xác suất, lấy mẫu vài nghìn điểm và sau đó tính toán các đặc tính riêng từ phân bố đó, nhưng tôi chắc chắn phải có cách để tìm ra các đặc tính riêng (ví dụ, bán chính và bán các trục nhỏ của hình elip gaussian) trực tiếp từ hình ảnh đó. Bất kỳ ý tưởng?

Thanks a lot :)

Trả lời

17

Chỉ cần lưu ý nhanh, có một số công cụ để phù hợp với gaussian thành hình ảnh. Điều duy nhất tôi có thể nghĩ ra khỏi đỉnh đầu là scikits.learn, không hoàn toàn định hướng hình ảnh, nhưng tôi biết có những thứ khác.

Để tính toán số liệu riêng của ma trận hiệp phương sai chính xác như bạn nghĩ trong đầu là rất tốn kém về mặt tính toán. Bạn phải liên kết từng pixel (hoặc một mẫu ngẫu nhiên lớn) với một điểm x, y.

Về cơ bản, bạn làm như sau:

import numpy as np 
    # grid is your image data, here... 
    grid = np.random.random((10,10)) 

    nrows, ncols = grid.shape 
    i,j = np.mgrid[:nrows, :ncols] 
    coords = np.vstack((i.reshape(-1), j.reshape(-1), grid.reshape(-1))).T 
    cov = np.cov(coords) 
    eigvals, eigvecs = np.linalg.eigh(cov) 

Bạn thay vì có thể tận dụng một thực tế rằng đó là một hình ảnh thường xuyên lấy mẫu và tính toán đó là những khoảnh khắc (hoặc "trục intertial") để thay thế. Điều này sẽ nhanh hơn đáng kể cho các hình ảnh lớn.

Như một ví dụ nhanh, (Tôi đang sử dụng một phần của một trong previous answers của tôi, trong trường hợp bạn tìm thấy nó hữu ích ...)

import numpy as np 
import matplotlib.pyplot as plt 

def main(): 
    data = generate_data() 
    xbar, ybar, cov = intertial_axis(data) 

    fig, ax = plt.subplots() 
    ax.imshow(data) 
    plot_bars(xbar, ybar, cov, ax) 
    plt.show() 

def generate_data(): 
    data = np.zeros((200, 200), dtype=np.float) 
    cov = np.array([[200, 100], [100, 200]]) 
    ij = np.random.multivariate_normal((100,100), cov, int(1e5)) 
    for i,j in ij: 
     data[int(i), int(j)] += 1 
    return data 

def raw_moment(data, iord, jord): 
    nrows, ncols = data.shape 
    y, x = np.mgrid[:nrows, :ncols] 
    data = data * x**iord * y**jord 
    return data.sum() 

def intertial_axis(data): 
    """Calculate the x-mean, y-mean, and cov matrix of an image.""" 
    data_sum = data.sum() 
    m10 = raw_moment(data, 1, 0) 
    m01 = raw_moment(data, 0, 1) 
    x_bar = m10/data_sum 
    y_bar = m01/data_sum 
    u11 = (raw_moment(data, 1, 1) - x_bar * m01)/data_sum 
    u20 = (raw_moment(data, 2, 0) - x_bar * m10)/data_sum 
    u02 = (raw_moment(data, 0, 2) - y_bar * m01)/data_sum 
    cov = np.array([[u20, u11], [u11, u02]]) 
    return x_bar, y_bar, cov 

def plot_bars(x_bar, y_bar, cov, ax): 
    """Plot bars with a length of 2 stddev along the principal axes.""" 
    def make_lines(eigvals, eigvecs, mean, i): 
     """Make lines a length of 2 stddev.""" 
     std = np.sqrt(eigvals[i]) 
     vec = 2 * std * eigvecs[:,i]/np.hypot(*eigvecs[:,i]) 
     x, y = np.vstack((mean-vec, mean, mean+vec)).T 
     return x, y 
    mean = np.array([x_bar, y_bar]) 
    eigvals, eigvecs = np.linalg.eigh(cov) 
    ax.plot(*make_lines(eigvals, eigvecs, mean, 0), marker='o', color='white') 
    ax.plot(*make_lines(eigvals, eigvecs, mean, -1), marker='o', color='red') 
    ax.axis('image') 

if __name__ == '__main__': 
    main() 

enter image description here

1

Bạn hãy thử Principal Component Analysis (PCA)? Có lẽ MDP package có thể thực hiện công việc với nỗ lực tối thiểu.

3

Lắp một Gaussian mạnh mẽ có thể khôn lanh. Có một bài viết vui về chủ đề này trong Tạp chí Xử lý tín hiệu IEEE:

Hongwei Guo, "Một thuật toán đơn giản cho một Lắp Chức năng Gaussian" IEEE Magazine Xử lý tín hiệu, tháng 9 năm 2011, pp 134--137

.

tôi cung cấp cho một thực hiện các trường hợp 1D đây:

http://scipy-central.org/item/28/2/fitting-a-gaussian-to-noisy-data-points

(Cuộn xuống để xem phù hợp kết quả)