--- /dev/null
+## Copyright (C) 2005 Søren Hauberg
+##
+## This program is free software; you can redistribute it and/or modify
+## it under the terms of the GNU General Public License as published by
+## the Free Software Foundation; either version 2 of the License, or
+## (at your option) any later version.
+##
+## This program is distributed in the hope that it will be useful,
+## but WITHOUT ANY WARRANTY; without even the implied warranty of
+## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+## GNU General Public License for more details.
+##
+## You should have received a copy of the GNU General Public License
+## along with this program; If not, see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{filter} = } fspecial(@var{type}, @var{arg1}, @var{arg2})
+## Create spatial filters for image processing.
+##
+## @var{type} determines the shape of the filter and can be
+## @table @t
+## @item "average"
+## Rectangular averaging filter. The optional argument @var{arg1} controls the
+## size of the filter. If @var{arg1} is an integer @var{N}, a @var{N} by @var{N}
+## filter is created. If it is a two-vector with elements @var{N} and @var{M}, the
+## resulting filter will be @var{N} by @var{M}. By default a 3 by 3 filter is
+## created.
+## @item "disk"
+## Circular averaging filter. The optional argument @var{arg1} controls the
+## radius of the filter. If @var{arg1} is an integer @var{N}, a 2 @var{N} + 1
+## filter is created. By default a radius of 5 is used.
+## @item "gaussian"
+## Gaussian filter. The optional argument @var{arg1} controls the size of the
+## filter. If @var{arg1} is an integer @var{N}, a @var{N} by @var{N}
+## filter is created. If it is a two-vector with elements @var{N} and @var{M}, the
+## resulting filter will be @var{N} by @var{M}. By default a 3 by 3 filter is
+## created. The optional argument @var{arg2} sets spread of the filter. By default
+## a spread of @math{0.5} is used.
+## @item "log"
+## Laplacian of Gaussian. The optional argument @var{arg1} controls the size of the
+## filter. If @var{arg1} is an integer @var{N}, a @var{N} by @var{N}
+## filter is created. If it is a two-vector with elements @var{N} and @var{M}, the
+## resulting filter will be @var{N} by @var{M}. By default a 5 by 5 filter is
+## created. The optional argument @var{arg2} sets spread of the filter. By default
+## a spread of @math{0.5} is used.
+## @item "laplacian"
+## 3x3 approximation of the laplacian. The filter is approximated as
+## @example
+## (4/(@var{alpha}+1))*[@var{alpha}/4, (1-@var{alpha})/4, @var{alpha}/4; ...
+## (1-@var{alpha})/4, -1, (1-@var{alpha})/4; ...
+## @var{alpha}/4, (1-@var{alpha})/4, @var{alpha}/4];
+## @end example
+## where @var{alpha} is a number between 0 and 1. This number can be controlled
+## via the optional input argument @var{arg1}. By default it is @math{0.2}.
+## @item "unsharp"
+## Sharpening filter. The following filter is returned
+## @example
+## (1/(@var{alpha}+1))*[-@var{alpha}, @var{alpha}-1, -@var{alpha}; ...
+## @var{alpha}-1, @var{alpha}+5, @var{alpha}-1; ...
+## -@var{alpha}, @var{alpha}-1, -@var{alpha}];
+## @end example
+## where @var{alpha} is a number between 0 and 1. This number can be controlled
+## via the optional input argument @var{arg1}. By default it is @math{0.2}.
+## @item "motion"
+## Moion blur filter of width 1 pixel. The optional input argument @var{arg1}
+## controls the length of the filter, which by default is 9. The argument @var{arg2}
+## controls the angle of the filter, which by default is 0 degrees.
+## @item "sobel"
+## Horizontal Sobel edge filter. The following filter is returned
+## @example
+## [ 1, 2, 1;
+## 0, 0, 0;
+## -1, -2, -1 ]
+## @end example
+## @item "prewitt"
+## Horizontal Prewitt edge filter. The following filter is returned
+## @example
+## [ 1, 1, 1;
+## 0, 0, 0;
+## -1, -1, -1 ]
+## @end example
+## @item "kirsch"
+## Horizontal Kirsch edge filter. The following filter is returned
+## @example
+## [ 3, 3, 3;
+## 3, 0, 3;
+## -5, -5, -5 ]
+## @end example
+## @end table
+## @end deftypefn
+
+## Remarks by Søren Hauberg (jan. 2nd 2007)
+## The motion filter and most of the documentation was taken from Peter Kovesi's
+## GPL'ed implementation of fspecial from
+## http://www.csse.uwa.edu.au/~pk/research/matlabfns/OctaveCode/fspecial.m
+
+function f = fspecial (type, arg1, arg2)
+ if (!ischar (type))
+ error ("fspecial: first argument must be a string");
+ endif
+
+ switch lower(type)
+ case "average"
+ ## Get filtersize
+ if (nargin > 1 && isreal (arg1) && length (arg1 (:)) <= 2)
+ fsize = arg1 (:);
+ else
+ fsize = 3;
+ endif
+ ## Create the filter
+ f = ones (fsize);
+ ## Normalize the filter to integral 1
+ f = f / sum (f (:));
+
+ case "disk"
+ ## Get the radius
+ if (nargin > 1 && isreal (arg1) && length (arg1 (:)) == 1)
+ radius = arg1;
+ else
+ radius = 5;
+ endif
+ ## Create the filter
+ [x, y] = meshgrid (-radius:radius, -radius:radius);
+ r = sqrt (x.^2 + y.^2);
+ f = (r <= radius);
+ ## Normalize the filter to integral 1
+ f = f / sum (f (:));
+
+ case "gaussian"
+ ## Get hsize
+ if (nargin > 1 && isreal (arg1))
+ if (length (arg1 (:)) == 1)
+ hsize = [arg1, arg1];
+ elseif (length (arg1 (:)) == 2)
+ hsize = arg1;
+ else
+ error ("fspecial: second argument must be a scalar or a vector of two scalars");
+ endif
+ else
+ hsize = [3, 3];
+ endif
+ ## Get sigma
+ if (nargin > 2 && isreal (arg2) && length (arg2 (:)) == 1)
+ sigma = arg2;
+ else
+ sigma = 0.5;
+ endif
+ h1 = hsize (1)-1; h2 = hsize (2)-1;
+ [x, y] = meshgrid(0:h2, 0:h1);
+ x = x-h2/2; y = y-h1/2;
+ gauss = exp( -( x.^2 + y.^2 ) / (2*sigma^2) );
+ f = gauss / sum (gauss (:));
+
+ case "laplacian"
+ ## Get alpha
+ if (nargin > 1 && isscalar (arg1))
+ alpha = arg1;
+ if (alpha < 0 || alpha > 1)
+ error ("fspecial: second argument must be between 0 and 1");
+ endif
+ else
+ alpha = 0.2;
+ endif
+ ## Compute filter
+ f = (4/(alpha+1))*[alpha/4, (1-alpha)/4, alpha/4; ...
+ (1-alpha)/4, -1, (1-alpha)/4; ...
+ alpha/4, (1-alpha)/4, alpha/4];
+ case "log"
+ ## Get hsize
+ if (nargin > 1 && isreal (arg1))
+ if (length (arg1 (:)) == 1)
+ hsize = [arg1, arg1];
+ elseif (length (arg1 (:)) == 2)
+ hsize = arg1;
+ else
+ error ("fspecial: second argument must be a scalar or a vector of two scalars");
+ endif
+ else
+ hsize = [5, 5];
+ endif
+ ## Get sigma
+ if (nargin > 2 && isreal (arg2) && length (arg2 (:)) == 1)
+ sigma = arg2;
+ else
+ sigma = 0.5;
+ endif
+ ## Compute the filter
+ h1 = hsize (1)-1; h2 = hsize (2)-1;
+ [x, y] = meshgrid(0:h2, 0:h1);
+ x = x-h2/2; y = y = y-h1/2;
+ gauss = exp( -( x.^2 + y.^2 ) / (2*sigma^2) );
+ f = ( (x.^2 + y.^2 - 2*sigma^2).*gauss )/( 2*pi*sigma^6*sum(gauss(:)) );
+
+ case "motion"
+ ## Taken (with some changes) from Peter Kovesis implementation
+ ## (http://www.csse.uwa.edu.au/~pk/research/matlabfns/OctaveCode/fspecial.m)
+ ## FIXME: The implementation is not quite matlab compatible.
+ if (nargin > 1 && isreal (arg1))
+ len = arg1;
+ else
+ len = 9;
+ endif
+ if (mod (len, 2) == 1)
+ sze = [len, len];
+ else
+ sze = [len+1, len+1];
+ end
+ if (nargin > 2 && isreal (arg2))
+ angle = arg2;
+ else
+ angle = 0;
+ endif
+
+ ## First generate a horizontal line across the middle
+ f = zeros (sze);
+ f (floor (len/2)+1, 1:len) = 1;
+
+ # Then rotate to specified angle
+ f = imrotate (f, angle, "bilinear", "loose");
+ f = f / sum (f (:));
+
+ case "prewitt"
+ ## The filter
+ f = [1, 1, 1; 0, 0, 0; -1, -1, -1];
+
+ case "sobel"
+ ## The filter
+ f = [1, 2, 1; 0, 0, 0; -1, -2, -1];
+
+ case "kirsch"
+ ## The filter
+ f = [3, 3, 3; 3, 0, 3; -5, -5, -5];
+
+ case "unsharp"
+ ## Get alpha
+ if (nargin > 1 && isscalar (arg1))
+ alpha = arg1;
+ if (alpha < 0 || alpha > 1)
+ error ("fspecial: second argument must be between 0 and 1");
+ endif
+ else
+ alpha = 0.2;
+ endif
+ ## Compute filter
+ f = (1/(alpha+1))*[-alpha, alpha-1, -alpha; ...
+ alpha-1, alpha+5, alpha-1; ...
+ -alpha, alpha-1, -alpha];
+
+ otherwise
+ error ("fspecial: filter type '%s' is not supported", type);
+ endswitch
+endfunction
git://git.creatis.insa-lyon.fr / CreaPhase.git / blobdiff