% FILTERGRID Generates grid for constructing frequency domain filters
%
% Usage: [radius, u1, u2] = filtergrid(rows, cols)
% [radius, u1, u2] = filtergrid([rows, cols])
%
% Arguments: rows, cols - Size of image/filter
%
% Returns: radius - Grid of size [rows cols] containing normalised
% radius values from 0 to 0.5. Grid is quadrant
% shifted so that 0 frequency is at radius(1,1)
% u1, u2 - Grids containing normalised frequency values
% ranging from -0.5 to 0.5 in x and y directions
% respectively. u1 and u2 are quadrant shifted.
%
% Used by PHASECONGMONO, PHASECONG3 etc etc
%
% See also: WAVENUMBERGRID
% Copyright (c) 1996-2017 Peter Kovesi
% Centre for Exploration Targeting
% The University of Western Australia
% peter.kovesi at uwa edu au
%
% Permission is hereby granted, free of charge, to any person obtaining a copy
% of this software and associated documentation files (the "Software"), to deal
% in the Software without restriction, subject to the following conditions:
%
% The above copyright notice and this permission notice shall be included in
% all copies or substantial portions of the Software.
%
% The Software is provided "as is", without warranty of any kind.
%
% May 2013
% September 2017 Correction to setting up matrices of frequency values for
% odd sized images.
function [radius, u1, u2] = filtergrid(rows, cols)
% Handle case where rows, cols has been supplied as a 2-vector
if nargin == 1 & length(rows) == 2
tmp = rows;
rows = tmp(1);
cols = tmp(2);
end
% Set up X and Y spatial frequency matrices, u1 and u2 The following code
% adjusts things appropriately for odd and even values of rows and columns
% so that the 0 frequency point is placed appropriately. See
% https://blogs.uoregon.edu/seis/wiki/unpacking-the-matlab-fft/
if mod(cols,2)
u1range = [-(cols-1)/2:(cols-1)/2]/cols;
else
u1range = [-cols/2:(cols/2-1)]/cols;
end
if mod(rows,2)
u2range = [-(rows-1)/2:(rows-1)/2]/rows;
else
u2range = [-rows/2:(rows/2-1)]/rows;
end
[u1,u2] = meshgrid(u1range, u2range);
% Quadrant shift so that filters are constructed with 0 frequency at
% the corners
u1 = ifftshift(u1);
u2 = ifftshift(u2);
% Construct spatial frequency values in terms of normalised radius from
% centre.
radius = sqrt(u1.^2 + u2.^2);