function [azim,pol,azimax,nav,x,y]=maxazi(data,varargin) % [azim,pol,azimax,nav,x,y]=MAXAZI(data) % [azim,pol,azimax,nav,x,y]=MAXAZI(data,azim) % [azim,pol,azimax,nav,x,y]=MAXAZI(data,azim,[cenx ceny]) % [azim,pol,azimax,nav,x,y]=MAXAZI(data,azim,[cenx ceny],[minx maxx]) % % OUTPUT: % % azim Azimuthal directions % pol Average coherence along these directions % azimax Direction of the maximum average coherence % nav Number of coherence samples in this average % x,y The indices through the coherence that contains the maximum % % May specify center of plot if it is not in the geometric center of the % indices. For FFT products hence symmetric about the origin! % If a minimum and maximum wavenumber is specified, will compare % the azimuthal curves with LOGLINEARLY spaced samples. In other words, % the mean won't be biased as badly towards the many short-wavelength % coherence values. % % Calculates mean 'pol' and maximum 'azimax' of mean azimuthal line integral % for a data set. Or just calculates the integral over the given azimuth % Mean calculated on the basis of 'nav' samples, from 0 to 180. % Suitable for kx, ky wavelength representation: limit to kx^2+ky^2 rings. % Also returns the x and y indices used to define the maximum azimuth. % % See also ISAV, RADAVG, MINAZI, HALFMAZI % % Note that for images, may need to do flipud before calculating. % % EXAMPLE: % % maxazi('demo'); % % Last modified by fjsimons-at-mit.edu, November 30th, 2005 if ~isstr(data) [ry,rx]=size(data); hrx=ceil(rx/2); hry=ceil(ry/2); if hrx~=hry error('Should input square data set') end resol=180; if nargin==1 | isempty(varargin{1}) azim=linspace(0,180,resol); else azim=varargin{1}; end if nargin>=3 senter=varargin{2}; else senter=[hrx hry]; end if nargin==4 minx=varargin{3}(1); maxx=varargin{3}(2); disp('Loglinear mean comparison!') else disp('Linear mean comparison!') end for index=1:length(azim) % Work with the minimum number of independent measurements % for a given azimuth; assume the center chosen is the symmetry center! % 'senter' and min(size(data)) dictate if indeed senter(1):end is a good % choice (see e.g. halfmazi) [x,y]=linazim(azim(index),senter,hrx,[],min(size(data))); polated=interp2(data,x,y,'linear'); if nargin==4 logpol=interp1(linspace(minx,maxx,length(polated(senter(1):end))),... polated(senter(1):end),logspace(log10(minx),log10(maxx),... length(polated((senter(1):end))))); % Compare % plot(linspace(minx,maxx,length(polated(senter(1):end))),... % polated(senter(1):end),'+-'); hold on % plot(logspace(log10(minx),log10(maxx),length(polated(senter(1):end))),... % logpol,'ro-') pol(index)=nanmean(logpol); nav(index)=length(logpol); else pol(index)=nanmean(polated(senter(1):end)); nav(index)=length(polated(senter(1):end)); end end [i,k]=max(pol); azimax=azim(k); [x,y]=linazim(azimax,senter,hrx,[],min(size(data))); elseif strmatch(data,'demo') [X,Y]=meshgrid(linspace(-10,10,200)); G1=gauss2(X,Y,3,6); [azim,pol]=maxazi(G1); subplot(221); imagesc(G1); axis image subplot(222) ; plot(azim,pol); xlim([0 180]); grid [Xrot,Yrot]=rotdom(linspace(-10,10,200),linspace(-10,10,200),pi/6); G2=gauss2(Xrot,Yrot,3,6); [azim,pol]=maxazi(G2); subplot(223); imagesc(G2); axis image subplot(224) ; plot(azim,pol); xlim([0 180]); grid end