THESE/codes/snake/basic_code/imgaussian.m

   1 function I=imgaussian(I,sigma,siz)\r
   2 % IMGAUSSIAN filters an 1D, 2D color/greyscale or 3D image with an \r
   3 % Gaussian filter. This function uses for filtering IMFILTER or if \r
   4 % compiled the fast  mex code imgaussian.c . Instead of using a \r
   5 % multidimensional gaussian kernel, it uses the fact that a Gaussian \r
   6 % filter can be separated in 1D gaussian kernels.\r
   7 %\r
   8 % J=IMGAUSSIAN(I,SIGMA,SIZE)\r
   9 %\r
  10 % inputs,\r
  11 %   I: The 1D, 2D greyscale/color, or 3D input image with \r
  12 %           data type Single or Double\r
  13 %   SIGMA: The sigma used for the Gaussian kernel\r
  14 %   SIZE: Kernel size (single value) (default: sigma*6)\r
  15 % \r
  16 % outputs,\r
  17 %   J: The gaussian filtered image\r
  18 %\r
  19 % note, compile the code with: mex imgaussian.c -v\r
  20 %\r
  21 % example,\r
  22 %   I = im2double(imread('peppers.png'));\r
  23 %   figure, imshow(imgaussian(I,10));\r
  24 % \r
  25 % Function is written by D.Kroon University of Twente (September 2009)\r
  26 \r
  27 if(~exist('siz','var')), siz=sigma*6; end\r
  28 \r
  29 if(sigma>0)\r
  30     % Make 1D Gaussian kernel\r
  31     x=-ceil(siz/2):ceil(siz/2);\r
  32     H = exp(-(x.^2/(2*sigma^2)));\r
  33     H = H/sum(H(:));\r
  34 \r
  35     % Filter each dimension with the 1D Gaussian kernels\\r
  36     if(ndims(I)==1)\r
  37         I=imfilter(I,H, 'same' ,'replicate');\r
  38     elseif(ndims(I)==2)\r
  39         Hx=reshape(H,[length(H) 1]);\r
  40         Hy=reshape(H,[1 length(H)]);\r
  41         I=imfilter(imfilter(I,Hx, 'same' ,'replicate'),Hy, 'same' ,'replicate');\r
  42     elseif(ndims(I)==3)\r
  43         if(size(I,3)<4) % Detect if 3D or color image\r
  44             Hx=reshape(H,[length(H) 1]);\r
  45             Hy=reshape(H,[1 length(H)]);\r
  46             for k=1:size(I,3)\r
  47                 I(:,:,k)=imfilter(imfilter(I(:,:,k),Hx, 'same' ,'replicate'),Hy, 'same' ,'replicate');\r
  48             end\r
  49         else\r
  50             Hx=reshape(H,[length(H) 1 1]);\r
  51             Hy=reshape(H,[1 length(H) 1]);\r
  52             Hz=reshape(H,[1 1 length(H)]);\r
  53             I=imfilter(imfilter(imfilter(I,Hx, 'same' ,'replicate'),Hy, 'same' ,'replicate'),Hz, 'same' ,'replicate');\r
  54         end\r
  55     else\r
  56         error('imgaussian:input','unsupported input dimension');\r
  57     end\r
  58 end