--- a/scripts/statistics/base/std.m
+++ b/scripts/statistics/base/std.m
@@ -19,6 +19,8 @@
 
 ## -*- texinfo -*-
 ## @deftypefn {Function File} {} std (@var{x})
+## @deftypefnx {Function File} {} std (@var{x}, @var{opt})
+## @deftypefnx {Function File} {} std (@var{x}, @var{opt}, @var{dim})
 ## If @var{x} is a vector, compute the standard deviation of the elements
 ## of @var{x}.
 ## @iftex
@@ -38,26 +40,50 @@
 ## @end ifinfo
 ## If @var{x} is a matrix, compute the standard deviation for
 ## each column and return them in a row vector.
+##
+## The argument @var{opt} determines the type of normalization to use. Valid values
+## are
+##
+## @table @asis 
+## @item 0:
+##   normalizes with N-1, provides the square root of best unbiased estimator of 
+##   the variance [default]
+## @item 1:
+##   normalizes with N, this provides the square root of the second moment around 
+##   the mean
+## @end table
+##
+## The third argument @var{dim} determines the dimension along which the standard
+## deviation is calculated.
 ## @end deftypefn
 ## @seealso{mean and median}
 
 ## Author: jwe
 
-function retval = std (a)
+function retval = std (a, opt, dim)
 
-  if (nargin != 1)
-    usage ("std (a)");
+  if (nargin < 1 || nargin > 3)
+    usage ("std (a, opt, dim)");
+  endif
+  if nargin < 3
+    dim = min(find(size(a)>1));
+    if isempty(dim), dim=1; endif;
+  endif
+  if ((nargin < 2) || isempty(opt))
+    opt = 0;
   endif
 
-  nr = rows (a);
-  nc = columns (a);
-  if (nc == 1 && nr == 1)
-    retval = 0;
-  elseif (nc == 1 || nr == 1)
-    n = length (a);
-    retval = sqrt (sumsq (a - mean (a)) / (n - 1));
-  elseif (nr > 1 && nc > 0)
-    retval = sqrt (sumsq (a - ones (nr, 1) * mean (a)) / (nr - 1));
+  sz = size(a);
+  if (sz (dim) == 1)
+    retval = zeros(sz);
+  elseif (numel (a) > 0)
+    rng = ones (1, length (sz));
+    rng (dim) = sz (dim);
+    if (opt == 0)
+      retval = sqrt (sumsq (a - repmat(mean (a, dim), rng), dim) / (sz(dim) - 1));
+    else
+      retval = sqrt (sumsq (a - repmat(mean (a, dim), rng), dim) / sz(dim));
+    endif
   else
     error ("std: invalid matrix argument");
   endif