# HG changeset patch
# User Markus Mützel <markus.muetzel@gmx.de>
# Date 1467835210 -7200
#      Wed Jul 06 22:00:10 2016 +0200
# Node ID f2c92d4a173e349c4201f16cd92e2363c82a9bb1
# Parent  ccf50f0299992c40caebb035135b4c23861110bd
Add function "smooth3" (patch #8864)

* scripts/plot/draw/smooth3.m: New function.
* scripts/plot/draw/module.mk: Update.
* __unimplemented__.m: Remove "smooth3" from list.
* doc/interpreter/plot.txi: Add to documentation.

diff -r ccf50f029999 -r f2c92d4a173e doc/interpreter/plot.txi
--- a/doc/interpreter/plot.txi	Wed Jul 06 14:57:54 2016 -0400
+++ b/doc/interpreter/plot.txi	Wed Jul 06 22:00:10 2016 +0200
@@ -384,6 +384,8 @@
 
 @DOCSTRING(isocolors)
 
+@DOCSTRING(smooth3)
+
 @DOCSTRING(shrinkfaces)
 
 @DOCSTRING(diffuse)
diff -r ccf50f029999 -r f2c92d4a173e scripts/help/__unimplemented__.m
--- a/scripts/help/__unimplemented__.m	Wed Jul 06 14:57:54 2016 -0400
+++ b/scripts/help/__unimplemented__.m	Wed Jul 06 22:00:10 2016 +0200
@@ -805,7 +805,6 @@
   "serial",
   "setpixelposition",
   "showplottool",
-  "smooth3",
   "snapnow",
   "splitapply",
   "ss2tf",
diff -r ccf50f029999 -r f2c92d4a173e scripts/plot/draw/module.mk
--- a/scripts/plot/draw/module.mk	Wed Jul 06 14:57:54 2016 -0400
+++ b/scripts/plot/draw/module.mk	Wed Jul 06 22:00:10 2016 +0200
@@ -82,6 +82,7 @@
   scripts/plot/draw/semilogy.m \
   scripts/plot/draw/shrinkfaces.m \
   scripts/plot/draw/slice.m \
+  scripts/plot/draw/smooth3.m \
   scripts/plot/draw/sombrero.m \
   scripts/plot/draw/sphere.m \
   scripts/plot/draw/stairs.m \
diff -r ccf50f029999 -r f2c92d4a173e scripts/plot/draw/smooth3.m
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/scripts/plot/draw/smooth3.m	Wed Jul 06 22:00:10 2016 +0200
@@ -0,0 +1,239 @@
+## Copyright (C) 2016 Markus Muetzel
+##
+## This file is part of Octave.
+##
+## Octave 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 3 of the License, or (at
+## your option) any later version.
+##
+## Octave 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 Octave; see the file COPYING.  If not, see
+## <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn  {Function File} {@var{smoothed_data} =} smooth3 (@var{data})
+## @deftypefnx {Function File} {@var{smoothed_data} =} smooth3 (@var{data}, @var{method})
+## @deftypefnx {Function File} {@var{smoothed_data} =} smooth3 (@var{data}, @var{method}, @var{sz})
+## @deftypefnx {Function File} {@var{smoothed_data} =} smooth3 (@var{data}, @var{method}, @var{sz}, @var{st_dev})
+## Smooth values of 3-dimensional matrix @var{data}.
+##
+## @var{data} must be a non-singleton 3-dimensional matrix.  The smoothed data
+## of this matrix is returned in @var{smoothed_data} which is of the same size
+## as @var{data}.
+##
+## Optionally, @var{method} determines which convolution kernel is used for the
+## smoothing process.  Choose one of:
+## @table @asis
+## @item @code{box}, @code{b} (default)
+## to use a convolution kernel with sharp edges.
+## @item @code{gaussian}, @code{g}
+## to use a convolution kernel that is represented by a non-correlated
+## trivariate normal distribution function.
+## @end table
+##
+## Optionally, @var{sz} can be either a vector of 3 elements representing the
+## size of the convolution kernel in x-, y- and z-direction or a scalar in
+## which case the same size is used in all three dimensions.
+## The default value is 3.
+##
+## When @var{method} is @code{gaussian}, @var{st_dev} can define the standard
+## deviation of the trivariate normal distribution function.  @var{st_dev} can
+## be either a vector of 3 elements representing the standard deviation of the
+## gaussian convolution kernel in x-, y- and z-direction or a scalar in which
+## case the same value is used in all three dimensions.  The default value is
+## 0.65.
+##
+## Example:
+## @example
+## @group
+## data = rand (10, 10, 10);
+## figure;
+## subplot (1, 2, 1)
+## patch (isosurface (data, .5), "FaceColor", "blue", "EdgeColor", "k");
+## title ("Original data")
+## view(3)
+## smoothed_data = smooth3 (data);
+## subplot (1, 2, 2)
+## patch (isosurface (smoothed_data, .5), "FaceColor", "blue", "EdgeColor", "k");
+## title ("Smoothed data")
+## view(3)
+## @end group
+## @end example
+##
+## @seealso{isosurface, isonormals, patch}
+## @end deftypefn
+
+## Author: mmuetzel
+
+function smoothed_data = smooth3 (data, method, sz, st_dev)
+
+  if (nargin < 1 || nargin > 4 || nargout > 1)
+    print_usage ();
+  endif
+
+  if (nargin < 2), method = "box"; endif
+  if (nargin < 3), sz = 3; endif
+  if (nargin < 4), st_dev = .65; endif
+
+  [data, conv_kernel, sz, st_dev] = ...
+                        __get_check_smooth3_args__ (data, method, sz, st_dev);
+
+  ## manually pad data by replicating the values at the edges
+  ## (convn would pad with zeros)
+  idx = cell (3, 1);
+  for i_dim = 1:3
+    sz_dim = size (data, i_dim);
+    pad_vec = ones (1, (sz(i_dim)-1)/2);
+    idx{i_dim} = [pad_vec 1:sz_dim sz_dim*pad_vec];
+  endfor
+  data_padded = data(idx{:});
+
+  ## actual smoothing
+  smoothed_data = convn (data_padded, conv_kernel, "valid");
+
+endfunction
+
+function [data, conv_kernel, sz, st_dev] = __get_check_smooth3_args__ (data, method, sz, st_dev);
+
+  ## check data
+  if (ndims (data) != 3)
+    error("smooth3: DATA must have 3 dimensions");
+  endif
+
+  ## check sz
+  if (isscalar (sz))
+    sz(1:3) = sz;
+  endif
+
+  if (numel (sz) != 3)
+    error (["smooth3: the size SZ of the convolution kernel must either " ...
+            "be a scalar or a vector of length 3"]);
+  endif
+
+  if (any (sz < 1) || any (rem (sz, 2) != 1))
+    error (["smooth3: the size SZ of the convolution kernel must consist " ...
+            "of positive odd integers"]);
+  endif
+
+  ## check method
+  switch lower(method)
+    case {"g", "gaussian"}
+      ## check st_dev
+      if (isscalar (st_dev))
+        st_dev(1:3) = st_dev;
+      endif
+
+      if (numel (st_dev) != 3)
+        error (["smooth3: the standard deviation of the gaussian " ...
+                "convolution kernel must either be a scalar or a vector " ...
+                "of length 3"]);
+      endif
+
+      conv_kernel = __smooth3_gaussian3__ (sz, st_dev);
+    case {"b", "box"}
+      conv_kernel = ones (sz) / prod (sz);
+    otherwise
+      error ("smooth3: METHOD '%s' unknown", method)
+  endswitch
+
+endfunction
+
+function gaussian3 = __smooth3_gaussian3__ (sz, st_dev)
+  ## trivariate non-correlated gaussian distribution function
+
+  x = (-(sz(2)-1)/2:(sz(2)-1)/2) / st_dev(2);
+  y = (-(sz(1)-1)/2:(sz(1)-1)/2) / st_dev(1);
+  z = (-(sz(3)-1)/2:(sz(3)-1)/2) / st_dev(3);
+
+  [xx, yy, zz] = meshgrid (x, y, z);
+
+  gaussian3 = exp (-(xx.*xx + yy.*yy + zz.*zz)/2);
+
+  gaussian3 = gaussian3 / sum (gaussian3(:)); ## normalize
+
+endfunction
+
+
+## one input argument (method: "box")
+%!test
+%! a = rand(10, 8, 7);
+%! b = smooth3(a);
+%! assert (size_equal (a, b), true);
+
+## two input argument (method: "gaussian")
+%!test
+%! a = rand(5, 8, 7);
+%! b = smooth3(a, "gaussian");
+%! assert (size_equal (a, b), true);
+
+## three input argument (method: "box")
+%!test
+%! a = rand(3, 8, 7);
+%! b = smooth3(a, "box", 5);
+%! assert (size_equal (a, b), true);
+
+## three input argument (method: "gaussian")
+%!test
+%! a = rand(3, 8, 7);
+%! b = smooth3(a, "gaussian", 7);
+%! assert (size_equal (a, b), true);
+
+## size of convolution kernel = 1: no smoothing (method: "box")
+%!test
+%! a = rand(9, 8, 7);
+%! b = smooth3(a, "box", 1);
+%! assert (a, b);
+
+## size of convolution kernel = 1: no smoothing (method: "gaussian")
+%!test
+%! a = rand(9, 8, 7);
+%! b = smooth3(a, "gaussian", 1);
+%! assert (a, b);
+
+## four input arguments (method: "gaussian")
+%!test
+%! a = rand(3, 8, 7);
+%! b = smooth3(a, "gaussian", 7, .5);
+%! assert (size_equal (a, b), true);
+
+## size of convolution kernel is different in x, y and z (method: "box")
+%!test
+%! a = rand(3, 8, 7);
+%! b = smooth3(a, "box", [5 3 7]);
+%! assert (size_equal (a, b), true);
+
+## size of convolution kernel is different in x, y and z (method: "gaussian")
+%!test
+%! a = rand(3, 8, 7);
+%! b = smooth3(a, "gaussian", [5 3 7]);
+%! assert (size_equal (a, b), true);
+
+## size and width of gaussian convolution kernel is different in x, y and z
+## (method: "gaussian")
+%!test
+%! a = rand(3, 8, 7);
+%! b = smooth3(a, "gaussian", [7 3 5], [.3 .5 .4]);
+%! assert (size_equal (a, b), true);
+
+## test for each error
+%!test
+%!error <DATA must have 3 dimensions> a = rand(3, 8); b = smooth3(a);
+%!error <the size SZ of the convolution kernel must either be a scalar or a vector of length 3>
+%! a = rand(3, 8, 3); b = smooth3(a, "box", [3 5]);
+%!error <the size SZ of the convolution kernel must consist of positive odd integers>
+%! a = rand(3, 8, 3); b = smooth3(a, "box", [3 2 5]);
+%!error <the size SZ of the convolution kernel must consist of positive odd integers>
+%! a = rand(3, 8, 3); b = smooth3(a, "box", [3 0 5]);
+%!error <the size SZ of the convolution kernel must consist of positive odd integers>
+%! a = rand(3, 8, 3); b = smooth3(a, "box", [3 2.5 5]);
+%!error <the standard deviation of the gaussian convolution kernel must either be a scalar or a vector of length 3>
+%! a = rand(3, 8, 3); b = smooth3(a, "gaussian", 3, [.3 .4]);
+%!error <METHOD 'other' unknown> a = rand(3, 8, 3); b = smooth3(a, "other");
+%!error <Invalid call to smooth3>
+%! a = rand(3, 8, 3); b = smooth3(a, "box", [7 3 5], .4, 1);
