octave_packages/nan-2.5.5/nanstd.m

   1 function [y] = nanstd(x,FLAG,DIM)
   2 % NANSTD same as STD but ignores NaN's.
   3 % NANSTD is OBSOLETE; use NaN/STD instead. NANSTD is included
   4 %    to fix a bug in alternative implementations and to
   5 %    provide some compatibility.
   6 %
   7 % Y = nanstd(x, FLAG, [,DIM])
   8 %
   9 % x     data
  10 % FLAG  0: [default] normalizes with (N-1), N = sample size
  11 % FLAG  1: normalizes with N, N = sample size
  12 % DIM   dimension
  13 %       1 sum of columns
  14 %       2 sum of rows
  15 %       default or []: first DIMENSION with more than 1 element
  16 % Y     resulting standard deviation
  17 %
  18 % see also: SUM, SUMSKIPNAN, NANSUM, STD
  19
  20 %    $Id: nanstd.m 9033 2011-11-08 20:58:07Z schloegl $
  21 %    Copyright (C) 2000-2003,2006,2008,2009,2010 by Alois Schloegl <alois.schloegl@gmail.com>
  22 %    This is part of the NaN-toolbox. For more details see
  23 %       http://pub.ist.ac.at/~schloegl/matlab/NaN/
  24
  25 %    This program is free software; you can redistribute it and/or modify
  26 %    it under the terms of the GNU General Public License as published by
  27 %    the Free Software Foundation; either version 3 of the License, or
  28 %    (at your option) any later version.
  29 %
  30 %    This program is distributed in the hope that it will be useful,
  31 %    but WITHOUT ANY WARRANTY; without even the implied warranty of
  32 %    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  33 %    GNU General Public License for more details.
  34 %
  35 %    You should have received a copy of the GNU General Public License
  36 %    along with this program; If not, see <http://www.gnu.org/licenses/>.
  37
  38 if nargin<2,
  39         FLAG = 0;
  40 end;
  41
  42 if nargin<3,
  43         DIM = [];
  44 end;
  45 if isempty(FLAG),
  46         FLAG = 0;
  47 end;
  48 if isempty(DIM),
  49         DIM = find(size(x)>1,1);
  50         if isempty(DIM), DIM=1; end;
  51 end;
  52
  53 [y,n,ssq] = sumskipnan(x,DIM);
  54 if all(ssq(:).*n(:) > 2*(y(:).^2)),
  55         %% rounding error is neglectable
  56         y = ssq - y.*y./n;
  57 else
  58         %% rounding error is not neglectable
  59         [y,n] = sumskipnan(center(x,DIM).^2,DIM);
  60 end;
  61
  62 if (FLAG==1)
  63         y = sqrt(y)./n; % normalize with N
  64 else
  65         % default method
  66         y = sqrt(y./max(n-1,0));        % normalize with N-1
  67 end;
  68
  69
  70 %!assert(nanstd(0),NaN)
  71