%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% To fit a 2-D gaussian
%% m = Image
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[cx,cy,sx,sy,PeakOD] = Gaussian2D(m,tol);

[sizey sizex] = size(m);
[x,y] = MeshGrid(1:sizex,1:sizey);
fit = abs(PeakOD)*(exp(-0.5*(x-cx).^2./(sx^2)-0.5*(y-cy).^2./(sy^2)));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% a function to fit a thermal cloud 2-D
function [cx,cy,sx,sy,PeakOD] = Gaussian2D(m,tol);
%% m = image
%% tol = fitting tolerance

options = optimset('Display','off','TolFun',tol,'LargeScale','off');

[sizey sizex] = size(m);
[cx,cy,sx,sy] = centerofmass(m);
pOD = max(max(m));

mx = m(round(cy),;
x1D = 1:sizex;
ip1D = [cx,sx,pOD];
fp1D = fminunc(@fitGaussian1D,ip1D,options,mx,x1D);

cx = fp1D(1);
sx = fp1D(2);
PeakOD = fp1D(3);

my = m(:,round(cx))';
y1D = 1:sizey;
ip1D = [cy,sy,pOD];
fp1D = fminunc(@fitGaussian1D,ip1D,options,my,y1D);

cy = fp1D(1);
sy = fp1D(2);
PeakOD = fp1D(3);
[X,Y] = meshgrid(1:sizex,1:sizey);

initpar = [cx,cy,sx,sy,PeakOD];
fp = fminunc(@fitGaussian2D,initpar,options,m,X,Y);
cx = fp(1);
cy = fp(2);
sx = fp(3);
sy = fp(4);
PeakOD = fp(5);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PURPOSE: find c of m of distribution
function [cx,cy,sx,sy] = centerofmass(m);

[sizey sizex] = size(m);
vx = sum(m);
vy = sum(m');

vx = vx.*(vx>0);
vy = vy.*(vy>0);

x = [1:sizex];
y = [1:sizey];

cx = sum(vx.*x)/sum(vx);
cy = sum(vy.*y)/sum(vy);

sx = sqrt(sum(vx.*(abs(x-cx).^2))/sum(vx));
sy = sqrt(sum(vy.*(abs(y-cy).^2))/sum(vy));


其中cx = sum(vx.*x)/sum(vx);
cy = sum(vy.*y)/sum(vy);

sx = sqrt(sum(vx.*(abs(x-cx).^2))/sum(vx));
sy = sqrt(sum(vy.*(abs(y-cy).^2))/sum(vy));是什么含义 返回小木虫查看更多