changeset 9726:b7b89061bd0e

reimplement median using nth_element
author Jaroslav Hajek <highegg@gmail.com>
date Thu, 15 Oct 2009 08:38:07 +0200
parents aea3a3a950e1
children 04386b72d3df
files scripts/ChangeLog scripts/statistics/base/median.m
diffstat 2 files changed, 11 insertions(+), 29 deletions(-) [+]
line wrap: on
line diff
--- a/scripts/ChangeLog
+++ b/scripts/ChangeLog
@@ -1,3 +1,7 @@
+2009-10-14  Jaroslav Hajek  <highegg@gmail.com>
+
+	* statistics/base/median.m: Rewrite using nth_element.
+
 2009-10-01  John W. Eaton  <jwe@octave.org>
 
 	* image/__img__.m: Adjust xlim and ylim correctly.
--- a/scripts/statistics/base/median.m
+++ b/scripts/statistics/base/median.m
@@ -1,5 +1,6 @@
 ## Copyright (C) 1996, 1997, 1998, 1999, 2000, 2004, 2005, 2006, 2007,
 ##               2008, 2009 John W. Eaton
+## Copyright (C) 2009 VZLU Prague
 ##
 ## This file is part of Octave.
 ##
@@ -53,39 +54,16 @@
     print_usage ();
   endif
   if (nargin < 2)
-    dim = find (size (a) > 1, 1);
-    if (isempty (dim))
-      dim = 1;
-    endif
+    dim = [find(size (a) != 1, 1), 1](1); # First non-singleton dim.
   endif
 
-  sz = size (a);
-  s = sort (a, dim);
   if (numel (a) > 0)
-    if (numel (a) == sz(dim))
-      if (rem (sz(dim), 2) == 0)
-	i = sz(dim) / 2;
-	retval = (s(i) + s(i+1)) / 2;
-      else
-	i = ceil (sz(dim) /2);
-	retval = s(i);
-      endif
+    n = size (a, dim);
+    k = floor ((n+1) / 2);
+    if (mod (n, 2) == 1)
+      retval = nth_element (a, k, dim);
     else
-      idx = cell ();
-      nd = length (sz);
-      for i = 1:nd
-	idx{i} = 1:sz(i);
-      endfor
-      if (rem (sz(dim), 2) == 0)
-	i = sz(dim) / 2;
-	idx{dim} = i;
-	retval = s(idx{:});
-	idx{dim} = i+1;
-	retval = (retval + s(idx{:})) / 2;
-      else
-	idx{dim} = ceil (sz(dim) / 2);
-	retval = s(idx{:});
-      endif
+      retval = mean (nth_element (a, k:k+1, dim), dim);
     endif
   else
     error ("median: invalid matrix argument");