Elipticalgauss2D_fit_act

%%%%%%%%
  To fit a 2D gaussian to the raster maps
  Use starting best guess values of 
  
  Expects input arrays from e.g. meshgrid XI, YI and ZI
   Code then makes into suitable vectors
  amp = amplitude
  sigmax = x-width
  sigmay = y-width
  x0 = offset in x
  y0 = offset in y0

 If no start points given, it will estimate from the data

 ACT 11/5/2012

%%%%%%%%%%%%%

0001 function [sf gof t] = Ellipticalgauss2D_fit_act(XI,YI,ZI,a1,w1,w2,x0,y0,c1,c2)
0002 
0003 
0004 %%%%%%%%%
0005 %  To fit a 2D gaussian to the raster maps
0006 %  Use starting best guess values of
0007 %
0008 %  Expects input arrays from e.g. meshgrid XI, YI and ZI
0009 %   Code then makes into suitable vectors
0010 %  amp = amplitude
0011 %  sigmax = x-width
0012 %  sigmay = y-width
0013 %  x0 = offset in x
0014 %  y0 = offset in y0
0015 %
0016 % If no start points given, it will estimate from the data
0017 %
0018 % ACT 11/5/2012
0019 %
0020 %%%%%%%%%%%%%%
0021 
0022 % Creat suitable vectors
0023 X = XI(:);
0024 Y = YI(:);
0025 Z = ZI(:);
0026 
0027 % Define a general Gaussian which has angle offset to allow for
0028 % ellipticity:
0029  a = (cos(t1)^2/2/w1^2 + sin(t1)^2/2/w2^2);
0030  b = (-sin(2*t1)/4/w1^2 + sin(2*t1)/4/w2^2) ;
0031  c = (sin(t1)^2/2/w1^2 + cos(t1)^2/2/w2^2);
0032  f = fittype('rat33')
0033     function
0034  ff(a,b,c,a1,x0,y0,t1,w1,w2) = a1*exp( - (a*(x-x0).^2 + 2*b*(x-x0).*(y-y0) + c*(y-y0).^2;
0035  
0036 gauss2=fittype( @(a, b, c, x, y)  a1*exp( - (a*(x-x0).^2 + 2*b*(x-x0).*(y-y0) + c*(y-y0).^2)),'independent', {'x', 'y'},'dependent', 'z');
0037 
0038 %gauss2 = fittype( 'a1*exp(-(((x-c1)*cosd(t1)+(y-c2)*sind(t1))/w1)^2-((-(x-c1)*sind(t1)+(y-c2)*cosd(t1))/w2)^2)', 'independent', {'x', 'y'}, 'dependent', 'z' );
0039 
0040 %gauss2 = fittype( 'a1*exp(-(x-x0).^2/(2*sigmax^2)-(y-y0).^2/(2*sigmay^2))','independent', {'x', 'y'},'dependent', 'z' );
0041 
0042 % Make sure we have a starting point
0043 
0044 if(nargin<8)
0045     
0046     a1 = max(Z(:)); % height, determine from image. may want to subtract background
0047     w1 = 0.4;%guess width
0048     w2 = 0.4;%guess width
0049     %sigmax = 0.4; % guess width
0050     %sigmay = 0.4; % guess width
0051     c1 = mean(X); % guess position (center seems a good place to start)
0052     c2 = mean(Y);
0053     t1=0.0; % Guess starting position
0054 end
0055 
0056 % Set some fit options
0057 options = fitoptions('Method','NonlinearLeastSquares');
0058 options.StartPoint = [a1,w1,w2,t1,c1,c2];
0059 
0060 %Check for Nans in the data
0061 aa = isnan(Z);
0062 options.Exclude = aa;
0063 
0064 % Compute fit and plot figure
0065 
0066 [sf gof t] = fit([X,Y],double(Z),gauss2,options);
0067 
0068 
0069 figure(20); clf; plot(sf,[X,Y],Z);
0070 
0071 function e_g
0072 %sf.x0 and sf.y0 is the center of gaussian.
0073 %sf.sigmax etc will get you the other parameters.