1 | # Fit d'une gaussienne 2D - source : http://www.scipy.org/Cookbook/FittingData |
---|
2 | |
---|
3 | from numpy import * |
---|
4 | from scipy import optimize |
---|
5 | |
---|
6 | def gaussian(height, center_x, center_y, width_x, width_y): |
---|
7 | """Returns a gaussian function with the given parameters""" |
---|
8 | width_x = float(width_x) |
---|
9 | width_y = float(width_y) |
---|
10 | return lambda x,y: height*exp( |
---|
11 | -(((center_x-x)/width_x)**2+((center_y-y)/width_y)**2)/2) |
---|
12 | |
---|
13 | def moments(data): |
---|
14 | """Returns (height, x, y, width_x, width_y) |
---|
15 | the gaussian parameters of a 2D distribution by calculating its |
---|
16 | moments """ |
---|
17 | total = data.sum() |
---|
18 | X, Y = indices(data.shape) |
---|
19 | x = (X*data).sum()/total |
---|
20 | y = (Y*data).sum()/total |
---|
21 | col = data[:, int(y)] |
---|
22 | width_x = sqrt(abs((arange(col.size)-y)**2*col).sum()/col.sum()) |
---|
23 | row = data[int(x), :] |
---|
24 | width_y = sqrt(abs((arange(row.size)-x)**2*row).sum()/row.sum()) |
---|
25 | height = data.max() |
---|
26 | return height, x, y, width_x, width_y |
---|
27 | |
---|
28 | def fitgaussian(data): |
---|
29 | """Returns (height, x, y, width_x, width_y) |
---|
30 | the gaussian parameters of a 2D distribution found by a fit""" |
---|
31 | params = moments(data) |
---|
32 | errorfunction = lambda p: ravel(gaussian(*p)(*indices(data.shape)) - |
---|
33 | data) |
---|
34 | p, success = optimize.leastsq(errorfunction, params) |
---|
35 | return p |
---|