--- a/liboctave/ChangeLog	Fri Mar 28 14:52:27 2008 -0400
+++ b/liboctave/ChangeLog	Fri Mar 28 15:53:09 2008 -0400
@@ -1,3 +1,7 @@
+2008-03-27  Jaroslav Hajek  <highegg@gmail.com>
+
+	* oct-lookup.h: New file.
+
 2008-03-26  David Bateman  <dbateman@feee.fr>
 
 	* Array.cc (assignN): Additional fix for vector assignments.

--- a/liboctave/Makefile.in	Fri Mar 28 14:52:27 2008 -0400
+++ b/liboctave/Makefile.in	Fri Mar 28 15:53:09 2008 -0400
@@ -81,8 +81,8 @@
 	lo-ieee.h lo-mappers.h lo-math.h lo-specfun.h \
 	lo-sysdep.h lo-utils.h mach-info.h md5.h oct-alloc.h oct-cmplx.h \
 	oct-env.h oct-fftw.h oct-getopt.h oct-group.h oct-inttypes.h \
-	oct-md5.h oct-passwd.h oct-rand.h oct-rl-edit.h oct-rl-hist.h \
-	oct-shlib.h oct-sort.h oct-spparms.h oct-syscalls.h \
+	oct-lookup.h oct-md5.h oct-passwd.h oct-rand.h oct-rl-edit.h \
+	oct-rl-hist.h oct-shlib.h oct-sort.h oct-spparms.h oct-syscalls.h \
 	oct-sparse.h oct-time.h oct-uname.h \
 	pathlen.h pathsearch.h prog-args.h \
 	randgamma.h randmtzig.h randpoisson.h \

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/liboctave/oct-lookup.h	Fri Mar 28 15:53:09 2008 -0400
@@ -0,0 +1,174 @@
+/*
+
+Copyright (C) 2008 VZLU Prague, a.s., Czech Republic
+
+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/>.
+
+*/
+
+// Author: Jaroslav Hajek <highegg@gmail.com>
+
+#if !defined (octave_oct_lookup)
+#define octave_oct_lookup 1
+
+#include <algorithm>
+#include <functional>
+
+#include "oct-types.h"
+
+// a simple binary lookup
+template<typename T, typename bpred>
+octave_idx_type
+bin_lookup (const T *table, octave_idx_type size, 
+            const T& val,
+            bpred comp)
+{
+  return std::upper_bound (table, table + size, val, comp) - table;
+}
+
+// version using < operator
+template<typename T>
+octave_idx_type
+bin_lookup (const T *table, octave_idx_type size,
+            const T& val)
+{
+  return std::upper_bound (table, table + size, val) - table;
+}
+
+// a unary functor that checks whether a value is outside [a,b) range
+template<class T, class bpred>
+class out_range : public std::unary_function<T, bool>
+{
+public:
+  out_range (const T& aa, const T& bb, const bpred& ccomp) 
+    : a(aa), b(bb), comp(ccomp) { }
+
+  bool operator() (const T& x) { return comp (x, a) || ! comp (x, b); }
+
+private:
+  T a;
+  T b;
+
+  bpred comp;
+};
+
+// conveniently constructs the above functor
+// NOTE: with SGI extensions, this can be written as
+// compose2 (logical_and<bool>(), 
+//           bind2nd (less<T>(), a),
+//           not1 (bind2nd (less<T>(), b)))
+template<class T, class bpred>
+out_range<T, bpred> 
+chk_out_range (const T& a, const T& b, bpred comp)
+{
+  return out_range<T, bpred> (a, b, comp);
+}
+
+template<typename T, typename bpred>
+void 
+seq_lookup (const T *table, octave_idx_type offset, octave_idx_type size,
+            const T *vals, octave_idx_type nvals,
+            octave_idx_type *idx, bpred comp)
+{
+  const T *begin = table + offset;
+
+  if (size == 0)
+    // the trivial case of empty table
+    std::fill_n (idx, nvals, offset);
+  else
+    {
+      const T *vcur = vals;
+      const T *vend = vals + nvals;
+
+      const T *cur = begin;
+      const T *end = begin + size;
+
+      while (vcur < vend)
+        {
+          // determine the enclosing interval for next value, trying
+          // ++cur as a special case;
+          if (cur == end || comp (*vcur, *cur))
+            cur = std::upper_bound (begin, cur, *vcur, comp);
+          else
+            {
+              ++cur;
+              if (cur < end && ! comp (*vcur, *cur))
+                cur = std::upper_bound (cur + 1, end, *vcur, comp);
+            }
+
+          // store index of the current interval.
+          *(idx++) = (cur - table);
+          ++vcur;
+
+          // find first value not in current subrange
+          const T *vnew;
+          if (cur < end)
+            if (cur > begin)
+              // inner interval
+              vnew = std::find_if (vcur, vend,
+                                   chk_out_range (*(cur-1), *cur, comp));
+
+            else
+              // special case: lowermost range (-Inf, min) 
+              vnew = std::find_if (vcur, vend,
+                                   not1 (bind2nd (comp, *cur)));
+          else
+            // special case: uppermost range [max, Inf)
+            vnew = std::find_if (vcur, vend,
+                                 bind2nd (comp, *(cur-1)));
+
+          // store index of the current interval.
+          idx = std::fill_n (idx, vnew - vcur, cur - table);
+          vcur = vnew;
+
+        }
+    }
+}
+
+// overload using < operator
+template<typename T, typename bpred>
+void 
+seq_lookup (const T *table, octave_idx_type offset, octave_idx_type size,
+            const T *vals, octave_idx_type nvals,
+            octave_idx_type *idx)
+{
+  seq_lookup (table, offset, size, vals, nvals, idx, std::less<T>());
+}
+
+// helper functions - determine whether an array is descending
+template<typename T>
+bool 
+is_descending (const T *table, octave_idx_type size)
+{
+  return size > 1 && table[size-1] < table[0];
+}
+
+template<typename T, typename bpred>
+bool 
+is_descending (const T *table, octave_idx_type size,
+                    bpred comp)
+{
+  return size > 1 && comp (table[size-1], table[0]);
+}
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/

--- a/scripts/ChangeLog	Fri Mar 28 14:52:27 2008 -0400
+++ b/scripts/ChangeLog	Fri Mar 28 15:53:09 2008 -0400
@@ -1,3 +1,11 @@
+2008-03-27  Jaroslav Hajek  <highegg@gmail.com>
+
+	* general/lookup.m: Remove (lookup moved to DLD-FUNCTIONS).
+	* general/Makefile.in (SOURCES): Delete lookup.m from the list.
+	* general/interp1.m, general/interp2.m, general/interpn.m,
+	polynomial/ppval.m: Fix buggy lookup calls.
+	* general/interp1.m: New test.
+
 2008-03-28  Thomas Weber  <thomas.weber.mail@gmail.com>
 
 	* miscellaneous/tempdir.m: Use correct function name in texinfo

--- a/scripts/general/Makefile.in	Fri Mar 28 14:52:27 2008 -0400
+++ b/scripts/general/Makefile.in	Fri Mar 28 15:53:09 2008 -0400
@@ -41,7 +41,7 @@
   interp1.m interp2.m interp3.m interpn.m interpft.m \
   is_duplicate_entry.m isa.m isdefinite.m isdir.m isequal.m \
   isequalwithequalnans.m isscalar.m issquare.m issymmetric.m \
-  isvector.m logical.m logspace.m lookup.m mod.m nargchk.m \
+  isvector.m logical.m logspace.m mod.m nargchk.m \
   nargoutchk.m nextpow2.m nthroot.m num2str.m perror.m pol2cart.m \
   polyarea.m postpad.m prepad.m quadl.m randperm.m rat.m rem.m \
   repmat.m rot90.m rotdim.m shift.m shiftdim.m sortrows.m \

--- a/scripts/general/interp1.m	Fri Mar 28 14:52:27 2008 -0400
+++ b/scripts/general/interp1.m	Fri Mar 28 15:53:09 2008 -0400
@@ -200,9 +200,7 @@
       yi = mkpp (x, [dy./dx, y(1:end-1)], szy(2:end));
     else
       ## find the interval containing the test point
-      idx = lookup (x(2:nx-1), xi)+1;
-				# 2:n-1 so that anything beyond the ends
-				# gets dumped into an interval
+      idx = lookup (x, xi, "lr");
       ## use the endpoints of the interval to define a line
       s = (xi - x(idx))./dx(idx);
       yi(range,:) = s(:,ones(1,nc)).*dy(idx,:) + y(idx,:);
@@ -243,7 +241,7 @@
     if (nx < 4 || ny < 4)
       error ("interp1: table too short");
     endif
-    idx = lookup (x(3:nx-2), xi) + 1;
+    idx = lookup (x(2:nx-1), xi, "lr");
 
     ## Construct cubic equations for each interval using divided
     ## differences (computation of c and d don't use divided differences
@@ -570,3 +568,5 @@
 %!assert (interp1(1:2:8,1:2:8,1.4,"*cubic"),1.4);
 %!error interp1(1:2,1:2,1, "*spline");
 %!assert (interp1(1:2:6,1:2:6,1.4,"*spline"),1.4);
+
+%!assert (interp1([3,2,1],[3,2,2],2.5),2.5)

--- a/scripts/general/interp2.m	Fri Mar 28 14:52:27 2008 -0400
+++ b/scripts/general/interp2.m	Fri Mar 28 15:53:09 2008 -0400
@@ -185,8 +185,8 @@
     XI = reshape (XI, 1, prod (shape));
     YI = reshape (YI, 1, prod (shape));
 
-    xidx = lookup (X(1, 2:end-1), XI) + 1;
-    yidx = lookup (Y(2:end-1, 1), YI) + 1;
+    xidx = lookup (X(1, :), XI, "lr");
+    yidx = lookup (Y(:, 1), YI, "lr");
 
     if (strcmp (method, "linear"))
       ## each quad satisfies the equation z(x,y)=a+b*x+c*y+d*xy

--- a/scripts/general/interpn.m	Fri Mar 28 14:52:27 2008 -0400
+++ b/scripts/general/interpn.m	Fri Mar 28 15:53:09 2008 -0400
@@ -144,7 +144,7 @@
     yidx = cell (1, nd);
     for i = 1 : nd
       y{i} = y{i}(:);
-      yidx{i} = lookup (x{i}(2:end-1), y{i}) + 1;
+      yidx{i} = lookup (x{i}, y{i}, "lr");
     endfor
     idx = cell (1,nd);
     for i = 1 : nd

--- a/scripts/general/lookup.m	Fri Mar 28 14:52:27 2008 -0400
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,100 +0,0 @@
-## Copyright (C) 2000, 2006, 2007 Paul Kienzle
-##
-## 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{idx} =} lookup (@var{table}, @var{y})
-## Lookup values in a sorted table.  Usually used as a prelude to
-## interpolation.
-##
-## If table is strictly increasing and @code{idx = lookup (table, y)}, then
-## @code{table(idx(i)) <= y(i) < table(idx(i+1))} for all @code{y(i)}
-## within the table.  If @code{y(i)} is before the table, then 
-## @code{idx(i)} is 0. If @code{y(i)} is after the table then
-## @code{idx(i)} is @code{table(n)}.
-##
-## If the table is strictly decreasing, then the tests are reversed.
-## There are no guarantees for tables which are non-monotonic or are not
-## strictly monotonic.
-##
-## To get an index value which lies within an interval of the table,
-## use:
-##
-## @example
-## idx = lookup (table(2:length(table)-1), y) + 1
-## @end example
-##
-## @noindent
-## This expression puts values before the table into the first
-## interval, and values after the table into the last interval.
-## @end deftypefn
-
-## Changed from binary search to sort.
-## Thanks to Kai Habel <kai.habel@gmx.de> for the suggestion.
-
-## TODO: sort-based lookup is significantly slower given a large table
-## TODO: and small lookup vector.  This shouldn't be a problem since
-## TODO: interpolation (the reason for the table lookup in the first
-## TODO: place) usually involves subsampling of an existing table.  The
-## TODO: other use of interpolation, looking up values one at a time, is
-## TODO: unfortunately significantly slower for large tables.  
-## TODO:    sort is order O((lt+lx)*log(lt+lx)) 
-## TODO:    search is order O(lx*log(lt))
-## TODO: Clearly, search is asymptotically better than sort, but sort
-## TODO: is compiled and search is not.  Could support both, or recode
-## TODO: search in C++, or assume things are good enough as they stand.
-
-function idx = lookup (table, xi)
-  if (nargin == 2)
-    if (isempty (table))
-      idx = zeros (size (xi));
-    elseif (isvector (table))
-      [nr, nc] = size (xi);
-      lt = length (table);
-      if (table(1) > table(lt))
-	## decreasing table
-	[v, p] = sort ([xi(:); table(:)]);
-	idx(p) = cumsum (p > nr*nc);
-	idx = lt - idx(1:nr*nc);
-      else
-	## increasing table
-	[v, p] = sort ([table(:); xi(:) ]);
-	idx(p) = cumsum (p <= lt);
-	idx = idx(lt+1:lt+nr*nc);
-      endif
-      idx = reshape (idx, nr, nc);
-    else
-      error ("lookup: table must be a vector");
-    endif
-  else
-    print_usage ();
-  endif
-endfunction
-  
-%!assert (lookup(1:3, 0.5), 0)     # value before table
-%!assert (lookup(1:3, 3.5), 3)     # value after table error
-%!assert (lookup(1:3, 1.5), 1)     # value within table error
-%!assert (lookup(1:3, [3,2,1]), [3,2,1])
-%!assert (lookup([1:4]', [1.2, 3.5]'), [1, 3]');
-%!assert (lookup([1:4], [1.2, 3.5]'), [1, 3]');
-%!assert (lookup([1:4]', [1.2, 3.5]), [1, 3]);
-%!assert (lookup([1:4], [1.2, 3.5]), [1, 3]);
-%!assert (lookup(1:3, [3, 2, 1]), [3, 2, 1]);
-%!assert (lookup([3:-1:1], [3.5, 3, 1.2, 2.5, 2.5]), [0, 1, 2, 1, 1])
-%!assert (isempty(lookup([1:3], [])))
-%!assert (isempty(lookup([1:3]', [])))
-%!assert (lookup(1:3, [1, 2; 3, 0.5]), [1, 2; 3, 0]);

--- a/scripts/polynomial/ppval.m	Fri Mar 28 14:52:27 2008 -0400
+++ b/scripts/polynomial/ppval.m	Fri Mar 28 15:53:09 2008 -0400
@@ -39,7 +39,7 @@
     transposed = (columns (xi) == 1);
     xi = xi(:);
     xn = length (xi);
-    idx = lookup (pp.x(2:pp.n), xi) + 1;
+    idx = lookup (pp.x, xi, "lr");
     dx = (xi - pp.x(idx)).';
     dx = reshape (dx(ones(1,prod(pp.d)),:),[pp.d,xn]);
     c = reshape (pp.P(:,1), pp.n, prod (pp.d));

--- a/src/ChangeLog	Fri Mar 28 14:52:27 2008 -0400
+++ b/src/ChangeLog	Fri Mar 28 15:53:09 2008 -0400
@@ -1,3 +1,8 @@
+2008-03-28  Jaroslav Hajek  <highegg@gmail.com>
+
+	* DLD-FUNCTIONS/lookup.cc: New file.
+	* Makefile.in (DLD_XSRC): Add it to the list.
+
 2008-03-28  David Bateman  <dbateman@free.fr>
 
 	* ov-complex.cc (SCALAR_MAPPER, CD_SCALAR_MAPPER): New macro for

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/DLD-FUNCTIONS/lookup.cc	Fri Mar 28 15:53:09 2008 -0400
@@ -0,0 +1,308 @@
+/*
+
+Copyright (C) 2008 VZLU Prague a.s., Czech Republic
+
+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/>.
+
+*/
+
+// Author: Jaroslav Hajek <highegg@gmail.com>
+
+#include <cctype>
+#include <functional>
+#include <algorithm>
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include "dNDArray.h"
+#include "CNDArray.h"
+#include "oct-lookup.h"
+
+#include "Cell.h"
+#include "defun-dld.h"
+#include "error.h"
+#include "gripes.h"
+#include "oct-obj.h"
+#include "ov.h"
+
+static
+bool
+contains_char (const std::string& str, char c)
+{
+  return (str.find (c) != std::string::npos 
+	  || str.find (std::toupper (c)) != std::string::npos);
+}
+
+// FIXME -- remove these one once octave_value supports octave_idx_type.
+static octave_value&
+assign (octave_value& ov, octave_idx_type idx)
+{
+  double tmp = idx;
+  ov = tmp;
+  return ov;
+}
+
+static octave_value&
+assign (octave_value& ov, const ArrayN<octave_idx_type>& ida)
+{
+  NDArray tmp (ida.dims ());
+  for (int i = 0; i < ida.numel (); i++)
+    tmp(i) = ida(i);
+  ov = tmp;
+  return ov;
+}
+
+// normal ascending comparator
+static bool
+ov_str_less (const octave_value& a, const octave_value& b)
+{
+  return a.string_value () < b.string_value ();
+}
+
+// normal descending comparator
+static bool
+ov_str_greater (const octave_value& a, const octave_value& b)
+{
+  return a.string_value () > b.string_value ();
+}
+
+// case-insensitive character comparison functors
+struct icmp_char_lt : public std::binary_function<char, char, bool>
+{
+  bool operator () (char x, char y) const
+    { return std::toupper (x) < std::toupper (y); }
+};
+
+struct icmp_char_gt : public std::binary_function<char, char, bool>
+{
+  bool operator () (char x, char y) const
+    { return std::toupper (x) > std::toupper (y); }
+};
+
+// case-insensitive ascending comparator
+static bool
+ov_stri_less (const octave_value& a, const octave_value& b)
+{
+  std::string as = a.string_value ();
+  std::string bs = b.string_value ();
+
+  return std::lexicographical_compare (as.begin (), as.end (), 
+                                       bs.begin (), bs.end (),
+                                       icmp_char_lt());
+}
+
+// case-insensitive descending comparator
+static bool
+ov_stri_greater (const octave_value& a, const octave_value& b)
+{
+  std::string as = a.string_value ();
+  std::string bs = b.string_value ();
+
+  return std::lexicographical_compare (as.begin (), as.end (), 
+                                       bs.begin (), bs.end (),
+                                       icmp_char_gt());
+}
+
+DEFUN_DLD (lookup, args, ,
+  "-*- texinfo -*-\n\
+@deftypefn {Function File} {@var{idx} =} lookup (@var{table}, @var{y}, @var{opt})\n\
+Lookup values in a sorted table.  Usually used as a prelude to\n\
+interpolation.\n\
+\n\
+If table is strictly increasing and @code{idx = lookup (table, y)}, then\n\
+@code{table(idx(i)) <= y(i) < table(idx(i+1))} for all @code{y(i)}\n\
+within the table.  If @code{y(i)} is before the table, then\n\
+@code{idx(i)} is 0. If @code{y(i)} is after the table then\n\
+@code{idx(i)} is @code{table(n)}.\n\
+\n\
+If the table is strictly decreasing, then the tests are reversed.\n\
+There are no guarantees for tables which are non-monotonic or are not\n\
+strictly monotonic.\n\
+\n\
+@var{table} and @var{y} can also be a cell array of strings\n\
+(or @var{y} can be a single string). In this case, string lookup\n\
+is performed using lexicographical comparison.\n\
+If @var{opts} is specified, it shall be a string with letters indicating\n\
+additional options.\n\
+\n\
+For numeric lookup, 'l' in @var{opts} indicates that\n\
+the leftmost subinterval shall be extended to infinity (i.e. all indices\n\
+at least 1), and 'r' indicates that the rightmost subinterval shall be\n\
+extended to infinity (i.e. all indices at most n-1).\n\
+\n\
+For string lookup, 'i' indicates case-insensitive comparison.\n\
+@end deftypefn") 
+{
+  octave_value_list retval;
+
+  int nargin = args.length ();
+
+  if (nargin < 2 || nargin > 3 || (nargin == 3 && ! args(2).is_string ()))
+    {
+      print_usage ();
+      return retval;
+    }
+
+  octave_value argtable = args(0), argy = args(1);
+  if (argtable.ndims () > 2 || (argtable.columns () > 1 && argtable.rows () > 1))
+    warning ("lookup: table is not a vector");
+
+  bool num_case = argtable.is_numeric_type () && argy.is_numeric_type ();
+  bool str_case = argtable.is_cell () && (argy.is_cell () || argy.is_string ());
+
+  if (num_case) 
+    {
+      bool left_inf = false;
+      bool right_inf = false;
+
+      if (nargin == 3)
+        {
+          std::string opt = args(2).string_value ();
+          left_inf = contains_char (opt, 'l');
+          right_inf = contains_char (opt, 'r');
+        }
+
+      // in the case of a complex array, absolute values will be used for compatibility
+      // (though it's not too meaningful).
+
+      NDArray table = (argtable.is_complex_type ()) 
+        ? argtable.complex_array_value ().abs ()
+        : argtable.array_value ();
+
+      NDArray y = (argy.is_complex_type ()) 
+        ? argy.complex_array_value ().abs ()
+        : argy.array_value ();
+
+      ArrayN<octave_idx_type> idx (y.dims ());
+
+      // determine whether the array is descending. 
+      bool desc = is_descending (table.data (), table.length ());
+      octave_idx_type offset = left_inf ? 1 : 0;
+      octave_idx_type size = table.length () - offset - (right_inf ? 1 : 0);
+      if (size < 0) 
+        size = 0;
+
+      if (desc)
+        seq_lookup (table.data (), offset, size, 
+                    y.data (), y.length (), idx.fortran_vec (),
+                    std::greater<double> ());
+      else
+        seq_lookup (table.data (), offset, size, 
+                    y.data (), y.length (), idx.fortran_vec (),
+                    std::less<double> ());
+
+
+      //retval(0) = idx;
+      assign (retval(0), idx);
+    }
+  else if (str_case)
+    {
+      Cell table = argtable.cell_value ();
+      
+      bool (*ov_str_comp) (const octave_value&, const octave_value&);
+
+      bool icase = false, desc;
+
+      // check for case-insensitive option
+      if (nargin == 3)
+        {
+          std::string opt = args(2).string_value ();
+          icase = contains_char (opt, 'i');
+        }
+
+      // pick the correct comparator
+      if (icase)
+        {
+          if (is_descending (table.data (), table.length (), ov_stri_less))
+            ov_str_comp = ov_stri_greater;
+          else
+            ov_str_comp = ov_stri_less;
+        }
+      else
+        {
+          if (is_descending (table.data (), table.length (), ov_str_less))
+            ov_str_comp = ov_str_greater;
+          else
+            ov_str_comp = ov_str_less;
+        }
+
+
+      // query just the first cell to verify it's a string
+      if (table(0).is_string ())
+        {
+          if (argy.is_cell ())
+            {
+              Cell y = argy.cell_value ();
+              ArrayN<octave_idx_type> idx (y.dims ());
+
+
+
+              for (int i = 0; i < y.numel (); i++)
+                  idx(i) = bin_lookup (table.data (), table.length (), y(i), 
+                                       std::ptr_fun (ov_str_comp));
+
+              //retval(0) = idx;
+              assign (retval(0), idx);
+            }
+          else
+            {
+              octave_idx_type idx;
+
+              idx = bin_lookup (table.data (), table.length (), argy, 
+                                std::ptr_fun (ov_str_comp));
+
+              //retval(0) = idx;
+              assign (retval(0), idx);
+            }
+        }
+      else
+        error("lookup: table is not a cell array of strings.");
+    }
+  else
+    print_usage ();
+
+  return retval;
+
+}  
+
+/*
+%!assert (real(lookup(1:3, 0.5)), 0)     # value before table
+%!assert (real(lookup(1:3, 3.5)), 3)     # value after table error
+%!assert (real(lookup(1:3, 1.5)), 1)     # value within table error
+%!assert (real(lookup(1:3, [3,2,1])), [3,2,1])
+%!assert (real(lookup([1:4]', [1.2, 3.5]')), [1, 3]');
+%!assert (real(lookup([1:4], [1.2, 3.5]')), [1, 3]');
+%!assert (real(lookup([1:4]', [1.2, 3.5])), [1, 3]);
+%!assert (real(lookup([1:4], [1.2, 3.5])), [1, 3]);
+%!assert (real(lookup(1:3, [3, 2, 1])), [3, 2, 1]);
+%!assert (real(lookup([3:-1:1], [3.5, 3, 1.2, 2.5, 2.5])), [0, 1, 2, 1, 1])
+%!assert (isempty(lookup([1:3], [])))
+%!assert (isempty(lookup([1:3]', [])))
+%!assert (real(lookup(1:3, [1, 2; 3, 0.5])), [1, 2; 3, 0]);
+%!
+%!assert (real(lookup({"apple","lemon","orange"}, {"banana","kiwi"; "ananas","mango"})), [1,1;0,2])
+%!assert (real(lookup({"apple","lemon","orange"}, "potato")), 3)
+%!assert (real(lookup({"orange","lemon","apple"}, "potato")), 0)
+*/
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/

--- a/src/Makefile.in	Fri Mar 28 14:52:27 2008 -0400
+++ b/src/Makefile.in	Fri Mar 28 15:53:09 2008 -0400
@@ -67,7 +67,7 @@
 	dasrt.cc dassl.cc det.cc dispatch.cc dlmread.cc dmperm.cc eig.cc \
 	expm.cc fft.cc fft2.cc fftn.cc fftw.cc filter.cc find.cc fsolve.cc \
 	gammainc.cc gcd.cc getgrent.cc getpwent.cc getrusage.cc \
-	givens.cc hess.cc hex2num.cc inv.cc kron.cc lsode.cc \
+	givens.cc hess.cc hex2num.cc inv.cc kron.cc lookup.cc lsode.cc \
 	lu.cc luinc.cc matrix_type.cc max.cc md5sum.cc pinv.cc qr.cc \
 	quad.cc qz.cc rand.cc regexp.cc schur.cc sparse.cc \
 	spparms.cc sqrtm.cc svd.cc syl.cc symrcm.cc symbfact.cc \