gaussian fit to 2D pcolor map

I have a 2D map using pcolor function. Now I want to fit gaussian to this 2D map P1. How can I do that?

%[xx,yy] = meshgrid(x,y);
xx = [    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323
    3.9202    7.6226   11.3250   15.0274   18.7298   22.4323];% 8x6 matrix
yy = [    5.0820    5.0820    5.0820    5.0820    5.0820    5.0820
    6.6220    6.6220    6.6220    6.6220    6.6220    6.6220
    8.1620    8.1620    8.1620    8.1620    8.1620    8.1620
    9.7020    9.7020    9.7020    9.7020    9.7020    9.7020
   11.2420   11.2420   11.2420   11.2420   11.2420   11.2420
   12.7820   12.7820   12.7820   12.7820   12.7820   12.7820
   14.3220   14.3220   14.3220   14.3220   14.3220   14.3220
   15.8620   15.8620   15.8620   15.8620   15.8620   15.8620];% 8x6 matrix
C = [   0    30   133   199   143    27     0   NaN
   123   614   765   832   810   590   100   NaN
   388   787   897   891   903   857   442   NaN
   125   570   744   737   782   659   176   NaN
     0    53   180   270   210    63     0   NaN
   NaN   NaN   NaN   NaN   NaN   NaN   NaN   NaN] ;% C is a 6x8 matrix of data
P1 = pcolor(xx',yy',C);

Answer 1

[Matt J]

You could try gaussfitn,

https://www.mathworks.com/matlabcentral/fileexchange/69116-gaussfitn

idx=~isnan(C);
xy=[xx(idx),yy(idx)];
params=gaussfitn(xy,C(idx));

[D,A,mu,sigma]=deal(params{:})

D =

190.8900

A =

835.3517

mu =

11.3936

11.1766

sigma =

151.2489 -130.1645

-130.1645 114.3164

[Ham Man]

after running the function I see this expression:

Local minimum possible.

lsqcurvefit stopped because the final change in the sum of squares relative to 
its initial value is less than the value of the function tolerance.

<stopping criteria details>

does it mean the function has worked fine? If so, how can I plot the gussian result in 2D or 3D?

I appreciate your time Matt!

[Matt J]

does it mean the function has worked fine?

The message means lsqcurvefit thinks it succeeded. No optimization code can ever be totally certain of its own success. You can call gaussfitn with additional outputs to get the same kind of additional diagnostic information that lsqcurvefit gives you,

[params,resnorm, residual,exitflag,output] = gaussfitn()

If so, how can I plot the gussian result in 2D or 3D?

You can use surf or similar surface plotting tools.

[Ham Man]

I tried to use surf:

Since the gaussfitn() outputs are not matrix, how can I use the outputs of this function as an input for surf??

[Matt J]

You don't. You use the parameters of the fit to calculate the Gaussian surface on a meshgrid of points. You would give that to surf().

[Ham Man]

Sorry, no luck to calculate the Gaussian surface!

[Matt J]

I don't know what to advise because I don't know what you've tried with surf(). The steps should have been virtually identical to this documentation example, except that instead of Z = sin(X) + cos(Y) you would use the formula for the Gaussian surface.

[Ham Man]

Thank you so much Matt for replying. what I'm doing is:

xx =[

   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000
   -1.2000   -1.0000   -0.6000   -0.2000    0.2000    0.6000    1.0000 ];

yy =[

   -1.1070   -1.1070   -1.1070   -1.1070   -1.1070   -1.1070   -1.1070
   -0.8556   -0.8556   -0.8556   -0.8556   -0.8556   -0.8556   -0.8556
   -0.6111   -0.6111   -0.6111   -0.6111   -0.6111   -0.6111   -0.6111
   -0.3667   -0.3667   -0.3667   -0.3667   -0.3667   -0.3667   -0.3667
   -0.1222   -0.1222   -0.1222   -0.1222   -0.1222   -0.1222   -0.1222
    0.1222    0.1222    0.1222    0.1222    0.1222    0.1222    0.1222
    0.3667    0.3667    0.3667    0.3667    0.3667    0.3667    0.3667
    0.6111    0.6111    0.6111    0.6111    0.6111    0.6111    0.6111
    0.8556    0.8556    0.8556    0.8556    0.8556    0.8556    0.8556 ];
    
C =  [ NaN   NaN   NaN   NaN   NaN   NaN   NaN   NaN   NaN
       NaN     0    30   133   199   143    27     0   NaN
       NaN   123   614   765   832   810   590   100   NaN
       NaN   388   787   897   891   903   857   442   NaN
       NaN   125   570   744   737   782   659   176   NaN
       NaN     0    53   180   270   210    63     0   NaN
       NaN   NaN   NaN   NaN   NaN   NaN   NaN   NaN   NaN ];

idx = ~isnan(C);
xy = [xx(idx),yy(idx)];
params = gaussfitn(xy,C(idx));
[D,A,mu,sigma] = deal(params{:});
% asuume mu=0,bucause mu from param is 2x1 matrix and the dim is
% not agree for z calculation:
z = D + A*exp( -0.5 * (xy-0) * inv(sigma) .*(xy-0) );
surf(z);

The first problem is: z is a function of x only,not x and y. I just put xy in z formula!!

The second problem is why I get two values for z? I know xy is 35x2 and I should expect z of 35x2 matrix. Does is make sence to have two values for z?

[Matt J]

Perhaps as follows.

idx = ~isnan(C);
xy = [xx(idx),yy(idx)];
params = gaussfitn(xy,C(idx));
[D,A,mu,sigma] = deal(params{:});
% asuume mu=0,bucause mu from param is 2x1 matrix and the dim is
% not agree for z calculation:

[X,Y]=meshgrid(mu(1)+linspace(-3,3)*sqrt(sigma(1)),...
               mu(2)+linspace(-3,3)*sqrt(sigma(end)),...);
               
dXY=[X(:),Y(:)]-mu(:)';

Z = D + A*exp( -0.5 *sum((sigma\dXY).*dXY.',1) );

Z=reshape(Z,size(X));

surf(X,Y,Z); hold on;
scatter(xy(:,1),xy(:,2)); hold off