changeset 15697:9811b32b645e

maint: merge in Faddeeva functions and relatives
author Jordi Gutiérrez Hermoso <jordigh@octave.org>
date Thu, 29 Nov 2012 10:40:32 -0500
parents 310db14aa605 (diff) 2fac72a256ce (current diff)
children ea2788117ee3
files
diffstat 12 files changed, 212 insertions(+), 134 deletions(-) [+]
line wrap: on
line diff
--- a/doc/interpreter/dynamic.txi
+++ b/doc/interpreter/dynamic.txi
@@ -1119,7 +1119,7 @@
 allocate any memory that is needed by the foreign code, with either the
 @nospell{fortran_vec} method or the @w{@code{OCTAVE_LOCAL_BUFFER}} macro.
 
-The Octave unwind_protect mechanism (@ref{The @code{unwind_protect} Statement})
+The Octave unwind_protect mechanism (@ref{The unwind_protect Statement})
 can also be used in oct-files.  In conjunction with the exception
 handling of Octave, it is important to enforce that certain code is run
 to allow variables, etc. to be restored even if an exception occurs.  An
--- a/doc/interpreter/errors.txi
+++ b/doc/interpreter/errors.txi
@@ -44,7 +44,7 @@
 Since an error can occur during the evaluation of a program, it is
 very convenient to be able to detect that an error occurred, so that
 the error can be fixed.  This is possible with the @code{try} statement
-described in @ref{The @code{try} Statement}.
+described in @ref{The try Statement}.
 
 @menu
 * Raising Errors::
@@ -136,7 +136,7 @@
 @subsection Catching Errors
 
 When an error occurs, it can be detected and handled using the
-@code{try} statement as described in @ref{The @code{try} Statement}.
+@code{try} statement as described in @ref{The try Statement}.
 As an example, the following piece of code counts the number of errors
 that occurs during a @code{for} loop.
 
--- a/doc/interpreter/octave.texi
+++ b/doc/interpreter/octave.texi
@@ -367,15 +367,15 @@
 
 Statements
 
-* The @code{if} Statement::            
-* The @code{switch} Statement::        
-* The @code{while} Statement::         
-* The @code{do-until} Statement::      
-* The @code{for} Statement::           
-* The @code{break} Statement::         
-* The @code{continue} Statement::      
-* The @code{unwind_protect} Statement::  
-* The @code{try} Statement::           
+* The if Statement::            
+* The switch Statement::        
+* The while Statement::         
+* The do-until Statement::      
+* The for Statement::           
+* The break Statement::         
+* The continue Statement::      
+* The unwind_protect Statement::  
+* The try Statement::           
 * Continuation Lines::          
 
 The @code{switch} Statement
--- a/doc/interpreter/stmt.txi
+++ b/doc/interpreter/stmt.txi
@@ -46,20 +46,20 @@
 @dfn{body} of a control statement.
 
 @menu
-* The @code{if} Statement::            
-* The @code{switch} Statement::        
-* The @code{while} Statement::         
-* The @code{do-until} Statement::      
-* The @code{for} Statement::           
-* The @code{break} Statement::         
-* The @code{continue} Statement::      
-* The @code{unwind_protect} Statement::  
-* The @code{try} Statement::           
+* The if Statement::            
+* The switch Statement::        
+* The while Statement::         
+* The do-until Statement::      
+* The for Statement::           
+* The break Statement::         
+* The continue Statement::      
+* The unwind_protect Statement::  
+* The try Statement::           
 * Continuation Lines::          
 @end menu
 
-@node The @code{if} Statement
-@section The @code{if} Statement
+@node The if Statement
+@section The if Statement
 @cindex @code{if} statement
 @cindex @code{else} statement
 @cindex @code{elseif} statement
@@ -207,8 +207,8 @@
 using the indentation to show how Octave groups the statements.
 @xref{Functions and Scripts}.
 
-@node The @code{switch} Statement
-@section The @code{switch} Statement
+@node The switch Statement
+@section The switch Statement
 @cindex @code{switch} statement
 @cindex @code{case} statement
 @cindex @code{otherwise} statement
@@ -385,8 +385,8 @@
 
 @end itemize
 
-@node The @code{while} Statement
-@section The @code{while} Statement
+@node The while Statement
+@section The while Statement
 @cindex @code{while} statement
 @cindex @code{endwhile} statement
 @cindex loop
@@ -454,8 +454,8 @@
 body; but using one makes the program clearer unless the body is very
 simple.
 
-@node The @code{do-until} Statement
-@section The @code{do-until} Statement
+@node The do-until Statement
+@section The do-until Statement
 @cindex @code{do-until} statement
 
 The @code{do-until} statement is similar to the @code{while} statement,
@@ -501,8 +501,8 @@
 body; but using one makes the program clearer unless the body is very
 simple.
 
-@node The @code{for} Statement
-@section The @code{for} Statement
+@node The for Statement
+@section The for Statement
 @cindex @code{for} statement
 @cindex @code{endfor} statement
 
@@ -658,8 +658,8 @@
 structure elements when the names of the elements do not need to be
 known.
 
-@node The @code{break} Statement
-@section The @code{break} Statement
+@node The break Statement
+@section The break Statement
 @cindex @code{break} statement
 
 The @code{break} statement jumps out of the innermost @code{for} or
@@ -713,8 +713,8 @@
 @end group
 @end example
 
-@node The @code{continue} Statement
-@section The @code{continue} Statement
+@node The continue Statement
+@section The continue Statement
 @cindex @code{continue} statement
 
 The @code{continue} statement, like @code{break}, is used only inside
@@ -762,8 +762,8 @@
 @end group
 @end example
 
-@node The @code{unwind_protect} Statement
-@section The @code{unwind_protect} Statement
+@node The unwind_protect Statement
+@section The unwind_protect Statement
 @cindex @code{unwind_protect} statement
 @cindex @code{unwind_protect_cleanup}
 @cindex @code{end_unwind_protect}
@@ -813,8 +813,8 @@
 point of the error and the statement to restore the value would not be
 executed.
 
-@node The @code{try} Statement
-@section The @code{try} Statement
+@node The try Statement
+@section The try Statement
 @cindex @code{try} statement
 @cindex @code{catch}
 @cindex @code{end_try_catch}
--- a/doc/interpreter/vectorize.txi
+++ b/doc/interpreter/vectorize.txi
@@ -635,7 +635,7 @@
 
 If you are using @code{eval} merely as an exception handling mechanism, and not
 because you need to execute some arbitrary text, use the @code{try}
-statement instead.  @xref{The @code{try} Statement}.
+statement instead.  @xref{The try Statement}.
 
 @item Use @code{ignore_function_time_stamp} when appropriate.
 If you are calling lots of functions, and none of them will need to change
--- a/libinterp/octave-value/ov-oncleanup.cc
+++ b/libinterp/octave-value/ov-oncleanup.cc
@@ -192,7 +192,7 @@
 elements) or returned from a function, @var{action} will be executed after\n\
 clearing the last copy of the object.  Note that if multiple local onCleanup\n\
 variables are created, the order in which they are called is unspecified.\n\
-For similar functionality @xref{The @code{unwind_protect} Statement}.\n\
+For similar functionality @xref{The unwind_protect Statement}.\n\
 @end deftypefn")
 {
   octave_value retval;
--- a/scripts/image/gray2ind.m
+++ b/scripts/image/gray2ind.m
@@ -28,16 +28,16 @@
 ## If not given @var{n} defaults to 64 for grayscale images or 2 for
 ## binary black and white images.
 ##
-## The output @var{img} is of class uint8 if @var{n} is less than or equal to
-## 256.  Otherwise, the output is of class uint16.
-## @seealso{ind2gray, rgb2ind} 
+## The output @var{img} is of class uint8 or uint 16 if @var{n} is less than or
+## equal to 256 or 65536 respectively.  Otherwise, the output is of class double.
+## @seealso{ind2gray, rgb2ind}
 ## @end deftypefn
 
 ## Author: Tony Richardson <arichard@stark.cc.oh.us>
 ## Created: July 1994
 ## Adapted-By: jwe
 
-function [img, map] = gray2ind (I, n = 64)
+function [I, map] = gray2ind (I, n = 64)
 
   if (nargin < 1 || nargin > 2)
     print_usage ();
@@ -45,55 +45,59 @@
     error ("gray2ind: I must be a grayscale or binary image");
   elseif (! isscalar (n) || n < 1 || n > 65536)
     error ("gray2ind: N must be a positive integer in the range [1, 65536]");
+  elseif (! ismatrix (I) || ndims (I) != 2)
+    error ("gray2ind: first input argument must be a gray scale image");
   endif
+
+  ## default n is different if image is logical
+  if (nargin == 1 && islogical (I))
+    n = 2;
+  endif
+
+  if (! isscalar (n) || n < 0)
+    error ("gray2ind: second input argument must be a positive integer");
+  endif
+
   cls = class (I);
   if (! any (isa (I, {"logical", "uint8", "uint16", "int16", "single", "double"})))
     error ("gray2ind: invalid data type '%s'", cls);
-  elseif (isfloat (I) && (any (I(:) < 0) || any (I(:) > 1)))
+  elseif (isfloat (I) && (min (I(:) < 0) || max (I(:) > 1)))
     error ("gray2ind: floating point images may only contain values between 0 and 1");
   endif
 
-  ## Create grayscale colormap
-  if (nargin == 1 && islogical (I))
-    n = 2;
-  endif
   map = gray (n);
 
   ## Set up scale factor
   if (isinteger (I))
-    low = intmin (cls);
-    if (low != 0)
-      I -= low;  # shift to zero-based indexing
-    endif
-    scale = double (intmax (cls)) - double (low);
+    low   = double (intmin (cls));
+    scale = double (intmax (cls)) - low;
+    I = double (I) - low;
   else
-    scale = 1;  # floating point doesn't need scaling
+    scale = 1;
   endif
-
-  ## Scale image
-  img = I * ((n-1)/scale);
-
-  ## Round and convert to appropriate output type.
   ## Note: no separate call to round () necessary because
   ##       type conversion does that automatically.
-  if (n <= 256)
-    img = uint8 (img) + 1;
+  I = I * ((n-1)/scale);
+  if (n < 256)
+    I = uint8 (I);
+  elseif (n < 65536)
+    I = uint16 (I);
   else
-    img = uint16 (img) + 1;
+    ## if uint16 is not enough, we return double in which case index
+    ## values should start at 1.
+    I = round (I) + 1;
   endif
-
 endfunction
 
-
-%!assert (gray2ind ([0 0.25 0.5 1]), uint8 ([1 17 33 64]))
-%!assert (gray2ind ([0 0.25 0.5 1], 400), uint16 ([1 101 201 400]))
-%!assert (gray2ind (logical ([1 0 0 1])), uint8 ([2 1 1 2]))
-%!assert (gray2ind (uint8 ([0 64 128 192 255])), uint8 ([1 17 33 48 64]))
+%!assert (gray2ind ([0 0.25 0.5 1]), uint8 ([0 16 32 63]))
+%!assert (gray2ind ([0 0.25 0.5 1], 400), uint16 ([0 100 200 399]))
+%!assert (gray2ind (logical ([1 0 0 1])), uint8 ([1 0 0 1]))
+%!assert (gray2ind (uint8 ([0 64 128 192 255])), uint8 ([0 16 32 47 63]))
 
 %!test
 %! i2g = ind2gray (1:100, gray (100));
 %! g2i = gray2ind (i2g, 100);
-%! assert (g2i, uint8 (1:100));
+%! assert (g2i, uint8 (0:99));
 
 %% Test input validation
 %!error gray2ind ()
--- a/scripts/image/ind2gray.m
+++ b/scripts/image/ind2gray.m
@@ -20,8 +20,10 @@
 ## @deftypefn  {Function File} {@var{I} =} ind2gray (@var{x})
 ## @deftypefnx {Function File} {@var{I} =} ind2gray (@var{x}, @var{map})
 ## Convert a color indexed image to a grayscale intensity image.
+##
 ## If @var{map} is omitted, the current colormap is used to determine the
-## intensities.
+## intensities.  If it doesn't contain enough colors, it is padded it with the
+## last color in the map.
 ##
 ## The output @var{I} is of the same class as the input @var{x} and may be
 ## one of @code{uint8}, @code{uint16}, @code{single}, or @code{double}.
@@ -38,32 +40,19 @@
 ## Created: July 1994
 ## Adapted-By: jwe
 
-function I = ind2gray (x, map)
+function I = ind2gray (x, map = colormap ())
 
   if (nargin < 1 || nargin > 2)
     print_usage ();
   endif
-
-  if (! isreal (x) || issparse (x)
-      || (isfloat (x) && (any (x(:) < 1 || any (x(:) != fix (x(:)))))))
-    error ("ind2gray: X must be an indexed image");
-  endif
-  cls = class (x);
-  if (! any (isa (x, {"logical", "uint8", "uint16", "single", "double"})))
-    error ("ind2gray: invalid data type '%s'", cls);
-  endif
-
-  if (nargin == 1)
-    map = colormap ();
-  elseif (! iscolormap (map))
-    error ("ind2gray: MAP must be a valid colormap");
-  endif
+  [x, map] = ind2x ("ind2rgb", x, map);
 
   ## Convert colormap to luminance intensity values 
-  map *= [0.299; 0.587; 0.114];
+  map *= [0.29894; 0.58704; 0.11402];
 
-  ## Convert colormap to same class as that of input so that
-  ## indexing in colormap will produce output of the same type as the input.
+  ## Convert colormap to same class as that of input so that when reshape
+  ## later, produces output of the same type as the input
+  cls = class (x);
   if (isinteger (x)) 
     map *= intmax (cls);
   elseif (strcmp (cls, "single"))
@@ -75,12 +64,10 @@
 
 endfunction
 
-
-%!test
+%!shared i2g
 %! i2g = ind2gray (1:100, gray (100));
-%! assert (i2g, 0:1/99:1, eps);
-%! g2i = gray2ind (i2g, 100);
-%! assert (g2i, uint8 (1:100));
+%!assert (i2g, 0:1/99:1, eps);
+%!assert (gray2ind (i2g, 100), uint8 (0:99));
 
 %%test input validation
 %!error ind2gray ()
--- a/scripts/image/ind2rgb.m
+++ b/scripts/image/ind2rgb.m
@@ -36,35 +36,12 @@
 ## Created: July 1994
 ## Adapted-By: jwe
 
-function [R, G, B] = ind2rgb (x, map)
+function [R, G, B] = ind2rgb (x, map = colormap ())
 
   if (nargin < 1 || nargin > 2)
     print_usage ();
   endif
-
-  if (! isreal (x) || issparse (x)
-      || (isfloat (x) && (any (x(:) < 1 || any (x(:) != fix (x(:)))))))
-    error ("ind2rgb: X must be an indexed image");
-  endif
-  cls = class (x);
-  if (! any (isa (x, {"logical", "uint8", "uint16", "single", "double"})))
-    error ("ind2rgb: invalid data type '%s'", cls);
-  endif
-
-  if (nargin == 1)
-    map = colormap ();
-  elseif (! iscolormap (map))
-    error ("ind2rgb: MAP must be a valid colormap");
-  endif
-
-  ## Do we have enough colors in the color map?
-  maxidx = max (x(:));
-  rm = rows (map);
-  if (rm < maxidx)
-    ## Pad with the last color in the map.
-    pad = repmat (map(end,:), maxidx-rm, 1);
-    map(end+1:maxidx, :) = pad;
-  endif
+  [x, map] = ind2x ("ind2rgb", x, map);
 
   ## Compute result
   [row, col] = size (x);
@@ -80,8 +57,27 @@
 
 endfunction
 
-
-%% FIXME: Need some functional tests or %!demo blocks
+%!shared img, map, ergb, rgb, r, g, b
+%! img = [2 4 5; 3 2 5; 1 2 4];
+%! map = [0.0  0.0  0.0
+%!        0.2  0.4  0.6
+%!        0.4  0.4  0.5
+%!        0.3  0.7  1.0
+%!        0.1  0.5  0.8];
+%! ergb(:,:,1) = [0.2 0.3 0.1; 0.4 0.2 0.1; 0.0 0.2 0.3];
+%! ergb(:,:,2) = [0.4 0.7 0.5; 0.4 0.4 0.5; 0.0 0.4 0.7];
+%! ergb(:,:,3) = [0.6 1.0 0.8; 0.5 0.6 0.8; 0.0 0.6 1.0];
+%! ## test basic usage with 1 and 3 outputs
+%! [rgb] = ind2rgb (img, map);
+%! [r, g, b] = ind2rgb (img, map);
+%!assert (ergb, rgb);
+%!assert (ergb, reshape ([r(:) g(:) b(:)], [size(img) 3]));
+%! ## test correction for integers
+%! img = uint8 (img -1);
+%! [rgb] = ind2rgb (img, map);
+%!assert (ergb, rgb);
+%! ## check it fails when image is a float with an index value of 0
+%!fail("[rgb] = ind2rgb (double(img), map)")
 
 %%test input validation
 %!error ind2rgb ()
@@ -97,4 +93,3 @@
 %!error <MAP must be a valid colormap> ind2rgb (1, ones (2,4))
 %!error <MAP must be a valid colormap> ind2rgb (1, [-1])
 %!error <MAP must be a valid colormap> ind2rgb (1, [2])
-
new file mode 100644
--- /dev/null
+++ b/scripts/image/private/ind2x.m
@@ -0,0 +1,56 @@
+## Copyright (C) 1994-2012 John W. Eaton
+## Copyright (C) 2012 Carnë Draug
+##
+## 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/>.
+
+## private function for the ind2something functions which have a lot of code
+## in common
+
+function [x, map] = ind2x (name, x, map)
+
+  ## Check if X is an indexed image.
+  if (ndims (x) != 2 || issparse (x) || (isfloat (x) && ! isindex (x)) ||
+      ! ismember (class (x), {"double", "single", "uint8", "uint16"}))
+    error ("%s: X must be an indexed image", name);
+  endif
+
+  ## Check the color map.
+  if (! iscolormap (map))
+    error ("%s: MAP must be a valid colormap", name);
+  endif
+
+  ## Do we have enough colors in the color map?
+  ## there's an offset of 1 when the indexed image is an integer class so we fix
+  ## it now and convert it to float only if really necessary and even then only
+  ## to single precision since its enough for both uint8 and uint16
+  maxidx = max (x(:));
+  if (isinteger (x))
+    if (maxidx == intmax (class (x)))
+      x = single (x);
+    endif
+    x      += 1;
+    maxidx += 1;
+  endif
+
+  rm = rows (map);
+  if (rm < maxidx)
+    ## Pad with the last color in the map.
+    pad = repmat (map(end,:), maxidx-rm, 1);
+    map(end+1:maxidx, :) = pad;
+  endif
+
+endfunction
--- a/scripts/image/rgb2ind.m
+++ b/scripts/image/rgb2ind.m
@@ -1,4 +1,5 @@
 ## Copyright (C) 1994-2012 John W. Eaton
+## Copyright (C) 2012 Carnë Draug
 ##
 ## This file is part of Octave.
 ##
@@ -23,7 +24,6 @@
 ## @seealso{ind2rgb, rgb2hsv, rgb2ntsc}
 ## @end deftypefn
 
-## Bugs: The color map may have duplicate entries.
 ## FIXME: This function has a very different syntax than the Matlab one of the same name.
 ##        Octave function does no support N, MAP, DITHER, or TOL arguments
 ## Author: Tony Richardson <arichard@stark.cc.oh.us>
@@ -51,15 +51,38 @@
 
   x = reshape (1:numel (R), size (R));
 
-  map = [R(:), G(:), B(:)];
+  map    = unique([R(:) G(:) B(:)], "rows");
+  [~, x] = ismember ([R(:) G(:) B(:)], map, "rows");
+  x      = reshape (x, size (R));
 
+  ## a colormap is of class double and values between 0 and 1
+  switch class (R)
+    case {"single", "double", "logical"}
+      ## do nothing, return the same
+    case {"uint8", "uint16"}
+      map = double (map) / double (intmax (class (R)));
+    case "int16"
+      map = (double (im) + 32768) / 65535;
+    otherwise
+      error ("unsupported image class %s", im_class);
+  endswitch
+
+  ## we convert to the smallest class necessary to encode the image. Matlab
+  ## documentation does not mention what it does when uint16 is not enough...
+  ## When an indexed image is of integer class, there's a -1 offset to the
+  ## colormap, hence the adjustment
+  if (rows (map) < 256)
+    x = uint8 (x - 1);
+  elseif (rows (map) < 65536)
+    x = uint16 (x - 1);
+  else
+    ## leave it as double
+  endif
 endfunction
 
-
 %% FIXME: Need some functional tests or %!demo blocks
 
 %% Test input validation
 %!error rgb2ind ()
 %!error rgb2ind (1,2)
 %!error rgb2ind (1,2,3,4)
-
--- a/scripts/set/powerset.m
+++ b/scripts/set/powerset.m
@@ -19,11 +19,17 @@
 ## -*- texinfo -*-
 ## @deftypefn  {Function File} {} powerset (@var{a})
 ## @deftypefnx {Function File} {} powerset (@var{a}, "rows")
+## Compute the powerset (all subsets) of the set @var{a}.
 ##
-## Return a cell array containing all subsets of the set @var{a}.
+## The set @var{a} must be a numerical matrix or a cell array of strings.  The
+## output will always be a cell array of either vectors or strings.
 ##
+## With the optional second argument @code{"rows"}, each row of the set @var{a}
+## is considered one element of the set.  As a result, @var{a} must then be a
+## numerical 2D matrix.
+##
+## @seealso{unique, union, setxor, setdiff, ismember}
 ## @end deftypefn
-## @seealso{unique, union, setxor, setdiff, ismember}
 
 function p = powerset (a, byrows_arg)
 
@@ -31,7 +37,7 @@
 
   if (nargin == 2)
     if (! strcmpi (byrows_arg, "rows"))
-      error ('powerset: expecting third argument to be "rows"');
+      error ('powerset: expecting second argument to be "rows"');
     elseif (iscell (a))
       warning ('powerset: "rows" not valid for cell arrays');
     else
@@ -40,6 +46,9 @@
   elseif (nargin != 1)
     print_usage ();
   endif
+  if (iscell (a) && ! iscellstr (a))
+    error ("powerset: cell arrays can only used for character strings");
+  endif
 
   if (byrows)
     a = unique (a, byrows_arg);
@@ -50,7 +59,7 @@
   endif
 
   if (n == 0)
-    p = [];
+    p = {};
   else
     if (n > 32)
       error ("powerset: not implemented for more than 32 elements");
@@ -77,8 +86,12 @@
 endfunction
 
 
-%!test
+%!shared c, p
 %! c = sort (cellstr ({ [], [1], [2], [3], [1, 2], [1, 3], [2, 3], [1, 2, 3]}));
 %! p = sort (cellstr (powerset ([1, 2, 3])));
-%! assert (p, c);
+%!assert (p, c);
+%! c = sort (cellstr ({ [], [1:3], [2:4], [3:5], [1:3; 2:4], [1:3; 3:5], [2:4; 3:5], [1:3; 2:4; 3:5]}));
+%! p = sort (cellstr (powerset ([1:3;2:4;3:5], "rows")));
+%!assert (p,c);
+%!assert (powerset([]), {});  # always return a cell array