# HG changeset patch
# User jwe
# Date 1069276999 0
# Node ID 79fe96966ca017751c22fac5b68ced95d8c31cf3
# Parent  d95dd6ae89fdc2dd7a6d83ab1e4578872e68ae5b
[project @ 2003-11-19 21:22:39 by jwe]

diff --git a/liboctave/CNDArray.cc b/liboctave/CNDArray.cc
--- a/liboctave/CNDArray.cc
+++ b/liboctave/CNDArray.cc
@@ -120,25 +120,15 @@
   MX_ND_COMPLEX_OP_REDUCTION (+= elem (iter_idx), Complex (0, 0));
 }
 
-Matrix 
+NDArray
 ComplexNDArray::abs (void) const
 {
-  Matrix retval;
-
-  if (dims () . length () == 2)
-    {
-      int nr = rows ();
-      int nc = cols ();
+  NDArray retval (dims ());
 
-      retval.resize (nr, nc);
-      
-      for (int j = 0; j < nc; j++)
-	for (int i = 0; i < nr; i++)
-	  retval (i, j) = ::abs (elem (i, j));
-    }
-  else
-    (*current_liboctave_error_handler)
-      ("abs is not yet implemented for N-d arrays");
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    retval(i) = ::abs (elem (i));
       
   return retval;
 }
diff --git a/liboctave/CNDArray.h b/liboctave/CNDArray.h
--- a/liboctave/CNDArray.h
+++ b/liboctave/CNDArray.h
@@ -83,7 +83,7 @@
   ComplexNDArray sum (int dim = -1) const;
   ComplexNDArray sumsq (int dim = -1) const;
 
-  Matrix abs (void) const;
+  NDArray abs (void) const;
 
   ComplexMatrix matrix_value (void) const;
 
@@ -108,7 +108,7 @@
 
 private:
 
-  ComplexNDArray (Complex *d, dim_vector& dv)
+  ComplexNDArray (Complex *d, const dim_vector& dv)
     : MArrayN<Complex> (d, dv) { }
 };
 
diff --git a/liboctave/ChangeLog b/liboctave/ChangeLog
--- a/liboctave/ChangeLog
+++ b/liboctave/ChangeLog
@@ -1,3 +1,21 @@
+2003-11-19  John W. Eaton  <jwe@bevo.che.wisc.edu>
+
+	* dNDArray.cc (NDArray::all_elements_are_int_or_inf_or_nan,
+	NDArray::any_element_is_inf_or_nan, NDArray::too_large_for_float):
+	New functions.
+	* dNDArray.h: Provide decls.
+
+	* dMatrix.h (Matrix::any_element_is_negative,
+	Matrix::any_element_is_inf_or_nan, Matrix::too_large_for_float,
+	Matrix::all_elements_are_int_or_inf_or_nan, Matrix::all_integers):
+	Simplify.
+
+	* dNDArray.cc (NDArray::abs): Make it work for N-d arrays.
+	* CNDArray.cc (ComplexNDArray::abs): Likewise.
+
+	* dNDArray.cc (real, imag): New functions.
+	* dNDArray.h: Provide decls.
+
 2003-11-18  John W. Eaton  <jwe@bevo.che.wisc.edu>
 
 	* Makefile.in (TEMPLATE_SRC): Move MArrayN.cc here from MATRIX_SRC.
diff --git a/liboctave/dMatrix.cc b/liboctave/dMatrix.cc
--- a/liboctave/dMatrix.cc
+++ b/liboctave/dMatrix.cc
@@ -2033,22 +2033,19 @@
 bool
 Matrix::any_element_is_negative (bool neg_zero) const
 {
-  int nr = rows ();
-  int nc = cols ();
+  int nel = nelem ();
 
   if (neg_zero)
     {
-      for (int j = 0; j < nc; j++)
-	for (int i = 0; i < nr; i++)
-	  if (lo_ieee_signbit (elem (i, j)))
-	    return true;
+      for (int i = 0; i < nel; i++)
+	if (lo_ieee_signbit (elem (i)))
+	  return true;
     }
   else
     {
-      for (int j = 0; j < nc; j++)
-	for (int i = 0; i < nr; i++)
-	  if (elem (i, j) < 0)
-	    return true;
+      for (int i = 0; i < nel; i++)
+	if (elem (i) < 0)
+	  return true;
     }
 
   return false;
@@ -2058,34 +2055,31 @@
 bool
 Matrix::any_element_is_inf_or_nan (void) const
 {
-  int nr = rows ();
-  int nc = cols ();
-
-  for (int j = 0; j < nc; j++)
-    for (int i = 0; i < nr; i++)
-      {
-	double val = elem (i, j);
-	if (xisinf (val) || xisnan (val))
-	  return 1;
-      }
-  return 0;
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+      if (xisinf (val) || xisnan (val))
+	return true;
+    }
+
+  return false;
 }
 
 bool
 Matrix::all_elements_are_int_or_inf_or_nan (void) const
 {
-  int nr = rows ();
-  int nc = cols ();
-
-  for (int j = 0; j < nc; j++)
-    for (int i = 0; i < nr; i++)
-      {
-	double val = elem (i, j);
-	if (xisnan (val) || D_NINT (val) == val)
-	  continue;
-	else
-	  return false;
-      }
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+      if (xisnan (val) || D_NINT (val) == val)
+	continue;
+      else
+	return false;
+    }
 
   return true;
 }
@@ -2096,31 +2090,29 @@
 bool
 Matrix::all_integers (double& max_val, double& min_val) const
 {
-  int nr = rows ();
-  int nc = cols ();
-
-  if (nr > 0 && nc > 0)
+  int nel = nelem ();
+
+  if (nel > 0)
     {
-      max_val = elem (0, 0);
-      min_val = elem (0, 0);
+      max_val = elem (0);
+      min_val = elem (0);
     }
   else
     return false;
 
-  for (int j = 0; j < nc; j++)
-    for (int i = 0; i < nr; i++)
-      {
-	double val = elem (i, j);
-
-	if (val > max_val)
-	  max_val = val;
-
-	if (val < min_val)
-	  min_val = val;
-
-	if (D_NINT (val) != val)
-	  return false;
-      }
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+
+      if (val > max_val)
+	max_val = val;
+
+      if (val < min_val)
+	min_val = val;
+
+      if (D_NINT (val) != val)
+	return false;
+    }
 
   return true;
 }
@@ -2128,17 +2120,15 @@
 bool
 Matrix::too_large_for_float (void) const
 {
-  int nr = rows ();
-  int nc = cols ();
-
-  for (int j = 0; j < nc; j++)
-    for (int i = 0; i < nr; i++)
-      {
-	double val = elem (i, j);
-
-	if (val > FLT_MAX || val < FLT_MIN)
-	  return true;
-      }
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+
+      if (val > FLT_MAX || val < FLT_MIN)
+	return true;
+    }
 
   return false;
 }
diff --git a/liboctave/dNDArray.cc b/liboctave/dNDArray.cc
--- a/liboctave/dNDArray.cc
+++ b/liboctave/dNDArray.cc
@@ -29,11 +29,14 @@
 #include <config.h>
 #endif
 
+#include <cfloat>
+
 #include "Array-util.h"
 #include "dNDArray.h"
 #include "mx-base.h"
 #include "lo-error.h"
 #include "lo-ieee.h"
+#include "lo-mappers.h"
 
 NDArray::NDArray (const boolNDArray& a)
   : MArrayN<double> (a.dims ())
@@ -62,9 +65,111 @@
   return b;
 }
 
+bool
+NDArray::any_element_is_negative (bool neg_zero) const
+{
+  int nel = nelem ();
+
+  if (neg_zero)
+    {
+      for (int i = 0; i < nel; i++)
+	if (lo_ieee_signbit (elem (i)))
+	  return true;
+    }
+  else
+    {
+      for (int i = 0; i < nel; i++)
+	if (elem (i) < 0)
+	  return true;
+    }
+
+  return false;
+}
+
+
+bool
+NDArray::any_element_is_inf_or_nan (void) const
+{
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+      if (xisinf (val) || xisnan (val))
+	return true;
+    }
+
+  return false;
+}
+
+bool
+NDArray::all_elements_are_int_or_inf_or_nan (void) const
+{
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+      if (xisnan (val) || D_NINT (val) == val)
+	continue;
+      else
+	return false;
+    }
+
+  return true;
+}
+
+// Return nonzero if any element of M is not an integer.  Also extract
+// the largest and smallest values and return them in MAX_VAL and MIN_VAL.
+
+bool
+NDArray::all_integers (double& max_val, double& min_val) const
+{
+  int nel = nelem ();
+
+  if (nel > 0)
+    {
+      max_val = elem (0);
+      min_val = elem (0);
+    }
+  else
+    return false;
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+
+      if (val > max_val)
+	max_val = val;
+
+      if (val < min_val)
+	min_val = val;
+
+      if (D_NINT (val) != val)
+	return false;
+    }
+
+  return true;
+}
+
+bool
+NDArray::too_large_for_float (void) const
+{
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      double val = elem (i);
+
+      if (val > FLT_MAX || val < FLT_MIN)
+	return true;
+    }
+
+  return false;
+}
+
 // XXX FIXME XXX -- this is not quite the right thing.
 
-
 boolNDArray
 NDArray::all (int dim) const
 {
@@ -107,30 +212,39 @@
   MX_ND_REAL_OP_REDUCTION (+= elem (iter_idx), 0);
 }
 
-Matrix
+NDArray
+real (const ComplexNDArray& a)
+{
+  int a_len = a.length ();
+  NDArray retval;
+  if (a_len > 0)
+    retval = NDArray (mx_inline_real_dup (a.data (), a_len), a.dims ());
+  return retval;
+}
+
+NDArray
+imag (const ComplexNDArray& a)
+{
+  int a_len = a.length ();
+  NDArray retval;
+  if (a_len > 0)
+    retval = NDArray (mx_inline_imag_dup (a.data (), a_len), a.dims ());
+  return retval;
+}
+
+NDArray
 NDArray::abs (void) const
 {
-  Matrix retval;
-
-  if (dims () . length () == 2)
-    {
-      int nr = rows ();
-      int nc = cols ();
+  NDArray retval (dims ());
 
-      retval.resize (nr, nc);
-      
-      for (int j = 0; j < nc; j++)
-	for (int i = 0; i < nr; i++)
-	  retval(i,j) = fabs (elem (i, j));
-    }
-  else
-    (*current_liboctave_error_handler)
-      ("abs is not yet implemented for N-d arrays");
+  int nel = nelem ();
+
+  for (int i = 0; i < nel; i++)
+    retval(i) = fabs (elem (i));
 
   return retval;
 }
 
-
 Matrix
 NDArray::matrix_value (void) const
 {
@@ -172,56 +286,6 @@
   return ::compute_index (ra_idx, dimensions);
 }
 
-bool
-NDArray::any_element_is_negative (bool neg_zero) const
-{
-  int n = length (); 
-  if (neg_zero)
-    {
-      for (int i = 0; i < n; i++)
-	if (lo_ieee_signbit (Array<double>::elem (i)))
-	  return true;
-    }
-  else
-    {
-      for (int i = 0; i < n; i++)
-	if (Array<double>::elem (i) < 0)
-	  return true;
-    }
- 
- return false;
-}
-
-bool
-NDArray::all_integers (double& max_val, double& min_val) const
-{
-  int n = length ();
-
-  if (n > 0)
-    {
-      max_val = Array<double>::elem (0);
-      min_val = Array<double>::elem (0);
-    }
-  else 
-    return false;
-
-  for (int i = 0; i < n; i++)
-    {
-      double val = Array<double>::elem (0);
-      
-      if (val > max_val)
-	max_val = val;
-
-      if (val < min_val)
-	min_val = val;
-
-      if (D_NINT (val) != val)
-	return false;
-    }
-
-  return true;
-}
-
 NDS_CMP_OPS(NDArray, , double, )
 NDS_BOOL_OPS(NDArray, double, 0.0)
 
diff --git a/liboctave/dNDArray.h b/liboctave/dNDArray.h
--- a/liboctave/dNDArray.h
+++ b/liboctave/dNDArray.h
@@ -70,6 +70,12 @@
 
   boolNDArray operator ! (void) const;
 
+  bool any_element_is_negative (bool = false) const;
+  bool any_element_is_inf_or_nan (void) const;
+  bool all_elements_are_int_or_inf_or_nan (void) const;
+  bool all_integers (double& max_val, double& min_val) const;
+  bool too_large_for_float (void) const;
+
   // XXX FIXME XXX -- this is not quite the right thing.
 
   boolNDArray all (int dim = -1) const;
@@ -81,7 +87,10 @@
   NDArray sum (int dim = -1) const;  
   NDArray sumsq (int dim = -1) const;
  
-  Matrix abs (void) const;
+  NDArray abs (void) const;
+
+  friend NDArray real (const ComplexNDArray& a);
+  friend NDArray imag (const ComplexNDArray& a);
 
   Matrix matrix_value (void) const;
 
@@ -101,12 +110,9 @@
 
   static double resize_fill_value (void) { return 0; }
 
-  bool any_element_is_negative (bool = false) const;
-  bool all_integers (double& max_val, double& min_val) const;
-
 private:
 
-  NDArray (double *d, dim_vector& dv) : MArrayN<double> (d, dv) { }
+  NDArray (double *d, const dim_vector& dv) : MArrayN<double> (d, dv) { }
 };
 
 NDS_CMP_OP_DECLS (NDArray, double)
diff --git a/src/DLD-FUNCTIONS/lsode.cc b/src/DLD-FUNCTIONS/lsode.cc
--- a/src/DLD-FUNCTIONS/lsode.cc
+++ b/src/DLD-FUNCTIONS/lsode.cc
@@ -63,8 +63,6 @@
 {
   ColumnVector retval;
 
-  int nstates = x.capacity ();
-
   octave_value_list args;
   args(1) = t;
   args(0) = x;
@@ -104,8 +102,6 @@
 {
   Matrix retval;
 
-  int nstates = x.capacity ();
-
   octave_value_list args;
   args(1) = t;
   args(0) = x;
diff --git a/src/ls-hdf5.cc b/src/ls-hdf5.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-hdf5.cc
@@ -0,0 +1,1241 @@
+/*
+
+Copyright (C) 1996, 1997 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+// Author: Steven G. Johnson <stevenj@alum.mit.edu>
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#if defined (HAVE_HDF5)
+
+#include <cfloat>
+#include <cstring>
+#include <cctype>
+
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#include <memory>
+#include <string>
+
+#include <hdf5.h>
+
+#include "byte-swap.h"
+#include "data-conv.h"
+#include "file-ops.h"
+#include "glob-match.h"
+#include "lo-mappers.h"
+#include "lo-sstream.h"
+#include "mach-info.h"
+#include "oct-env.h"
+#include "oct-time.h"
+#include "quit.h"
+#include "str-vec.h"
+
+#include "Cell.h"
+#include "defun.h"
+#include "error.h"
+#include "gripes.h"
+#include "load-save.h"
+#include "oct-obj.h"
+#include "oct-map.h"
+#include "ov-cell.h"
+#include "pager.h"
+#include "pt-exp.h"
+#include "symtab.h"
+#include "sysdep.h"
+#include "unwind-prot.h"
+#include "utils.h"
+#include "variables.h"
+#include "version.h"
+#include "dMatrix.h"
+
+#include "ls-utils.h"
+#include "ls-hdf5.h"
+
+// this is only used for HDF5 import
+// try to convert s into a valid identifier, replacing invalid chars with "_":
+
+static std::string
+make_valid_identifier (const std::string& nm)
+{
+  std::string retval;
+
+  size_t nm_len = nm.length ();
+
+  if (nm_len > 0)
+    {
+      if (! isalpha (nm[0]))
+	retval += '_';
+
+      for (size_t i = 0; i < nm_len; i++)
+	{
+	  char c = nm[i];
+	  retval += (isalnum (c) || c == '_') ? c : '_';
+	}
+    }
+
+  return retval;
+}
+
+static bool
+ident_is_all_digits (const std::string& id)
+{
+  bool retval = true;
+
+  size_t len = 0;
+
+  for (size_t i = 0; i < len; i++)
+    {
+      if (! isdigit (id[i]))
+	{
+	  retval = false;
+	  break;
+	}
+    }
+
+  return retval;
+}
+
+// Define this to 1 if/when HDF5 supports automatic conversion between
+// integer and floating-point binary data:
+#define HAVE_HDF5_INT2FLOAT_CONVERSIONS 0
+
+// Given two compound types t1 and t2, determine whether they 
+// are compatible for reading/writing.  This function only
+// works for non-nested types composed of simple elements (ints, floats...),
+// which is all we need it for
+
+bool
+hdf5_types_compatible (hid_t t1, hid_t t2)
+{
+  int n;
+  if ((n = H5Tget_nmembers (t1)) != H5Tget_nmembers (t2))
+    return false;
+
+  for (int i = 0; i < n; ++i)
+    {
+      hid_t mt1 = H5Tget_member_type (t1, i);
+      hid_t mt2 = H5Tget_member_type (t2, i);
+
+      if (H5Tget_class (mt1) != H5Tget_class (mt2))
+	return false;
+
+      H5Tclose (mt2);
+      H5Tclose (mt1);
+    }
+
+  return true;
+}
+
+// Import a multidimensional (rank >= 3) dataset whose id is data_id, into tc.
+// This works by calling itself recursively, building up lists of lists
+//  of lists ... of 2d matrices.  rank and dims are the rank and dimensions
+//  of the dataset.  type_id is the datatype to read into.  If it is
+//  H5T_NATIVE_DOUBLE, we are reading a real matrix.  Otherwise, type_id
+//  is assumed to be a complex type for reading a complex matrix.
+//
+//  Upon entry, we should have curdim = rank - 1, start = an array
+//  of length rank = all zeros, and count = an array of length rank =
+//  all ones except for the first two dimensions which equal the corresponding
+//  entries in dims[]. 
+//
+//  Note that we process the dimensions in reverse order, reflecting
+//  the fact that Octave is uses column-major (Fortran-order) data while
+//  HDF5 is row-major.  This means that the HDF5 file is read
+//  non-contiguously, but on the other hand means that for a 3d array
+//  we get a list of xy-plane slices, which seems nice.  We could change
+//  this behavior without much trouble; what is the best thing to do?
+//
+//  Returns a positive value upon success.
+
+static herr_t
+hdf5_import_multidim (hid_t data_id, hid_t space_id, hsize_t rank,
+		      const hsize_t *dims, hsize_t curdim,
+		      hssize_t *start, const hsize_t *count,
+		      hid_t type_id, octave_value &tc)
+{
+  herr_t retval = 1;
+
+  if (rank < 3 || curdim < 1 || curdim >= rank)
+    return -1;
+
+  if (curdim == 1)
+    {
+      // import 2d dataset for 1st 2 dims directly as a matrix
+      int nr, nc;    // rows and columns
+      nc = dims[0];  // octave uses column-major & HDF5 uses row-major
+      nr = dims[1];
+
+      hid_t mem_space_id = H5Screate_simple (2, dims, 0);
+
+      if (mem_space_id < 0)
+	return -1;
+
+      if (H5Sselect_all (mem_space_id) < 0)
+	return -1;
+    
+      if (H5Sselect_hyperslab (space_id, H5S_SELECT_SET,
+			       start, 0, count, 0) < 0)
+	{
+	  H5Sclose (mem_space_id);
+	  return -1;
+	}
+    
+      if (type_id == H5T_NATIVE_DOUBLE)
+	{
+	  // real matrix
+	  Matrix m (nr, nc);
+	  double *re = m.fortran_vec ();
+	  if (H5Dread (data_id, type_id, mem_space_id, space_id,
+		       H5P_DEFAULT, (void *) re) < 0)
+	    retval = -1;  // error
+	  else
+	    tc = m;
+	}
+      else
+	{
+	  // assume that we are using complex numbers
+	  // complex matrix
+	  ComplexMatrix m (nr, nc);
+	  Complex *reim = m.fortran_vec ();
+	  if (H5Dread (data_id, type_id, mem_space_id, space_id,
+		       H5P_DEFAULT, (void *) X_CAST (double *, reim)) < 0)
+	    retval = -1;  // error
+	  else
+	    tc = m;
+	}
+    
+      H5Sclose (mem_space_id);
+
+    }
+  else
+    {
+      octave_value_list lst;
+
+      for (hsize_t i = 0; i < dims[curdim]; ++i)
+	{
+	  octave_value slice;
+	  start[curdim] = i;
+	  retval = hdf5_import_multidim (data_id, space_id, rank,
+					 dims, curdim-1, start, count,
+					 type_id, slice);
+	  if (retval < 0)
+	    break;
+	  lst.append (slice);
+	}
+
+      if (retval > 0)
+	tc = octave_value (lst);
+    }
+
+  return retval;
+}
+
+// Return true if loc_id has the attribute named attr_name, and false
+// otherwise.
+
+bool
+hdf5_check_attr (hid_t loc_id, const char *attr_name)
+{
+  bool retval = false;
+
+  // we have to pull some shenanigans here to make sure
+  // HDF5 doesn't print out all sorts of error messages if we
+  // call H5Aopen for a non-existing attribute
+
+  H5E_auto_t err_func;
+  void *err_func_data;
+
+  // turn off error reporting temporarily, but save the error
+  // reporting function:
+
+  H5Eget_auto (&err_func, &err_func_data);
+  H5Eset_auto (0, 0);
+
+  hid_t attr_id = H5Aopen_name (loc_id, attr_name);
+
+  if (attr_id >= 0)
+    {
+      // successful
+      retval = 1;
+      H5Aclose (attr_id);
+    }
+
+  // restore error reporting:
+  H5Eset_auto (err_func, err_func_data);
+
+  return retval;
+}
+
+// The following two subroutines create HDF5 representations of the way
+// we will store Octave complex and range types (pairs and triplets of
+// floating-point numbers, respectively).  NUM_TYPE is the HDF5 numeric
+// type to use for storage (e.g. H5T_NATIVE_DOUBLE to save as 'double').
+// Note that any necessary conversions are handled automatically by HDF5.
+
+static hid_t
+hdf5_make_complex_type (hid_t num_type)
+{
+  hid_t type_id = H5Tcreate (H5T_COMPOUND, sizeof (double) * 2);
+
+  H5Tinsert (type_id, "real", 0 * sizeof (double), num_type);
+  H5Tinsert (type_id, "imag", 1 * sizeof (double), num_type);
+
+  return type_id;
+}
+
+static hid_t
+hdf5_make_range_type (hid_t num_type)
+{
+  hid_t type_id = H5Tcreate (H5T_COMPOUND, sizeof (double) * 3);
+
+  H5Tinsert (type_id, "base", 0 * sizeof (double), num_type);
+  H5Tinsert (type_id, "limit", 1 * sizeof (double), num_type);
+  H5Tinsert (type_id, "increment", 2 * sizeof (double), num_type);
+
+  return type_id;
+}
+
+// Callback data structure for passing data to hdf5_read_next_data, below.
+
+struct
+hdf5_callback_data
+{
+  hdf5_callback_data (void)
+    : name (), global (false), tc (), doc (),
+      complex_type (hdf5_make_complex_type (H5T_NATIVE_DOUBLE)),
+      range_type (hdf5_make_range_type (H5T_NATIVE_DOUBLE)),
+      import (false) { }
+
+  // the following fields are set by hdf5_read_data on successful return:
+
+  // the name of the variable
+  std::string name;
+
+  // whether it is global
+  bool global;
+
+  // the value of the variable, in Octave form
+  octave_value tc;
+
+  // a documentation string (NULL if none)
+  std::string doc;
+
+  // the following fields are input to hdf5_read_data:
+
+  // HDF5 rep's of complex and range type
+  hid_t complex_type, range_type;
+
+  // whether to try extra hard to import "foreign" data
+  bool import;
+};
+
+// This variable, set in read_hdf5_data(), tells whether we are using
+// a version of HDF5 with a buggy H5Giterate (i.e. which neglects to
+// increment the index parameter to the next unread item).
+static bool have_h5giterate_bug = false;
+
+// This function is designed to be passed to H5Giterate, which calls it
+// on each data item in an HDF5 file.  For the item whose name is NAME in
+// the group GROUP_ID, this function sets dv->tc to an Octave representation
+// of that item.  (dv must be a pointer to hdf5_callback_data.)  (It also
+// sets the other fields of dv).
+//
+// It returns 1 on success (in which case H5Giterate stops and returns),
+// -1 on error, and 0 to tell H5Giterate to continue on to the next item
+// (e.g. if NAME was a data type we don't recognize).
+
+static herr_t
+hdf5_read_next_data (hid_t group_id, const char *name, void *dv)
+{
+  hdf5_callback_data *d = static_cast <hdf5_callback_data *> (dv);
+
+  H5G_stat_t info;
+  herr_t retval = 0;
+  bool ident_valid = valid_identifier (name);
+
+  std::string vname = name;
+
+  // Allow identifiers as all digits so we can load lists saved by
+  // earlier versions of Octave.
+
+  if (! ident_valid && (d->import || ident_is_all_digits (vname)))
+    {
+      // fix the identifier, replacing invalid chars with underscores
+      vname = make_valid_identifier (vname);
+
+      // check again (in case vname was null, empty, or some such thing):
+      ident_valid = valid_identifier (vname); 
+    }
+
+  H5Gget_objinfo (group_id, name, 1, &info);
+
+  if (info.type == H5G_DATASET && ident_valid)
+    {
+      retval = 1;
+
+      hid_t data_id = H5Dopen (group_id, name);
+
+      if (data_id < 0)
+	{
+	  retval = data_id;
+
+	  goto done;
+	}
+
+      hid_t type_id = H5Dget_type (data_id);
+
+      hid_t type_class_id = H5Tget_class (type_id);
+
+#if HAVE_HDF5_INT2FLOAT_CONVERSIONS
+      if (type_class_id == H5T_INTEGER || type_class_id == H5T_FLOAT)
+	{
+#else
+      // hdf5 doesn't (yet) support automatic float/integer conversions
+      if (type_class_id == H5T_FLOAT)
+	{
+#endif
+	  // read real matrix or scalar variable
+
+	  hid_t space_id = H5Dget_space (data_id);
+
+	  hsize_t rank = H5Sget_simple_extent_ndims (space_id);
+
+	  if (rank == 0)
+	    {
+	      // real scalar:
+	      double dtmp;
+	      if (H5Dread (data_id, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, 
+			   H5P_DEFAULT, (void *) &dtmp) < 0)
+		retval = -1;  // error
+	      else
+		d->tc = dtmp;
+	    }
+	  else if (rank > 0 && rank <= 2)
+	    {
+	      // real matrix
+	      OCTAVE_LOCAL_BUFFER (hsize_t, dims, rank);
+	      OCTAVE_LOCAL_BUFFER (hsize_t, maxdims, rank);
+
+	      H5Sget_simple_extent_dims (space_id, dims, maxdims);
+
+	      int nr, nc;  // rows and columns
+	      // octave uses column-major & HDF5 uses row-major
+	      nc = dims[0];
+	      nr = rank > 1 ? dims[1] : 1;
+	      Matrix m (nr, nc);
+	      double *re = m.fortran_vec ();
+	      if (H5Dread (data_id, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, 
+			   H5P_DEFAULT, (void *) re) < 0)
+		retval = -1;  // error
+	      else
+		d->tc = m;
+	    }
+	  else if (rank >= 3 && d->import)
+	    {
+	      OCTAVE_LOCAL_BUFFER (hsize_t, dims, rank);
+	      OCTAVE_LOCAL_BUFFER (hsize_t, maxdims, rank);
+
+	      H5Sget_simple_extent_dims (space_id, dims, maxdims);
+
+	      OCTAVE_LOCAL_BUFFER (hssize_t, start, rank);
+	      OCTAVE_LOCAL_BUFFER (hsize_t, count, rank);
+
+	      for (hsize_t i = 0; i < rank; ++i)
+		{
+		  start[i] = 0;
+		  count[i] = 1;
+		}
+	      count[0] = dims[0];
+	      count[1] = dims[1];
+	      retval = hdf5_import_multidim (data_id, space_id,
+					     rank, dims, rank-1,
+					     start, count,
+					     H5T_NATIVE_DOUBLE, d->tc);
+	    }
+	  else
+	    {
+	      warning ("load: can't read %d-dim. hdf5 dataset %s",
+		       rank, name);
+	      retval = 0;  // skip; we can't read 3+ dimensional datasets
+	    }
+
+	  H5Sclose (space_id);
+	}
+      else if (type_class_id == H5T_STRING)
+	{
+	  // read string variable
+	  hid_t space_id = H5Dget_space (data_id);
+	  hsize_t rank = H5Sget_simple_extent_ndims (space_id);
+
+	  if (rank == 0)
+	    {
+	      // a single string:
+	      int slen = H5Tget_size (type_id);
+	      if (slen < 0)
+		retval = -1;  // error
+	      else
+		{
+		  OCTAVE_LOCAL_BUFFER (char, s, slen);
+		  // create datatype for (null-terminated) string
+		  // to read into:
+		  hid_t st_id = H5Tcopy (H5T_C_S1);
+		  H5Tset_size (st_id, slen);
+		  if (H5Dread (data_id, st_id, H5S_ALL, H5S_ALL, 
+			       H5P_DEFAULT, (void *) s) < 0)
+		    {
+		      retval = -1;  // error
+		    }
+		  else
+		    d->tc = s;
+
+		  H5Tclose (st_id);
+		}
+	    }
+	  else if (rank == 1)
+	    {
+	      // string vector
+	      hsize_t elements, maxdim;
+	      H5Sget_simple_extent_dims (space_id, &elements, &maxdim);
+	      int slen = H5Tget_size (type_id);
+	      if (slen < 0)
+		retval = -1;  // error
+	      else
+		{
+		  // hdf5 string arrays store strings of all the
+		  // same physical length (I think), which is
+		  // slightly wasteful, but oh well.
+
+		  OCTAVE_LOCAL_BUFFER (char, s, elements * slen);
+
+		  // create datatype for (null-terminated) string
+		  // to read into:
+		  hid_t st_id = H5Tcopy (H5T_C_S1);
+		  H5Tset_size (st_id, slen);
+
+		  if (H5Dread (data_id, st_id, H5S_ALL, H5S_ALL, 
+			       H5P_DEFAULT, (void *) s) < 0)
+		    retval = -1;  // error
+		  else
+		    {
+		      charMatrix chm (elements, slen - 1);
+		      for (hsize_t i = 0; i < elements; ++i)
+			{
+			  chm.insert (s + i*slen, i, 0);
+			}
+		      d->tc = octave_value (chm, true);
+		    }
+
+		  H5Tclose (st_id);
+		}
+	    }
+	  else
+	    {
+	      warning ("load: can't read %d-dim. hdf5 string vector %s",
+		       rank, name); 
+	      // skip; we can't read higher-dimensional string vectors
+	      retval = 0;
+	    }
+	}
+      else if (type_class_id == H5T_COMPOUND)
+	{
+	  // check for complex or range data:
+
+	  if (hdf5_types_compatible (type_id, d->complex_type))
+	    {
+	      // read complex matrix or scalar variable
+
+	      hid_t space_id = H5Dget_space (data_id);
+	      hsize_t rank = H5Sget_simple_extent_ndims (space_id);
+
+	      if (rank == 0)
+		{
+		  // complex scalar:
+		  Complex ctmp;
+		  if (H5Dread (data_id, d->complex_type, H5S_ALL,
+			       H5S_ALL, H5P_DEFAULT,
+			       (void *) X_CAST (double *, &ctmp)) < 0)
+		    retval = -1;  // error
+		  else
+		    d->tc = ctmp;
+		}
+	      else if (rank > 0 && rank <= 2)
+		{
+		  // complex matrix
+		  OCTAVE_LOCAL_BUFFER (hsize_t, dims, rank);
+		  OCTAVE_LOCAL_BUFFER (hsize_t, maxdims, rank);
+		  H5Sget_simple_extent_dims (space_id, dims, maxdims);
+		  int nr, nc;  // rows and columns
+		  // octave uses column-major & HDF5 uses row-major
+		  nc = dims[0];
+		  nr = rank > 1 ? dims[1] : 1;
+		  ComplexMatrix m (nr, nc);
+		  Complex *reim = m.fortran_vec ();
+		  if (H5Dread (data_id, d->complex_type, H5S_ALL,
+			       H5S_ALL, H5P_DEFAULT,
+			       (void *) X_CAST (double *, reim)) < 0)
+		    retval = -1;  // error
+		  else
+		    d->tc = m;
+		}
+	      else if (rank >= 3 && d->import)
+		{
+		  OCTAVE_LOCAL_BUFFER (hsize_t, dims, rank);
+		  OCTAVE_LOCAL_BUFFER (hsize_t, maxdims, rank);
+		  H5Sget_simple_extent_dims (space_id, dims, maxdims);
+		  OCTAVE_LOCAL_BUFFER (hssize_t, start, rank);
+		  OCTAVE_LOCAL_BUFFER (hsize_t, count, rank);
+		  for (hsize_t i = 0; i < rank; ++i)
+		    {
+		      start[i] = 0;
+		      count[i] = 1;
+		    }
+		  count[0] = dims[0];
+		  count[1] = dims[1];
+		  retval = hdf5_import_multidim (data_id, space_id,
+						 rank, dims, rank-1,
+						 start, count,
+						 d->complex_type,
+						 d->tc);
+		}
+	      else
+		{
+		  warning ("load: can't read %d-dim. hdf5 dataset %s",
+			   rank, name);
+		  // skip; we can't read 3+ dimensional datasets
+		  retval = 0;
+		}
+	      H5Sclose (space_id);
+	    }
+	  else if (hdf5_types_compatible (type_id, d->range_type))
+	    {
+	      // read range variable:
+	      hid_t space_id = H5Dget_space (data_id);
+	      hsize_t rank = H5Sget_simple_extent_ndims (space_id);
+
+	      if (rank == 0)
+		{
+		  double rangevals[3];
+		  if (H5Dread (data_id, d->range_type, H5S_ALL, H5S_ALL, 
+			       H5P_DEFAULT, (void *) rangevals) < 0)
+		    retval = -1;  // error
+		  else
+		    {
+		      Range r (rangevals[0], rangevals[1], rangevals[2]);
+		      d->tc = r;
+		    }
+		}
+	      else
+		{
+		  warning ("load: can't read range array `%s' in hdf5 file",
+			   name);
+		  // skip; we can't read arrays of range variables
+		  retval = 0;
+		}
+
+	      H5Sclose (space_id);
+	    }
+	  else
+	    {
+	      warning ("load: can't read `%s' (unknown compound datatype)",
+		       name);
+	      retval = 0; // unknown datatype; skip.
+	    }
+	}
+      else
+	{
+	  warning ("load: can't read `%s' (unknown datatype)", name);
+	  retval = 0; // unknown datatype; skip
+	}
+
+      H5Tclose (type_id);
+
+      // check for OCTAVE_GLOBAL attribute:
+      d->global = hdf5_check_attr (data_id, "OCTAVE_GLOBAL");
+
+      H5Dclose (data_id);
+    }
+  else if (info.type == H5G_GROUP && ident_valid)
+    {
+      // read in group as a list or a structure
+      retval = 1;
+
+      hid_t subgroup_id = H5Gopen (group_id, name);
+
+      if (subgroup_id < 0)
+	{
+	  retval = subgroup_id;
+	  goto done;
+	}
+
+      // an HDF5 group is treated as an octave structure by
+      // default (since that preserves name information), and an
+      // octave list otherwise.
+
+      bool is_list = hdf5_check_attr (subgroup_id, "OCTAVE_LIST");
+
+      hdf5_callback_data dsub;
+
+      dsub.complex_type = d->complex_type;
+      dsub.range_type = d->range_type;
+      dsub.import = d->import;
+
+      herr_t retval2;
+      octave_value_list lst;
+      Octave_map m;
+      int current_item = 0;
+      while ((retval2 = H5Giterate (group_id, name, &current_item,
+				    hdf5_read_next_data, &dsub)) > 0)
+	{
+	  if (is_list)
+	    lst.append (dsub.tc);
+	  else
+	    {
+	      octave_value ov = dsub.tc;
+
+	      if (ov.is_list ())
+		m [dsub.name] = ov.list_value ();
+	      else
+		m [dsub.name] = ov;
+	    }
+
+	  if (have_h5giterate_bug)
+	    current_item++;  // H5Giterate returned the last index processed
+	}
+
+      if (retval2 < 0)
+	retval = retval2;
+      else
+	{
+	  d->global = hdf5_check_attr (group_id, "OCTAVE_GLOBAL");
+
+	  if (is_list)
+	    d->tc = octave_value (lst);
+	  else
+	    d->tc = m;
+	}
+
+      H5Gclose (subgroup_id);
+    }
+  else if (! ident_valid)
+    {
+      // should we attempt to handle invalid identifiers by converting
+      // bad characters to '_', say?
+      warning ("load: skipping invalid identifier `%s' in hdf5 file",
+	       name);
+    }
+
+ done:
+
+  if (retval < 0)
+    error ("load: error while reading hdf5 item %s", name);
+
+  if (retval > 0)
+    {
+      // get documentation string, if any:
+      int comment_length = H5Gget_comment (group_id, name, 0, 0);
+
+      if (comment_length > 1)
+	{
+	  OCTAVE_LOCAL_BUFFER (char, tdoc, comment_length);
+	  H5Gget_comment (group_id, name, comment_length, tdoc);
+	  d->doc = tdoc;
+	}
+      else if (vname != name)
+	{
+	  // the name was changed by import; store the original name
+	  // as the documentation string:
+	  d->doc = name;
+	}
+
+      // copy name (actually, vname):
+      d->name = vname;
+    }
+
+  return retval;
+}
+
+// Read the next Octave variable from the stream IS, which must really be
+// an hdf5_ifstream.  Return the variable value in tc, its doc string
+// in doc, and whether it is global in global.  The return value is
+// the name of the variable, or NULL if none were found or there was
+// and error.  If import is true, we try extra hard to import "foreign"
+// datasets (not created by Octave), although we usually do a reasonable
+// job anyway.  (c.f. load -import documentation.)
+std::string
+read_hdf5_data (std::istream& is,
+		const std::string& filename, bool& global,
+		octave_value& tc, std::string& doc, bool import)
+{
+  std::string retval;
+
+  doc.resize (0);
+
+  hdf5_ifstream& hs = (hdf5_ifstream&) is;
+  hdf5_callback_data d;
+
+  d.import = import;
+
+  // Versions of HDF5 prior to 1.2.2 had a bug in H5Giterate where it
+  // would return the index of the last item processed instead of the
+  // next item to be processed, forcing us to increment the index manually.
+
+  unsigned int vers_major, vers_minor, vers_release;
+
+  H5get_libversion (&vers_major, &vers_minor, &vers_release);
+
+  // XXX FIXME XXX -- this test looks wrong.
+  have_h5giterate_bug
+    = (vers_major < 1
+       || (vers_major == 1 && (vers_minor < 2
+			       || (vers_minor == 2 && vers_release < 2))));
+
+  herr_t H5Giterate_retval = H5Giterate (hs.file_id, "/", &hs.current_item,
+					 hdf5_read_next_data, &d);
+
+  if (have_h5giterate_bug)
+    {
+      // H5Giterate sets current_item to the last item processed; we want
+      // the index of the next item (for the next call to read_hdf5_data)
+
+      hs.current_item++;
+    }
+
+  if (H5Giterate_retval > 0)
+    {
+      global = d.global;
+      tc = d.tc;
+      doc = d.doc;
+    }
+  else
+    {
+      // an error occurred (H5Giterate_retval < 0) or there are no
+      // more datasets print an error message if retval < 0?
+      // hdf5_read_next_data already printed one, probably.
+    }
+
+  H5Tclose (d.complex_type);
+  H5Tclose (d.range_type);
+
+  if (! d.name.empty ())
+    retval = d.name;
+
+  return retval;
+}
+
+// Add an attribute named attr_name to loc_id (a simple scalar
+// attribute with value 1).  Return value is >= 0 on success.
+static herr_t
+hdf5_add_attr (hid_t loc_id, const char *attr_name)
+{
+  herr_t retval = 0;
+
+  hid_t as_id = H5Screate (H5S_SCALAR);
+
+  if (as_id >= 0)
+    {
+      hid_t a_id = H5Acreate (loc_id, attr_name,
+			      H5T_NATIVE_UCHAR, as_id, H5P_DEFAULT);
+
+      if (a_id >= 0)
+	{
+	  unsigned char attr_val = 1;
+
+	  retval = H5Awrite (a_id, H5T_NATIVE_UCHAR, (void*) &attr_val);
+
+	  H5Aclose (a_id);
+	}
+      else
+	retval = a_id;
+
+      H5Sclose (as_id);
+    }
+  else
+    retval = as_id;
+
+  return retval;
+}
+
+
+// save_type_to_hdf5 is not currently used, since hdf5 doesn't yet support
+// automatic float<->integer conversions:
+
+#if HAVE_HDF5_INT2FLOAT_CONVERSIONS
+
+// return the HDF5 type id corresponding to the Octave save_type
+
+static hid_t
+save_type_to_hdf5 (save_type st)
+{
+  switch (st)
+    {
+    case LS_U_CHAR:
+      return H5T_NATIVE_UCHAR;
+
+    case LS_U_SHORT:
+      return H5T_NATIVE_USHORT;
+
+    case LS_U_INT:
+      return H5T_NATIVE_UINT;
+
+    case LS_CHAR:
+      return H5T_NATIVE_CHAR;
+
+    case LS_SHORT:
+      return H5T_NATIVE_SHORT;
+
+    case LS_INT:
+      return H5T_NATIVE_INT;
+
+    case LS_FLOAT:
+      return H5T_NATIVE_FLOAT;
+
+    case LS_DOUBLE:
+    default:
+      return H5T_NATIVE_DOUBLE;
+    }
+}
+#endif /* HAVE_HDF5_INT2FLOAT_CONVERSIONS */
+
+// Add the data from TC to the HDF5 location loc_id, which could
+// be either a file or a group within a file.  Return true if
+// successful.  This function calls itself recursively for lists
+// (stored as HDF5 groups).
+
+static bool
+add_hdf5_data (hid_t loc_id, const octave_value& tc,
+	       const std::string& name, const std::string& doc,
+	       bool mark_as_global, bool save_as_floats)
+{
+  hsize_t dims[3];
+  hid_t type_id = -1, space_id = -1, data_id = -1;
+  bool data_is_group = 0;
+  bool retval = 0;
+
+  if (tc.is_string ())
+    {
+      int nr = tc.rows ();
+      charMatrix chm = tc.char_matrix_value ();
+      int nc = chm.cols ();
+
+      // create datatype for (null-terminated) string to write from:
+      type_id = H5Tcopy (H5T_C_S1); H5Tset_size (type_id, nc + 1);
+      if (type_id < 0)
+	goto error_cleanup;
+
+      dims[0] = nr;
+      space_id = H5Screate_simple (nr > 0 ? 1 : 0, dims, (hsize_t*) 0);
+      if (space_id < 0)
+	goto error_cleanup;
+
+      data_id = H5Dcreate (loc_id, name.c_str (), 
+			   type_id, space_id, H5P_DEFAULT);
+      if (data_id < 0)
+	goto error_cleanup;
+
+      OCTAVE_LOCAL_BUFFER (char, s, nr * (nc + 1));
+
+      for (int i = 0; i < nr; ++i)
+	{
+	  std::string tstr = chm.row_as_string (i);
+	  strcpy (s + i * (nc+1), tstr.c_str ());
+	}
+
+      if (H5Dwrite (data_id, type_id, H5S_ALL, H5S_ALL, H5P_DEFAULT,
+		    (void*) s) < 0) {
+	goto error_cleanup;
+      }
+    }
+  else if (tc.is_range ())
+    {
+      space_id = H5Screate_simple (0, dims, (hsize_t*) 0);
+      if (space_id < 0)
+	goto error_cleanup;
+
+      type_id = hdf5_make_range_type (H5T_NATIVE_DOUBLE);
+      if (type_id < 0)
+	goto error_cleanup;
+
+      data_id = H5Dcreate (loc_id, name.c_str (), 
+			   type_id, space_id, H5P_DEFAULT);
+      if (data_id < 0)
+	goto error_cleanup;
+    
+      Range r = tc.range_value ();
+      double range_vals[3];
+      range_vals[0] = r.base ();
+      range_vals[1] = r.limit ();
+      range_vals[2] = r.inc ();
+
+      if (H5Dwrite (data_id, type_id, H5S_ALL, H5S_ALL, H5P_DEFAULT,
+		    (void*) range_vals) < 0)
+	goto error_cleanup;
+    }
+  else if (tc.is_real_scalar ())
+    {
+      space_id = H5Screate_simple (0, dims, (hsize_t*) 0);
+      if (space_id < 0) goto error_cleanup;
+
+      data_id = H5Dcreate (loc_id, name.c_str (), 
+			   H5T_NATIVE_DOUBLE, space_id, H5P_DEFAULT);
+      if (data_id < 0)
+	goto error_cleanup;
+
+      double tmp = tc.double_value ();
+      if (H5Dwrite (data_id, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL,
+		    H5P_DEFAULT, (void*) &tmp) < 0)
+	goto error_cleanup;
+    }
+  else if (tc.is_real_matrix ())
+    {
+      Matrix m = tc.matrix_value ();
+      dims[1] = m.rows ();    // Octave uses column-major, while
+      dims[0] = m.columns (); // HDF5 uses row-major ordering
+
+      space_id = H5Screate_simple (dims[1] > 1 ?2:1, dims, (hsize_t*) 0);
+      if (space_id < 0)
+	goto error_cleanup;
+
+      hid_t save_type_id = H5T_NATIVE_DOUBLE;
+
+      if (save_as_floats)
+	{
+	  if (m.too_large_for_float ())
+	    {
+	      warning ("save: some values too large to save as floats --");
+	      warning ("save: saving as doubles instead");
+	    }
+	  else
+	    save_type_id = H5T_NATIVE_FLOAT;
+	}
+#if HAVE_HDF5_INT2FLOAT_CONVERSIONS
+      // hdf5 currently doesn't support float/integer conversions
+      else
+	{
+	  double max_val, min_val;
+
+	  if (m.all_integers (max_val, min_val))
+	    save_type_id
+	      = save_type_to_hdf5 (get_save_type (max_val, min_val));
+	}
+#endif /* HAVE_HDF5_INT2FLOAT_CONVERSIONS */
+
+      data_id = H5Dcreate (loc_id, name.c_str (), 
+			   save_type_id, space_id, H5P_DEFAULT);
+      if (data_id < 0)
+	goto error_cleanup;
+    
+      double *mtmp = m.fortran_vec ();
+      if (H5Dwrite (data_id, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL,
+		    H5P_DEFAULT, (void*) mtmp) < 0)
+	goto error_cleanup;
+    }
+  else if (tc.is_complex_scalar ())
+    {
+      space_id = H5Screate_simple (0, dims, (hsize_t*) 0);
+      if (space_id < 0)
+	goto error_cleanup;
+
+      type_id = hdf5_make_complex_type (H5T_NATIVE_DOUBLE);
+      if (type_id < 0)
+	goto error_cleanup;
+
+      data_id = H5Dcreate (loc_id, name.c_str (), 
+			   type_id, space_id, H5P_DEFAULT);
+      if (data_id < 0)
+	goto error_cleanup;
+
+      Complex tmp = tc.complex_value ();
+      if (H5Dwrite (data_id, type_id, H5S_ALL, H5S_ALL, H5P_DEFAULT,
+		    (void*) X_CAST (double*, &tmp)) < 0)
+	goto error_cleanup;
+    }
+  else if (tc.is_complex_matrix ())
+    {
+      ComplexMatrix m = tc.complex_matrix_value ();
+
+      dims[1] = m.rows ();    // Octave uses column-major, while
+      dims[0] = m.columns (); // HDF5 uses row-major ordering
+
+      space_id = H5Screate_simple (dims[1] > 1 ?2:1, dims, (hsize_t*) 0);
+      if (space_id < 0)
+	goto error_cleanup;
+
+      hid_t save_type_id = H5T_NATIVE_DOUBLE;
+
+      if (save_as_floats)
+	{
+	  if (m.too_large_for_float ())
+	    {
+	      warning ("save: some values too large to save as floats --");
+	      warning ("save: saving as doubles instead");
+	    }
+	  else
+	    save_type_id = H5T_NATIVE_FLOAT;
+	}
+#if HAVE_HDF5_INT2FLOAT_CONVERSIONS
+      // hdf5 currently doesn't support float/integer conversions
+      else
+	{
+	  double max_val, min_val;
+
+	  if (m.all_integers (max_val, min_val))
+	    save_type_id
+	      = save_type_to_hdf5 (get_save_type (max_val, min_val));
+	}
+#endif /* HAVE_HDF5_INT2FLOAT_CONVERSIONS */
+
+      type_id = hdf5_make_complex_type (save_type_id);
+      if (type_id < 0) goto error_cleanup;
+
+      data_id = H5Dcreate (loc_id, name.c_str (), 
+			   type_id, space_id, H5P_DEFAULT);
+      if (data_id < 0)
+	goto error_cleanup;
+    
+      hid_t complex_type_id = hdf5_make_complex_type (H5T_NATIVE_DOUBLE);
+      if (complex_type_id < 0)
+	goto error_cleanup;
+
+      Complex *mtmp = m.fortran_vec ();
+      if (H5Dwrite (data_id, complex_type_id, H5S_ALL, H5S_ALL, H5P_DEFAULT,
+		    (void*) X_CAST (double *, mtmp)) < 0)
+	{
+	  H5Tclose (complex_type_id);
+	  goto error_cleanup;
+	}
+
+      H5Tclose (complex_type_id);
+    }
+  else if (tc.is_list ())
+    {
+      data_id = H5Gcreate (loc_id, name.c_str (), 0);
+      if (data_id < 0)
+	goto error_cleanup;
+
+      data_is_group = 1;
+
+      // recursively add each element of the list to this group
+      octave_value_list lst = tc.list_value ();
+
+      for (int i = 0; i < lst.length (); ++i)
+	{
+	  // should we use lst.name_tags () to label the elements?
+	  char s[20];
+	  sprintf (s, "_%d", i);
+	  bool retval2 = add_hdf5_data (data_id, lst (i), s, "",
+					false, save_as_floats);
+	  if (! retval2)
+	    goto error_cleanup;
+	}
+
+      // mark with an attribute "OCTAVE_LIST" with value 1
+      // to distinguish from structures (also stored as HDF5 groups):
+      if (hdf5_add_attr (data_id, "OCTAVE_LIST") < 0)
+	goto error_cleanup;
+    }
+  else if (tc.is_map ())
+    {
+      // an Octave structure
+      data_id = H5Gcreate (loc_id, name.c_str (), 0);
+      if (data_id < 0)
+	goto error_cleanup;
+
+      data_is_group = 1;
+
+      // recursively add each element of the structure to this group
+      Octave_map m = tc.map_value ();
+      Octave_map::iterator i = m.begin ();
+      while (i != m.end ())
+	{
+	  // XXX FIXME XXX -- if the length of the structure array is
+	  // 1, should we really create a list object?
+	  bool retval2 = add_hdf5_data (data_id, octave_value (m.contents (i)),
+					m.key (i), "", false, save_as_floats);
+	  if (! retval2)
+	    goto error_cleanup;
+
+	  i++;
+	}
+    }
+  else
+    {
+      gripe_wrong_type_arg ("save", tc, false);
+      goto error_cleanup;
+    }
+
+  // attach doc string as comment:
+  if (doc.length () > 0
+      && H5Gset_comment (loc_id, name.c_str (), doc.c_str ()) < 0)
+    goto error_cleanup;
+
+  retval = 1;
+
+  // if it's global, add an attribute "OCTAVE_GLOBAL" with value 1
+  if (mark_as_global)
+    retval = hdf5_add_attr (data_id, "OCTAVE_GLOBAL") >= 0;
+
+ error_cleanup:
+
+  if (! retval)
+    error ("save: error while writing `%s' to hdf5 file", name.c_str ());
+
+  if (data_id >= 0)
+    {
+      if (data_is_group)
+	H5Gclose (data_id);
+      else
+	H5Dclose (data_id);
+    }
+
+  if (space_id >= 0)
+    H5Sclose (space_id);
+
+  if (type_id >= 0)
+    H5Tclose (type_id);
+
+  return retval;
+}
+
+// Write data from TC in HDF5 (binary) format to the stream OS,
+// which must be an hdf5_ofstream, returning true on success.
+
+bool
+save_hdf5_data (std::ostream& os, const octave_value& tc,
+		const std::string& name, const std::string& doc,
+		bool mark_as_global, bool save_as_floats)
+{
+  hdf5_ofstream& hs = (hdf5_ofstream&) os;
+
+  return add_hdf5_data (hs.file_id, tc, name, doc,
+			mark_as_global, save_as_floats);
+}
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
diff --git a/src/ls-hdf5.h b/src/ls-hdf5.h
new file mode 100644
--- /dev/null
+++ b/src/ls-hdf5.h
@@ -0,0 +1,128 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_hdf5_h)
+#define octave_ls_hdf5_h 1
+
+// first, we need to define our own dummy stream subclass, since
+// HDF5 needs to do its own file i/o
+
+// hdf5_fstreambase is used for both input and output streams, modeled
+// on the fstreambase class in <fstream.h>
+
+class hdf5_fstreambase : virtual public std::ios
+{
+public:
+
+  // HDF5 uses an "id" to refer to an open file
+  hid_t file_id;
+
+  // keep track of current item index in the file
+  int current_item;
+
+  hdf5_fstreambase () { file_id = -1; }
+
+  hdf5_fstreambase (const char *name, int mode, int prot = 0)
+    {
+      if (mode == std::ios::in)
+	file_id = H5Fopen (name, H5F_ACC_RDONLY, H5P_DEFAULT);
+      else if (mode == std::ios::out)
+	file_id = H5Fcreate (name, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
+
+      if (file_id < 0)
+	std::ios::setstate (std::ios::badbit);
+
+      current_item = 0;
+    }
+
+  void close ()
+    { 
+      if (file_id >= 0)
+	{
+	  if (H5Fclose (file_id) < 0)
+	    std::ios::setstate (std::ios::badbit);
+	  file_id = -1;
+	}
+    }
+
+  void open (const char *name, int mode, int prot = 0)
+    {
+      clear ();
+
+      if (mode == std::ios::in)
+	file_id = H5Fopen (name, H5F_ACC_RDONLY, H5P_DEFAULT);
+      else if (mode == std::ios::out)
+	file_id = H5Fcreate (name, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
+
+      if (file_id < 0)
+	std::ios::setstate (std::ios::badbit);
+
+      current_item = 0;
+    }
+};
+
+// input and output streams, subclassing istream and ostream
+// so that we can pass them for stream parameters in the functions below.
+
+class hdf5_ifstream : public hdf5_fstreambase, public std::istream
+{
+public:
+
+  hdf5_ifstream () : hdf5_fstreambase (), std::istream (0) { }
+
+  hdf5_ifstream (const char *name, int mode = std::ios::in, int prot = 0)
+    : hdf5_fstreambase (name, mode, prot), std::istream (0) { }
+
+  void open (const char *name, int mode = std::ios::in, int prot = 0)
+    { hdf5_fstreambase::open (name, mode, prot); }
+};
+
+class hdf5_ofstream : public hdf5_fstreambase, public std::ostream
+{
+public:
+
+  hdf5_ofstream () : hdf5_fstreambase (), std::ostream (0) { }
+
+  hdf5_ofstream (const char *name, int mode = std::ios::out, int prot = 0)
+    : hdf5_fstreambase (name, mode, prot), std::ostream (0) { }
+
+  void open (const char *name, int mode = std::ios::out, int prot = 0)
+    { hdf5_fstreambase::open (name, mode, prot); }
+};
+
+extern std::string
+read_hdf5_data (std::istream& is,
+		const std::string& filename, bool& global,
+		octave_value& tc, std::string& doc, bool import);
+
+extern bool
+save_hdf5_data (std::ostream& os, const octave_value& tc,
+		const std::string& name, const std::string& doc,
+		bool mark_as_global, bool save_as_floats);
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
diff --git a/src/ls-mat-ascii.cc b/src/ls-mat-ascii.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-mat-ascii.cc
@@ -0,0 +1,339 @@
+/*
+
+Copyright (C) 1996, 1997 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <cfloat>
+#include <cstring>
+#include <cctype>
+
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#include <memory>
+#include <string>
+
+#include "byte-swap.h"
+#include "data-conv.h"
+#include "file-ops.h"
+#include "glob-match.h"
+#include "lo-mappers.h"
+#include "lo-sstream.h"
+#include "mach-info.h"
+#include "oct-env.h"
+#include "oct-time.h"
+#include "quit.h"
+#include "str-vec.h"
+
+#include "Cell.h"
+#include "defun.h"
+#include "error.h"
+#include "gripes.h"
+#include "load-save.h"
+#include "oct-obj.h"
+#include "oct-map.h"
+#include "ov-cell.h"
+#include "pager.h"
+#include "pt-exp.h"
+#include "symtab.h"
+#include "sysdep.h"
+#include "unwind-prot.h"
+#include "utils.h"
+#include "variables.h"
+#include "version.h"
+#include "dMatrix.h"
+
+#include "ls-mat-ascii.h"
+
+static std::string
+get_mat_data_input_line (std::istream& is)
+{
+  std::string retval;
+
+  bool have_data = false;
+
+  do
+    {
+      retval = "";
+
+      char c;
+      while (is.get (c))
+	{
+	  if (c == '\n' || c == '\r')
+	    break;
+
+	  if (c == '%' || c == '#')
+	    {
+	      // skip to end of line
+	      while (is.get (c))
+		if (c == '\n' || c == '\r')
+		  break;
+
+	      break;
+	    }
+
+	  if (! is.eof ())
+	    {
+	      if (! have_data && c != ' ' && c != '\t')
+		have_data = true;
+
+	      retval += c;
+	    }
+	}
+    }
+  while (! (have_data || is.eof ()));
+
+  return retval;
+}
+
+static void
+get_lines_and_columns (std::istream& is, const std::string& filename, int& nr, int& nc)
+{
+  std::streampos pos = is.tellg ();
+
+  int file_line_number = 0;
+
+  nr = 0;
+  nc = 0;
+
+  while (is && ! error_state)
+    {
+      OCTAVE_QUIT;
+
+      std::string buf = get_mat_data_input_line (is);
+
+      file_line_number++;
+
+      size_t beg = buf.find_first_not_of (", \t");
+
+      // If we see a CR as the last character in the buffer, we had a
+      // CRLF pair as the line separator.  Any other CR in the text
+      // will not be considered as whitespace.
+
+      if (beg != NPOS && buf[beg] == '\r' && beg == buf.length () - 1)
+	{
+	  // We had a blank line ending with a CRLF.  Handle it the
+	  // same as an empty line.
+	  beg = NPOS;
+	}
+
+      int tmp_nc = 0;
+
+      while (beg != NPOS)
+	{
+	  tmp_nc++;
+
+	  size_t end = buf.find_first_of (", \t", beg);
+
+	  if (end != NPOS)
+	    {
+	      beg = buf.find_first_not_of (", \t", end);
+
+	      if (buf[beg] == '\r' && beg == buf.length () - 1)
+		{
+		  // We had a line with trailing spaces and
+		  // ending with a CRLF, so this should look like EOL,
+		  // not a new colum.
+		  break;
+		}
+	    }
+	  else
+	    break;
+	}
+
+      if (tmp_nc > 0)
+	{
+	  if (nc == 0)
+	    {
+	      nc = tmp_nc;
+	      nr++;
+	    }
+	  else if (nc == tmp_nc)
+	    nr++;
+	  else
+	    error ("load: %s: inconsistent number of columns near line %d",
+		   filename.c_str (), file_line_number);
+	}
+    }
+
+  if (nr == 0 || nc == 0)
+    error ("load: file `%s' seems to be empty!", filename.c_str ());
+
+  is.clear ();
+  is.seekg (pos, std::ios::beg);
+}
+
+// Extract a matrix from a file of numbers only.
+//
+// Comments are not allowed.  The file should only have numeric values.
+//
+// Reads the file twice.  Once to find the number of rows and columns,
+// and once to extract the matrix.
+//
+// FILENAME is used for error messages.
+//
+// This format provides no way to tag the data as global.
+
+std::string
+read_mat_ascii_data (std::istream& is, const std::string& filename,
+		     octave_value& tc)
+{
+  std::string retval;
+
+  std::string varname;
+
+  size_t pos = filename.rfind ('/');
+
+  if (pos != NPOS)
+    varname = filename.substr (pos+1);
+  else
+    varname = filename;
+
+  pos = varname.find ('.');
+
+  if (pos != NPOS)
+    varname = varname.substr (0, pos);
+
+  size_t len = varname.length ();
+  for (size_t i = 0; i < len; i++)
+    {
+      char c = varname[i];
+      if (! (isalnum (c) || c == '_'))
+	varname[i] = '_';
+    }
+
+  if (! isalpha (varname[0]))
+    varname.insert (0, "X");
+
+  if (valid_identifier (varname))
+    {
+      int nr = 0;
+      int nc = 0;
+
+      int total_count = 0;
+
+      get_lines_and_columns (is, filename, nr, nc);
+
+      OCTAVE_QUIT;
+
+      if (! error_state && nr > 0 && nc > 0)
+	{
+	  Matrix tmp (nr, nc);
+
+	  if (nr < 1 || nc < 1)
+	    is.clear (std::ios::badbit);
+	  else
+	    {
+	      double d;
+	      for (int i = 0; i < nr; i++)
+		{
+		  std::string buf = get_mat_data_input_line (is);
+
+#ifdef HAVE_SSTREAM
+		  std::istringstream tmp_stream (buf);
+#else
+		  std::istrstream tmp_stream (buf.c_str ());
+#endif
+
+		  for (int j = 0; j < nc; j++)
+		    {
+		      OCTAVE_QUIT;
+
+		      d = octave_read_double (tmp_stream);
+
+		      if (tmp_stream || tmp_stream.eof ())
+			{
+			  tmp.elem (i, j) = d;
+			  total_count++;
+
+			  // Skip whitespace and commas.
+			  char c;
+			  while (1)
+			    {
+			      tmp_stream >> c;
+
+			      if (! tmp_stream)
+				break;
+
+			      if (! (c == ' ' || c == '\t' || c == ','))
+				{
+				  tmp_stream.putback (c);
+				  break;
+				}
+			    }
+
+			  if (tmp_stream.eof ())
+			    break;
+			}
+		      else
+			{
+			  error ("load: failed to read matrix from file `%s'",
+				 filename.c_str ());
+
+			  return retval;
+			}
+
+		    }
+		}
+	    }
+
+	  if (is || is.eof ())
+	    {
+	      // XXX FIXME XXX -- not sure this is best, but it works.
+
+	      if (is.eof ())
+		is.clear ();
+
+	      int expected = nr * nc;
+
+	      if (expected == total_count)
+		{
+		  tc = tmp;
+		  retval = varname;
+		}
+	      else
+		error ("load: expected %d elements, found %d",
+		       expected, total_count);
+	    }
+	  else
+	    error ("load: failed to read matrix from file `%s'",
+		   filename.c_str ());
+	}
+      else
+	error ("load: unable to extract matrix size from file `%s'",
+	       filename.c_str ());
+    }
+  else
+    error ("load: unable to convert filename `%s' to valid identifier",
+	   filename.c_str ());
+
+  return retval;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-mat-ascii.h b/src/ls-mat-ascii.h
new file mode 100644
--- /dev/null
+++ b/src/ls-mat-ascii.h
@@ -0,0 +1,37 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_mat_ascii_h)
+#define octave_ls_mat_ascii_h 1
+
+extern std::string
+read_mat_ascii_data (std::istream& is, const std::string& filename,
+		     octave_value& tc);
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-mat4.cc b/src/ls-mat4.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-mat4.cc
@@ -0,0 +1,490 @@
+/*
+
+Copyright (C) 1996, 1997 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <cfloat>
+#include <cstring>
+#include <cctype>
+
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#include <memory>
+#include <string>
+
+#include "byte-swap.h"
+#include "data-conv.h"
+#include "file-ops.h"
+#include "glob-match.h"
+#include "lo-mappers.h"
+#include "lo-sstream.h"
+#include "mach-info.h"
+#include "oct-env.h"
+#include "oct-time.h"
+#include "quit.h"
+#include "str-vec.h"
+
+#include "Cell.h"
+#include "defun.h"
+#include "error.h"
+#include "gripes.h"
+#include "load-save.h"
+#include "oct-obj.h"
+#include "oct-map.h"
+#include "ov-cell.h"
+#include "pager.h"
+#include "pt-exp.h"
+#include "symtab.h"
+#include "sysdep.h"
+#include "unwind-prot.h"
+#include "utils.h"
+#include "variables.h"
+#include "version.h"
+#include "dMatrix.h"
+
+#include "ls-mat4.h"
+
+// Read LEN elements of data from IS in the format specified by
+// PRECISION, placing the result in DATA.  If SWAP is TRUE, swap
+// the bytes of each element before copying to DATA.  FLT_FMT
+// specifies the format of the data if we are reading floating point
+// numbers.
+
+static void
+read_mat_binary_data (std::istream& is, double *data, int precision,
+		      int len, bool swap,
+		      oct_mach_info::float_format flt_fmt)
+{
+  switch (precision)
+    {
+    case 0:
+      read_doubles (is, data, LS_DOUBLE, len, swap, flt_fmt);
+      break;
+
+    case 1:
+      read_doubles (is, data, LS_FLOAT, len, swap, flt_fmt);
+      break;
+
+    case 2:
+      read_doubles (is, data, LS_INT, len, swap, flt_fmt);
+      break;
+
+    case 3:
+      read_doubles (is, data, LS_SHORT, len, swap, flt_fmt);
+      break;
+
+    case 4:
+      read_doubles (is, data, LS_U_SHORT, len, swap, flt_fmt);
+      break;
+
+    case 5:
+      read_doubles (is, data, LS_U_CHAR, len, swap, flt_fmt);
+      break;
+
+    default:
+      break;
+    }
+}
+
+int
+read_mat_file_header (std::istream& is, bool& swap, FOUR_BYTE_INT& mopt, 
+		      FOUR_BYTE_INT& nr, FOUR_BYTE_INT& nc,
+		      FOUR_BYTE_INT& imag, FOUR_BYTE_INT& len,
+		      int quiet)
+{
+  swap = false;
+
+  // We expect to fail here, at the beginning of a record, so not
+  // being able to read another mopt value should not result in an
+  // error.
+
+  is.read (X_CAST (char *, &mopt), 4);
+  if (! is)
+    return 1;
+
+  if (! is.read (X_CAST (char *, &nr), 4))
+    goto data_read_error;
+
+  if (! is.read (X_CAST (char *, &nc), 4))
+    goto data_read_error;
+
+  if (! is.read (X_CAST (char *, &imag), 4))
+    goto data_read_error;
+
+  if (! is.read (X_CAST (char *, &len), 4))
+    goto data_read_error;
+
+// If mopt is nonzero and the byte order is swapped, mopt will be
+// bigger than we expect, so we swap bytes.
+//
+// If mopt is zero, it means the file was written on a little endian
+// machine, and we only need to swap if we are running on a big endian
+// machine.
+//
+// Gag me.
+
+  if (oct_mach_info::words_big_endian () && mopt == 0)
+    swap = true;
+
+  // mopt is signed, therefore byte swap may result in negative value.
+
+  if (mopt > 9999 || mopt < 0)
+    swap = true;
+
+  if (swap)
+    {
+      swap_4_bytes (X_CAST (char *, &mopt));
+      swap_4_bytes (X_CAST (char *, &nr));
+      swap_4_bytes (X_CAST (char *, &nc));
+      swap_4_bytes (X_CAST (char *, &imag));
+      swap_4_bytes (X_CAST (char *, &len));
+    }
+
+  if (mopt > 9999 || mopt < 0 || imag > 1 || imag < 0)
+    {
+      if (! quiet)
+	error ("load: can't read binary file");
+      return -1;
+    }
+
+  return 0;
+
+ data_read_error:
+  return -1;
+}
+
+// We don't just use a cast here, because we need to be able to detect
+// possible errors.
+
+oct_mach_info::float_format
+mopt_digit_to_float_format (int mach)
+{
+  oct_mach_info::float_format flt_fmt = oct_mach_info::flt_fmt_unknown;
+
+  switch (mach)
+    {
+    case 0:
+      flt_fmt = oct_mach_info::flt_fmt_ieee_little_endian;
+      break;
+
+    case 1:
+      flt_fmt = oct_mach_info::flt_fmt_ieee_big_endian;
+      break;
+
+    case 2:
+      flt_fmt = oct_mach_info::flt_fmt_vax_d;
+      break;
+
+    case 3:
+      flt_fmt = oct_mach_info::flt_fmt_vax_g;
+      break;
+
+    case 4:
+      flt_fmt = oct_mach_info::flt_fmt_cray;
+      break;
+
+    default:
+      flt_fmt = oct_mach_info::flt_fmt_unknown;
+      break;
+    }
+
+  return flt_fmt;
+}
+
+int
+float_format_to_mopt_digit (oct_mach_info::float_format flt_fmt)
+{
+  int retval = -1;
+
+  switch (flt_fmt)
+    {
+    case oct_mach_info::flt_fmt_ieee_little_endian:
+      retval = 0;
+      break;
+
+    case oct_mach_info::flt_fmt_ieee_big_endian:
+      retval = 1;
+      break;
+
+    case oct_mach_info::flt_fmt_vax_d:
+      retval = 2;
+      break;
+
+    case oct_mach_info::flt_fmt_vax_g:
+      retval = 3;
+      break;
+
+    case oct_mach_info::flt_fmt_cray:
+      retval = 4;
+      break;
+
+    default:
+      break;
+    }
+
+  return retval;
+}
+
+// Extract one value (scalar, matrix, string, etc.) from stream IS and
+// place it in TC, returning the name of the variable.
+//
+// The data is expected to be in Matlab version 4 .mat format, though
+// not all the features of that format are supported.
+//
+// FILENAME is used for error messages.
+//
+// This format provides no way to tag the data as global.
+
+std::string
+read_mat_binary_data (std::istream& is, const std::string& filename,
+		      octave_value& tc)
+{
+  std::string retval;
+
+  // These are initialized here instead of closer to where they are
+  // first used to avoid errors from gcc about goto crossing
+  // initialization of variable.
+
+  Matrix re;
+  oct_mach_info::float_format flt_fmt = oct_mach_info::flt_fmt_unknown;
+  bool swap = false;
+  int type = 0;
+  int prec = 0;
+  int order = 0;
+  int mach = 0;
+  int dlen = 0;
+
+  FOUR_BYTE_INT mopt, nr, nc, imag, len;
+
+  int err = read_mat_file_header (is, swap, mopt, nr, nc, imag, len);
+  if (err)
+    {
+      if (err < 0)
+	goto data_read_error;
+      else
+	return retval;
+    }
+
+  type = mopt % 10;  // Full, sparse, etc.
+  mopt /= 10;        // Eliminate first digit.
+  prec = mopt % 10;  // double, float, int, etc.
+  mopt /= 10;        // Eliminate second digit.
+  order = mopt % 10; // Row or column major ordering.
+  mopt /= 10;        // Eliminate third digit.
+  mach = mopt % 10;  // IEEE, VAX, etc.
+
+  flt_fmt = mopt_digit_to_float_format (mach);
+
+  if (flt_fmt == oct_mach_info::flt_fmt_unknown)
+    {
+      error ("load: unrecognized binary format!");
+      return retval;
+    }
+
+  if (type != 0 && type != 1)
+    {
+      error ("load: can't read sparse matrices");
+      return retval;
+    }
+
+  if (imag && type == 1)
+    {
+      error ("load: encountered complex matrix with string flag set!");
+      return retval;
+    }
+
+  // LEN includes the terminating character, and the file is also
+  // supposed to include it, but apparently not all files do.  Either
+  // way, I think this should work.
+
+  {
+    OCTAVE_LOCAL_BUFFER (char, name, len+1);
+    name[len] = '\0';
+    if (! is.read (X_CAST (char *, name), len))
+      goto data_read_error;
+    retval = name;
+
+    dlen = nr * nc;
+    if (dlen < 0)
+      goto data_read_error;
+
+    if (order)
+      {
+	int tmp = nr;
+	nr = nc;
+	nc = tmp;
+      }
+
+      re.resize (nr, nc);
+
+      read_mat_binary_data (is, re.fortran_vec (), prec, dlen, swap, flt_fmt);
+
+      if (! is || error_state)
+	{
+	  error ("load: reading matrix data for `%s'", name);
+	  goto data_read_error;
+	}
+
+      if (imag)
+	{
+	  Matrix im (nr, nc);
+
+	  read_mat_binary_data (is, im.fortran_vec (), prec, dlen, swap,
+				flt_fmt);
+
+	  if (! is || error_state)
+	    {
+	      error ("load: reading imaginary matrix data for `%s'", name);
+	      goto data_read_error;
+	    }
+
+	  ComplexMatrix ctmp (nr, nc);
+
+	  for (int j = 0; j < nc; j++)
+	    for (int i = 0; i < nr; i++)
+	      ctmp (i, j) = Complex (re (i, j), im (i, j));
+
+	  tc = order ? ctmp.transpose () : ctmp;
+	}
+      else
+	tc = order ? re.transpose () : re;
+
+      if (type == 1)
+	tc = tc.convert_to_str ();
+
+      return retval;
+    }
+
+ data_read_error:
+  error ("load: trouble reading binary file `%s'", filename.c_str ());
+  return retval;
+}
+
+// Save the data from TC along with the corresponding NAME on stream OS 
+// in the MatLab version 4 binary format.
+
+bool
+save_mat_binary_data (std::ostream& os, const octave_value& tc,
+		      const std::string& name) 
+{
+  FOUR_BYTE_INT mopt = 0;
+
+  mopt += tc.is_string () ? 1 : 0;
+
+  oct_mach_info::float_format flt_fmt =
+    oct_mach_info::native_float_format ();;
+
+  mopt += 1000 * float_format_to_mopt_digit (flt_fmt);
+
+  os.write (X_CAST (char *, &mopt), 4);
+  
+  FOUR_BYTE_INT nr = tc.rows ();
+  os.write (X_CAST (char *, &nr), 4);
+
+  FOUR_BYTE_INT nc = tc.columns ();
+  os.write (X_CAST (char *, &nc), 4);
+
+  int len = nr * nc;
+
+  FOUR_BYTE_INT imag = tc.is_complex_type () ? 1 : 0;
+  os.write (X_CAST (char *, &imag), 4);
+
+  // LEN includes the terminating character, and the file is also
+  // supposed to include it.
+
+  FOUR_BYTE_INT name_len = name.length () + 1;
+
+  os.write (X_CAST (char *, &name_len), 4);
+  os << name << '\0';
+
+  if (tc.is_string ())
+    {
+      unwind_protect::begin_frame ("save_mat_binary_data");
+
+      charMatrix chm = tc.char_matrix_value ();
+
+      int nrow = chm.rows ();
+      int ncol = chm.cols ();
+	
+      OCTAVE_LOCAL_BUFFER (double, buf, ncol*nrow);
+	
+      for (int i = 0; i < nrow; i++)
+      	{
+	  std::string tstr = chm.row_as_string (i);
+	  const char *s = tstr.data ();
+	  
+	  for (int j = 0; j < ncol; j++)
+	    buf[j*nrow+i] = static_cast<double> (*s++ & 0x00FF);
+       	}
+      os.write ((char *)buf, nrow*ncol*sizeof(double));
+      
+      unwind_protect::run_frame ("save_mat_binary_data");
+    }
+  else if (tc.is_range ())
+    {
+      Range r = tc.range_value ();
+      double base = r.base ();
+      double inc = r.inc ();
+      int nel = r.nelem ();
+      for (int i = 0; i < nel; i++)
+	{
+	  double x = base + i * inc;
+	  os.write (X_CAST (char *, &x), 8);
+	}
+    }
+  else if (tc.is_real_scalar ())
+    {
+      double tmp = tc.double_value ();
+      os.write (X_CAST (char *, &tmp), 8);
+    }
+  else if (tc.is_real_matrix ())
+    {
+      Matrix m = tc.matrix_value ();
+      os.write (X_CAST (char *, m.data ()), 8 * len);
+    }
+  else if (tc.is_complex_scalar ())
+    {
+      Complex tmp = tc.complex_value ();
+      os.write (X_CAST (char *, &tmp), 16);
+    }
+  else if (tc.is_complex_matrix ())
+    {
+      ComplexMatrix m_cmplx = tc.complex_matrix_value ();
+      Matrix m = ::real (m_cmplx);
+      os.write (X_CAST (char *, m.data ()), 8 * len);
+      m = ::imag (m_cmplx);
+      os.write (X_CAST (char *, m.data ()), 8 * len);
+    }
+  else
+    gripe_wrong_type_arg ("save", tc, false);
+
+  return os;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
diff --git a/src/ls-mat4.h b/src/ls-mat4.h
new file mode 100644
--- /dev/null
+++ b/src/ls-mat4.h
@@ -0,0 +1,52 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_mat4_h)
+#define octave_ls_mat4 1
+
+extern oct_mach_info::float_format
+mopt_digit_to_float_format (int mach);
+
+extern int
+float_format_to_mopt_digit (oct_mach_info::float_format flt_fmt);
+
+extern int
+read_mat_file_header (std::istream& is, bool& swap, FOUR_BYTE_INT& mopt, 
+		      FOUR_BYTE_INT& nr, FOUR_BYTE_INT& nc,
+		      FOUR_BYTE_INT& imag, FOUR_BYTE_INT& len,
+		      int quiet = 0);
+
+extern std::string
+read_mat_binary_data (std::istream& is, const std::string& filename,
+		      octave_value& tc);
+
+extern bool
+save_mat_binary_data (std::ostream& os, const octave_value& tc,
+		      const std::string& name) ;
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
diff --git a/src/ls-mat5.cc b/src/ls-mat5.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-mat5.cc
@@ -0,0 +1,909 @@
+/*
+
+Copyright (C) 1996, 1997 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+// Author: James R. Van Zandt <jrv@vanzandt.mv.com>
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <cfloat>
+#include <cstring>
+#include <cctype>
+
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#include <memory>
+#include <string>
+
+#include "byte-swap.h"
+#include "data-conv.h"
+#include "file-ops.h"
+#include "glob-match.h"
+#include "lo-mappers.h"
+#include "lo-sstream.h"
+#include "mach-info.h"
+#include "oct-env.h"
+#include "oct-time.h"
+#include "quit.h"
+#include "str-vec.h"
+
+#include "Cell.h"
+#include "defun.h"
+#include "error.h"
+#include "gripes.h"
+#include "load-save.h"
+#include "oct-obj.h"
+#include "oct-map.h"
+#include "ov-cell.h"
+#include "pager.h"
+#include "pt-exp.h"
+#include "symtab.h"
+#include "sysdep.h"
+#include "unwind-prot.h"
+#include "utils.h"
+#include "variables.h"
+#include "version.h"
+#include "dMatrix.h"
+
+#include "ls-utils.h"
+#include "ls-mat5.h"
+
+#define PAD(l) (((l)<=4)?4:(((l)+7)/8)*8)
+#define TAGLENGTH(l) ((l)<=4?4:8)
+
+enum arrayclasstype
+  {
+    mxCELL_CLASS=1,		// cell array
+    mxSTRUCT_CLASS,		// structure
+    mxOBJECT_CLASS,		// object
+    mxCHAR_CLASS,		// character array
+    mxSPARSE_CLASS,		// sparse array
+    mxDOUBLE_CLASS,		// double precision array
+    mxSINGLE_CLASS,		// single precision floating point
+    mxINT8_CLASS,		// 8 bit signed integer
+    mxUINT8_CLASS,		// 8 bit unsigned integer
+    mxINT16_CLASS,		// 16 bit signed integer
+    mxUINT16_CLASS,		// 16 bit unsigned integer
+    mxINT32_CLASS,		// 32 bit signed integer
+    mxUINT32_CLASS		// 32 bit unsigned integer
+  };
+
+// Read COUNT elements of data from IS in the format specified by TYPE,
+// placing the result in DATA.  If SWAP is TRUE, swap the bytes of
+// each element before copying to DATA.  FLT_FMT specifies the format
+// of the data if we are reading floating point numbers.
+
+static void
+read_mat5_binary_data (std::istream& is, double *data,
+		       int count, bool swap, mat5_data_type type,
+		       oct_mach_info::float_format flt_fmt)
+{
+  
+  switch (type)
+    {
+    case miINT8:
+      read_doubles (is, data, LS_CHAR, count, swap, flt_fmt);
+      break;
+
+    case miUINT8:
+      read_doubles (is, data, LS_U_CHAR, count, swap, flt_fmt);
+      break;
+
+    case miINT16:
+      read_doubles (is, data, LS_SHORT, count, swap, flt_fmt);
+      break;
+
+    case miUINT16:
+      read_doubles (is, data, LS_U_SHORT, count, swap, flt_fmt);
+      break;
+
+    case miINT32:
+      read_doubles (is, data, LS_INT, count, swap, flt_fmt);
+      break;
+
+    case miUINT32:
+      read_doubles (is, data, LS_U_INT, count, swap, flt_fmt);
+      break;
+
+    case miSINGLE:
+      read_doubles (is, data, LS_FLOAT, count, swap, flt_fmt);
+      break;
+
+    case miRESERVE1:
+      break;
+
+    case miDOUBLE:
+      read_doubles (is, data, LS_DOUBLE, count, swap, flt_fmt);
+      break;
+
+    case miRESERVE2:
+    case miRESERVE3:
+      break;
+
+    case miINT64:
+#ifdef EIGHT_BYTE_INT
+      read_doubles (is, data, LS_LONG, count, swap, flt_fmt);
+#endif
+      break;
+
+    case miUINT64:
+#ifdef EIGHT_BYTE_INT
+      read_doubles (is, data, LS_U_LONG, count, swap, flt_fmt);
+#endif
+      break;
+
+    case miMATRIX:
+    default:
+      break;
+    }
+}
+
+// Read one element tag from stream IS, 
+// place the type code in TYPE and the byte count in BYTES
+// return nonzero on error
+static int
+read_mat5_tag (std::istream& is, bool swap, int& type, int& bytes)
+{
+  unsigned int upper;
+  FOUR_BYTE_INT temp;
+
+  if (! is.read (X_CAST (char *, &temp), 4 ))
+    goto data_read_error;
+
+  if (swap)
+    swap_4_bytes ((char *)&temp);
+
+  upper = (temp >> 16) & 0xffff;
+  type = temp & 0xffff;
+
+  if (upper)
+    {
+      // "compressed" format
+      bytes = upper;
+    }
+  else
+    {
+      if (! is.read (X_CAST (char *, &temp), 4 ))
+	goto data_read_error;
+      if (swap)
+	swap_4_bytes ((char *)&temp);
+      bytes = temp;
+    }
+
+  return 0;
+
+ data_read_error:
+  return 1;
+}
+
+// Extract one data element (scalar, matrix, string, etc.) from stream
+// IS and place it in TC, returning the name of the variable.
+//
+// The data is expected to be in Matlab's "Version 5" .mat format,
+// though not all the features of that format are supported.
+//
+// FILENAME is used for error messages.
+
+std::string
+read_mat5_binary_element (std::istream& is, const std::string& filename,
+			  bool swap, bool& global, octave_value& tc)
+{
+  std::string retval;
+
+  // These are initialized here instead of closer to where they are
+  // first used to avoid errors from gcc about goto crossing
+  // initialization of variable.
+
+  NDArray re;
+  oct_mach_info::float_format flt_fmt = oct_mach_info::flt_fmt_unknown;
+  int type = 0;
+  bool imag;
+  bool logicalvar;
+  enum arrayclasstype arrayclass;
+  FOUR_BYTE_INT junk;
+  FOUR_BYTE_INT flags;
+  FOUR_BYTE_INT dimension_length;
+  dim_vector dims;
+  int len;
+  int element_length;
+  std::streampos pos;
+  TWO_BYTE_INT number;
+  number = *(TWO_BYTE_INT *)"\x00\x01";
+
+  global = false;
+
+  // MAT files always use IEEE floating point
+  if ((number == 1) ^ swap)
+    flt_fmt = oct_mach_info::flt_fmt_ieee_big_endian;
+  else
+    flt_fmt = oct_mach_info::flt_fmt_ieee_little_endian;
+
+  // element type and length
+  if (read_mat5_tag (is, swap, type, element_length))
+    return retval;			// EOF
+
+  if (type != miMATRIX)
+    {
+      error ("load: invalid element type");
+      goto early_read_error;
+    }
+
+  if (element_length == 0)
+    {
+      tc = Matrix ();
+      return retval;
+    }
+
+  pos = is.tellg ();
+
+  // array flags subelement
+  if (read_mat5_tag (is, swap, type, len) || type != miUINT32 || len != 8)
+    {
+      error ("load: invalid array flags subelement");
+      goto early_read_error;
+    }
+
+  read_int (is, swap, flags);
+  imag = (flags & 0x0800) != 0;	// has an imaginary part?
+  global = (flags & 0x0400) != 0; // global variable?
+  logicalvar = (flags & 0x0200) != 0; // we don't use this yet
+  arrayclass = (arrayclasstype)(flags & 0xff);
+  read_int (is, swap, junk);	// an "undefined" entry
+  
+  // dimensions array subelement
+  {
+    if (read_mat5_tag (is, swap, type, dimension_length) || type != miINT32)
+      {
+	error ("load: invalid dimensions array subelement");
+	goto early_read_error;
+      }
+
+    int ndims = dimension_length / 4;
+    dims.resize (ndims);
+    for (int i = 0; i < ndims; i++)
+      {
+	FOUR_BYTE_INT n;
+	read_int (is, swap, n);
+	dims(i) = n;
+      }
+    re.resize (dims);
+
+    std::streampos tmp_pos = is.tellg ();
+    is.seekg (tmp_pos + static_cast<std::streamoff> (PAD (dimension_length) - dimension_length));
+  }
+
+  if (read_mat5_tag (is, swap, type, len) || type != miINT8)
+    {
+      error ("load: invalid array name subelement");
+      goto early_read_error;
+    }
+
+  {
+    OCTAVE_LOCAL_BUFFER (char, name, len+1);
+
+    // Structure field subelements have zero-length array name subelements.
+
+    std::streampos tmp_pos = is.tellg ();
+
+    if (len)
+      {
+	if (! is.read (X_CAST (char *, name), len ))
+	  goto data_read_error;
+	
+	is.seekg (tmp_pos + static_cast<std::streamoff> (PAD (len)));
+      }
+
+    name[len] = '\0';
+    retval = name;
+  }
+
+  switch (arrayclass)
+    {
+    case mxCELL_CLASS:
+      {
+	Cell cell_array (dims);
+
+	int n = cell_array.length ();
+
+	for (int i = 0; i < n; i++)
+	  {
+	    octave_value tc2;
+
+	    std::string nm
+	      = read_mat5_binary_element (is, filename, swap, global, tc2);
+
+	    if (! is || error_state)
+	      {
+		error ("load: reading cell data for `%s'", nm.c_str ());
+		goto data_read_error;
+	      }
+
+	    cell_array(i) = tc2;
+	  }
+
+	tc = cell_array;
+      }
+      break;
+
+    case mxOBJECT_CLASS:
+      warning ("load: objects are not implemented");
+      goto skip_ahead;
+
+    case mxSPARSE_CLASS:
+      warning ("load: sparse arrays are not implemented");
+      goto skip_ahead;
+
+    case mxSTRUCT_CLASS:
+      {
+	Octave_map m;
+	FOUR_BYTE_INT fn_type;
+	FOUR_BYTE_INT fn_len;
+	FOUR_BYTE_INT field_name_length;
+	int i;
+
+	// field name length subelement -- actually the maximum length
+	// of a field name.  The Matlab docs promise this will always
+	// be 32.  We read and use the actual value, on the theory
+	// that eventually someone will recognize that's a waste of
+	// space.
+	if (read_mat5_tag (is, swap, fn_type, fn_len) || fn_type != miINT32)
+	  {
+	    error ("load: invalid field name subelement");
+	    goto data_read_error;
+	  }
+
+	if (! is.read (X_CAST (char *, &field_name_length), fn_len ))
+	  goto data_read_error;
+
+	if (swap)
+	  swap_4_bytes ((char *)&field_name_length);
+
+	// field name subelement.  The length of this subelement tells
+	// us how many fields there are.
+	if (read_mat5_tag (is, swap, fn_type, fn_len) || fn_type != miINT8)
+	  {
+	    error ("load: invalid field name subelement");
+	    goto data_read_error;
+	  }
+
+	int n_fields = fn_len/field_name_length;
+
+	fn_len = PAD (fn_len);
+
+	OCTAVE_LOCAL_BUFFER (char, elname, fn_len);
+
+	if (! is.read (elname, fn_len))
+	  goto data_read_error;
+
+	int n;
+	if (dims(0) == 1)
+	  n = dims(1);
+	else if (dims(1) == 1)
+	  n = dims(0);
+	else
+	  {
+	    error ("load: can only handle one-dimensional structure arrays");
+	    goto data_read_error;
+	  }
+
+	Cell field_elts (n_fields, n);
+
+	// fields subelements
+	for (int j = 0; j < n; j++)
+	  {
+	    for (i = 0; i < n_fields; i++)
+	      {
+		octave_value fieldtc;
+		read_mat5_binary_element (is, filename, swap, global, fieldtc);
+		field_elts(i,j) = fieldtc;
+	      }
+	  }
+
+	for (int j = n_fields-1; j >= 0; j--)
+	  {
+	    const char *key = elname + j*field_name_length;
+
+	    for (int k = n-1; k >=0; k--)
+	      m[key](k) = field_elts(j,k);
+	  }
+
+	tc = m;
+      }
+      break;
+
+    case mxCHAR_CLASS:
+      // handle as a numerical array to start with
+
+    case mxDOUBLE_CLASS:
+    case mxSINGLE_CLASS:
+    case mxINT8_CLASS:
+    case mxUINT8_CLASS:
+    case mxINT16_CLASS:
+    case mxUINT16_CLASS:
+    case mxINT32_CLASS:
+    case mxUINT32_CLASS:
+    default:
+      // handle all these numerical formats as double arrays
+      
+      // real data subelement
+      {
+	std::streampos tmp_pos;
+	
+	if (read_mat5_tag (is, swap, type, len))
+	  {
+	    error ("load: reading matrix data for `%s'", retval.c_str ());
+	    goto data_read_error;
+	  }
+
+	int n = re.length ();
+	tmp_pos = is.tellg ();
+	read_mat5_binary_data (is, re.fortran_vec (), n, swap,
+			       (enum mat5_data_type) type, flt_fmt);
+
+	if (! is || error_state)
+	  {
+	    error ("load: reading matrix data for `%s'", retval.c_str ());
+	    goto data_read_error;
+	  }
+
+	is.seekg (tmp_pos + static_cast<std::streamoff> (PAD (len)));
+      }
+      
+      // imaginary data subelement
+      if (imag)
+	{
+	  NDArray im (dims);
+	  
+	  if (read_mat5_tag (is, swap, type, len))
+	    {
+	      error ("load: reading matrix data for `%s'", retval.c_str ());
+	      goto data_read_error;
+	    }
+
+	  int n = im.length ();
+	  read_mat5_binary_data (is, im.fortran_vec (), n, swap,
+				 (enum mat5_data_type) type, flt_fmt);
+
+	  if (! is || error_state)
+	    {
+	      error ("load: reading imaginary matrix data for `%s'",
+		     retval.c_str ());
+	      goto data_read_error;
+	    }
+
+	  ComplexNDArray ctmp (dims);
+
+	  for (int i = 0; i < n; i++)
+	    ctmp(i) = Complex (re(i), im(i));
+
+	  tc = ctmp;
+	}
+      else
+	tc = re;
+
+      if (arrayclass == mxCHAR_CLASS)
+	tc = tc.convert_to_str ();
+    }
+
+  is.seekg (pos + static_cast<std::streamoff> (element_length));
+
+  if (is.eof ())
+    is.clear ();
+
+  return retval;
+
+ data_read_error:
+ early_read_error:
+  error ("load: trouble reading binary file `%s'", filename.c_str ());
+  return std::string ();
+
+ skip_ahead:
+  warning ("skipping over `%s'", retval.c_str ());
+  is.seekg (pos + static_cast<std::streamoff> (element_length));
+  return read_mat5_binary_element (is, filename, swap, global, tc);
+}
+
+int
+read_mat5_binary_file_header (std::istream& is, bool& swap, bool quiet)
+{
+  TWO_BYTE_INT version=0, magic=0;
+
+  is.seekg (124, std::ios::beg);
+  is.read (X_CAST (char *, &version), 2);
+  is.read (X_CAST (char *, &magic), 2);
+
+  if (magic == 0x4d49)
+    swap = 0;
+  else if (magic == 0x494d)
+    swap = 1;
+  else
+    {
+      if (! quiet)
+	error ("load: can't read binary file");
+      return -1;
+    }
+
+  if (! swap)			// version number is inverse swapped!
+    version = ((version >> 8) & 0xff) + ((version & 0xff) << 8);
+
+  if (version != 1 && !quiet)
+    warning ("load: found version %d binary MAT file, "
+	     "but only prepared for version 1", version);
+
+  return 0;
+}
+
+static int 
+write_mat5_tag (std::ostream& is, int type, int bytes)
+{
+  FOUR_BYTE_INT temp;
+
+  if (bytes <= 4)
+    temp = (bytes << 16) + type;
+  else
+    {
+      temp = type;
+      if (! is.write ((char *)&temp, 4))
+	goto data_write_error;
+      temp = bytes;
+    }
+
+  if (! is.write ((char *)&temp, 4))
+    goto data_write_error;
+
+  return 0;
+
+ data_write_error:
+  return 1;
+}
+
+// write out the numeric values in M to OS,
+// preceded by the appropriate tag.
+static void 
+write_mat5_array (std::ostream& os, const NDArray& m, bool save_as_floats)
+{
+  int nel = m.nelem ();
+  double max_val, min_val;
+  save_type st = LS_DOUBLE;
+  mat5_data_type mst;
+  int size;
+  unsigned len;
+  const double *data = m.data ();
+
+// Have to use copy here to avoid writing over data accessed via
+// Matrix::data().
+
+#define MAT5_DO_WRITE(TYPE, data, count, stream)			\
+  do									\
+    {									\
+      OCTAVE_LOCAL_BUFFER (TYPE, ptr, count);                           \
+      for (int i = 0; i < count; i++)					\
+        ptr[i] = X_CAST (TYPE, data[i]);				\
+      stream.write (X_CAST (char *, ptr), count * sizeof (TYPE));	\
+    }									\
+  while (0)
+
+  if (save_as_floats)
+    {
+      if (m.too_large_for_float ())
+	{
+	  warning ("save: some values too large to save as floats --");
+	  warning ("save: saving as doubles instead");
+	}
+      else
+	st = LS_FLOAT;
+    }
+
+  if (m.all_integers (max_val, min_val))
+    st = get_save_type (max_val, min_val);
+
+  switch (st)
+    {
+    default:
+    case LS_DOUBLE:  mst = miDOUBLE; size = 8; break;
+    case LS_FLOAT:   mst = miSINGLE; size = 4; break;
+    case LS_U_CHAR:  mst = miUINT8;  size = 1; break;
+    case LS_U_SHORT: mst = miUINT16; size = 2; break;
+    case LS_U_INT:   mst = miUINT32; size = 4; break;
+    case LS_CHAR:    mst = miINT8;   size = 1; break;
+    case LS_SHORT:   mst = miINT16;  size = 2; break;
+    case LS_INT:     mst = miINT32;  size = 4; break;
+    }
+
+  len = nel*size;
+  write_mat5_tag (os, mst, len);
+
+  {
+    switch (st)
+      {
+      case LS_U_CHAR:
+	MAT5_DO_WRITE (unsigned char, data, nel, os);
+	break;
+	
+      case LS_U_SHORT:
+	MAT5_DO_WRITE (unsigned TWO_BYTE_INT, data, nel, os);
+	break;
+	
+      case LS_U_INT:
+	MAT5_DO_WRITE (unsigned FOUR_BYTE_INT, data, nel, os);
+	break;
+	
+	// provide for 64 bit ints, even though get_save_type does
+	// not yet implement them
+#ifdef EIGHT_BYTE_INT
+      case LS_U_LONG:
+	MAT5_DO_WRITE (unsigned EIGHT_BYTE_INT, data, nel, os);
+	break;
+#endif
+
+      case LS_CHAR:
+	MAT5_DO_WRITE (signed char, data, nel, os);
+	break;
+	
+      case LS_SHORT:
+	MAT5_DO_WRITE (TWO_BYTE_INT, data, nel, os);
+	break;
+
+      case LS_INT:
+	MAT5_DO_WRITE (FOUR_BYTE_INT, data, nel, os);
+	break;
+
+#ifdef EIGHT_BYTE_INT
+      case LS_LONG:
+	MAT5_DO_WRITE (EIGHT_BYTE_INT, data, nel, os);
+	break;
+#endif
+
+      case LS_FLOAT:
+	MAT5_DO_WRITE (float, data, nel, os);
+	break;
+
+      case LS_DOUBLE: // No conversion necessary.
+	os.write (X_CAST (char *, data), len);
+	break;
+
+      default:
+	(*current_liboctave_error_handler)
+	  ("unrecognized data format requested");
+	break;
+      }
+  }
+  if (PAD (len) > len)
+    {
+      static char buf[9]="\x00\x00\x00\x00\x00\x00\x00\x00";
+      os.write (buf, PAD (len) - len);
+    }
+}
+
+// Write out cell element values in the cell array to OS, preceded by
+// the appropriate tag.
+
+static bool 
+write_mat5_cell_array (std::ostream& os, Cell& cell, bool mark_as_global,
+		       bool save_as_floats)
+{
+  int nel = cell.nelem ();
+
+  for (int i = 0; i < nel; i++)
+    {
+      octave_value ov = cell(i);
+
+      if (! save_mat5_binary_element (os, ov, "", mark_as_global,
+				      save_as_floats))
+	return false;
+    }
+
+  return true;
+}
+
+// save the data from TC along with the corresponding NAME on stream
+// OS in the MatLab version 5 binary format.  Return true on success.
+
+bool
+save_mat5_binary_element (std::ostream& os,
+			  const octave_value& tc, const std::string& name,
+			  bool mark_as_global, bool save_as_floats) 
+{
+  FOUR_BYTE_INT flags=0;
+  FOUR_BYTE_INT junk=0;
+  std::streampos fixup, contin;
+
+  // element type and length
+  fixup = os.tellp ();
+  write_mat5_tag (os, miMATRIX, 99); // we don't know the real length yet
+  
+  // array flags subelement
+  write_mat5_tag (os, miUINT32, 8);
+
+  if (mark_as_global)
+    flags |= 0x0400;
+
+  if (tc.is_complex_scalar () || tc.is_complex_matrix ())
+    flags |= 0x0800;
+
+  if (tc.is_string ())
+    flags |= mxCHAR_CLASS;
+  else if (tc.is_real_scalar ())
+    flags |= mxDOUBLE_CLASS;
+  else if (tc.is_real_matrix () || tc.is_range ())
+    flags |= mxDOUBLE_CLASS;
+  else if (tc.is_complex_scalar ())
+    flags |= mxDOUBLE_CLASS;
+  else if (tc.is_complex_matrix ())
+    flags |= mxDOUBLE_CLASS;
+  else if (tc.is_map ()) 
+    flags |= mxSTRUCT_CLASS;
+  else if (tc.is_cell ())
+    flags |= mxCELL_CLASS;
+  else
+    {
+      gripe_wrong_type_arg ("save", tc, false);
+      goto error_cleanup;
+    }
+
+  os.write ((char *)&flags, 4);
+  os.write ((char *)&junk, 4);
+
+  {
+    dim_vector dv = tc.dims ();
+    int nd = tc.ndims ();
+
+    write_mat5_tag (os, miINT32, 8);
+
+    for (int i = 0; i < nd; i++)
+      {
+	int n = dv(i);
+	os.write ((char *)&n, 4);
+      }
+  }
+
+  // array name subelement
+  {
+    int namelen = name.length ();
+
+    if (namelen > 31)
+      namelen = 31; // only 31 char names permitted in mat file
+
+    int paddedlength = PAD (namelen);
+
+    write_mat5_tag (os, miINT8, namelen);
+    OCTAVE_LOCAL_BUFFER (char, paddedname, paddedlength);
+    memset (paddedname, 0, paddedlength);
+    strncpy (paddedname, name.c_str (), namelen);
+    os.write (paddedname, paddedlength);
+  }
+
+  // data element
+  if (tc.is_string ())
+    {
+      charMatrix chm = tc.char_matrix_value ();
+      int nr = chm.rows ();
+      int nc = chm.cols ();
+      int len = nr*nc*2;
+      int paddedlength = PAD (nr*nc*2);
+
+      OCTAVE_LOCAL_BUFFER (TWO_BYTE_INT, buf, nc*nr+3);
+      write_mat5_tag (os, miUINT16, len);
+
+      for (int i = 0; i < nr; i++)
+	{
+	  std::string tstr = chm.row_as_string (i);
+	  const char *s = tstr.data ();
+
+	  for (int j = 0; j < nc; j++)
+	    buf[j*nr+i] = *s++ & 0x00FF;
+	}
+      os.write ((char *)buf, nr*nc*2);
+      
+      if (paddedlength > len)
+	os.write ((char *)buf, paddedlength - len);
+    }
+  else if (tc.is_real_scalar () || tc.is_real_matrix () || tc.is_range ())
+    {
+      NDArray m = tc.array_value ();
+
+      write_mat5_array (os, m, save_as_floats);
+    }
+  else if (tc.is_cell ())
+    {
+      Cell cell = tc.cell_value ();
+
+      if (! write_mat5_cell_array (os, cell, mark_as_global, save_as_floats))
+	goto error_cleanup;
+    }
+  else if (tc.is_complex_scalar () || tc.is_complex_matrix ()) 
+    {
+      ComplexNDArray m_cmplx = tc.complex_matrix_value ();
+
+      write_mat5_array (os, ::real (m_cmplx), save_as_floats);
+      write_mat5_array (os, ::imag (m_cmplx), save_as_floats);
+    }
+  else if (tc.is_map ()) 
+    {
+      // an Octave structure */
+      // recursively write each element of the structure
+      Octave_map m = tc.map_value ();
+
+      {
+	char buf[32];
+	FOUR_BYTE_INT maxfieldnamelength = 32;
+	int fieldcnt = 0;
+
+	for (Octave_map::iterator i = m.begin (); i != m.end (); i++)
+	  fieldcnt++;
+
+	write_mat5_tag (os, miINT32, 4);
+	os.write ((char *)&maxfieldnamelength, 4);
+	write_mat5_tag (os, miINT8, fieldcnt*32);
+
+	for (Octave_map::iterator i = m.begin (); i != m.end (); i++)
+	  {
+	    // write the name of each element
+	    std::string tstr = m.key (i);
+	    memset (buf, 0, 32);
+	    strncpy (buf, tstr.c_str (), 31); // only 31 char names permitted
+	    os.write (buf, 32);
+	  }
+
+	int len = m.numel ();
+
+	for (Octave_map::iterator i = m.begin (); i != m.end (); i++)
+	  {
+	    // write the data of each element
+	    Cell elts = m.contents (i);
+
+	    for (int j = 0; j < len; j++)
+	      {
+		bool retval2 = save_mat5_binary_element (os, elts(j), "",
+							 mark_as_global,
+							 save_as_floats);
+		if (! retval2)
+		  goto error_cleanup;
+	      }
+	  }
+      }
+    }
+  else
+    gripe_wrong_type_arg ("save", tc, false);
+
+  contin = os.tellp ();
+  os.seekp (fixup);
+  write_mat5_tag (os, miMATRIX, 
+                  static_cast<int>(contin - fixup) - 8); // the actual length
+  os.seekp (contin);
+
+  return true;
+
+ error_cleanup:
+  error ("save: error while writing `%s' to MAT file", name.c_str ());
+
+  return false;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-mat5.h b/src/ls-mat5.h
new file mode 100644
--- /dev/null
+++ b/src/ls-mat5.h
@@ -0,0 +1,62 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_mat5_h)
+#define octave_ls_mat5_h 1
+
+enum mat5_data_type
+  {
+    miINT8 = 1,			// 8 bit signed
+    miUINT8,			// 8 bit unsigned
+    miINT16,			// 16 bit signed
+    miUINT16,			// 16 bit unsigned
+    miINT32,			// 32 bit signed
+    miUINT32,			// 32 bit unsigned
+    miSINGLE,			// IEEE 754 single precision float
+    miRESERVE1,
+    miDOUBLE,			// IEEE 754 double precision float
+    miRESERVE2,
+    miRESERVE3,
+    miINT64,			// 64 bit signed
+    miUINT64,			// 64 bit unsigned
+    miMATRIX			// MATLAB array
+  };
+
+extern int
+read_mat5_binary_file_header (std::istream& is, bool& swap,
+			      bool quiet = false);
+extern std::string
+read_mat5_binary_element (std::istream& is, const std::string& filename,
+			  bool swap, bool& global, octave_value& tc);
+extern bool
+save_mat5_binary_element (std::ostream& os,
+			  const octave_value& tc, const std::string& name,
+			  bool mark_as_global, bool save_as_floats);
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-oct-ascii.cc b/src/ls-oct-ascii.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-oct-ascii.cc
@@ -0,0 +1,972 @@
+/*
+
+Copyright (C) 1996, 1997 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+// Author: John W. Eaton.
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <cfloat>
+#include <cstring>
+#include <cctype>
+
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#include <memory>
+#include <string>
+
+#include "byte-swap.h"
+#include "data-conv.h"
+#include "file-ops.h"
+#include "glob-match.h"
+#include "lo-mappers.h"
+#include "lo-sstream.h"
+#include "mach-info.h"
+#include "oct-env.h"
+#include "oct-time.h"
+#include "quit.h"
+#include "str-vec.h"
+
+#include "Cell.h"
+#include "defun.h"
+#include "error.h"
+#include "gripes.h"
+#include "load-save.h"
+#include "oct-obj.h"
+#include "oct-map.h"
+#include "ov-cell.h"
+#include "pager.h"
+#include "pt-exp.h"
+#include "symtab.h"
+#include "sysdep.h"
+#include "unwind-prot.h"
+#include "utils.h"
+#include "variables.h"
+#include "version.h"
+#include "dMatrix.h"
+
+#include "ls-oct-ascii.h"
+
+// The number of decimal digits to use when writing ascii data.
+static int Vsave_precision;
+
+#define CELL_ELT_TAG "<cell-element>"
+
+// Used when converting Inf to something that gnuplot can read.
+
+#ifndef OCT_RBV
+#define OCT_RBV DBL_MAX / 100.0
+#endif
+
+// Functions for reading ascii data.
+
+static void
+ascii_save_type (std::ostream& os, const char *type, bool mark_as_global)
+{
+  if (mark_as_global)
+    os << "# type: global ";
+  else
+    os << "# type: ";
+
+  os << type << "\n";
+}
+
+static Matrix
+strip_infnan (const Matrix& m)
+{
+  int nr = m.rows ();
+  int nc = m.columns ();
+
+  Matrix retval (nr, nc);
+
+  int k = 0;
+  for (int i = 0; i < nr; i++)
+    {
+      for (int j = 0; j < nc; j++)
+	{
+	  double d = m (i, j);
+	  if (xisnan (d))
+	    goto next_row;
+	  else
+	    retval (k, j) = xisinf (d) ? (d > 0 ? OCT_RBV : -OCT_RBV) : d;
+	}
+      k++;
+
+    next_row:
+      continue;
+    }
+
+  if (k > 0)
+    retval.resize (k, nc);
+
+  return retval;
+}
+
+static ComplexMatrix
+strip_infnan (const ComplexMatrix& m)
+{
+  int nr = m.rows ();
+  int nc = m.columns ();
+
+  ComplexMatrix retval (nr, nc);
+
+  int k = 0;
+  for (int i = 0; i < nr; i++)
+    {
+      for (int j = 0; j < nc; j++)
+	{
+	  Complex c = m (i, j);
+	  if (xisnan (c))
+	    goto next_row;
+	  else
+	    {
+	      double re = real (c);
+	      double im = imag (c);
+
+	      re = xisinf (re) ? (re > 0 ? OCT_RBV : -OCT_RBV) : re;
+	      im = xisinf (im) ? (im > 0 ? OCT_RBV : -OCT_RBV) : im;
+
+	      retval (k, j) = Complex (re, im);
+	    }
+	}
+      k++;
+
+    next_row:
+      continue;
+    }
+
+  if (k > 0)
+    retval.resize (k, nc);
+
+  return retval;
+}
+
+// Skip white space and comments on stream IS.
+
+static void
+skip_comments (std::istream& is)
+{
+  char c = '\0';
+  while (is.get (c))
+    {
+      if (c == ' ' || c == '\t' || c == '\n')
+	; // Skip whitespace on way to beginning of next line.
+      else
+	break;
+    }
+
+  for (;;)
+    {
+      if (is && (c == '%' || c == '#'))
+	while (is.get (c) && c != '\n')
+	  ; // Skip to beginning of next line, ignoring everything.
+      else
+	break;
+    }
+}
+
+// Extract a KEYWORD and its value from stream IS, returning the
+// associated value in a new string.
+//
+// Input should look something like:
+//
+//  [%#][ \t]*keyword[ \t]*:[ \t]*string-value[ \t]*\n
+
+std::string
+extract_keyword (std::istream& is, const char *keyword)
+{
+  std::string retval;
+
+  char c;
+  while (is.get (c))
+    {
+      if (c == '%' || c == '#')
+	{
+	  OSSTREAM buf;
+	
+	  while (is.get (c) && (c == ' ' || c == '\t' || c == '%' || c == '#'))
+	    ; // Skip whitespace and comment characters.
+
+	  if (isalpha (c))
+	    buf << c;
+
+	  while (is.get (c) && isalpha (c))
+	    buf << c;
+
+	  buf << OSSTREAM_ENDS;
+	  const char *tmp = OSSTREAM_C_STR (buf);
+	  OSSTREAM_FREEZE (buf);
+	  int match = (strncmp (tmp, keyword, strlen (keyword)) == 0);
+
+	  if (match)
+	    {
+	      OSSTREAM value;
+	      while (is.get (c) && (c == ' ' || c == '\t' || c == ':'))
+		; // Skip whitespace and the colon.
+
+	      if (c != '\n')
+		{
+		  value << c;
+		  while (is.get (c) && c != '\n')
+		    value << c;
+		}
+	      value << OSSTREAM_ENDS;
+	      retval = OSSTREAM_STR (value);
+	      OSSTREAM_FREEZE (value);
+	      break;
+	    }
+	}
+    }
+
+  int len = retval.length ();
+
+  if (len > 0)
+    {
+      while (len)
+	{
+	  c = retval[len-1];
+
+	  if (c == ' ' || c == '\t')
+	    len--;
+	  else
+	    {
+	      retval.resize (len);
+	      break;
+	    }
+	}
+    }
+
+  return retval;
+}
+
+// Match KEYWORD on stream IS, placing the associated value in VALUE,
+// returning TRUE if successful and FALSE otherwise.
+//
+// Input should look something like:
+//
+//  [%#][ \t]*keyword[ \t]*int-value.*\n
+
+static bool
+extract_keyword (std::istream& is, const char *keyword, int& value)
+{
+  bool status = false;
+  value = 0;
+
+  char c;
+  while (is.get (c))
+    {
+      if (c == '%' || c == '#')
+	{
+	  OSSTREAM buf;
+
+	  while (is.get (c) && (c == ' ' || c == '\t' || c == '%' || c == '#'))
+	    ; // Skip whitespace and comment characters.
+
+	  if (isalpha (c))
+	    buf << c;
+
+	  while (is.get (c) && isalpha (c))
+	    buf << c;
+
+	  buf << OSSTREAM_ENDS;
+	  const char *tmp = OSSTREAM_C_STR (buf);
+	  int match = (strncmp (tmp, keyword, strlen (keyword)) == 0);
+	  OSSTREAM_FREEZE (buf);
+
+	  if (match)
+	    {
+	      while (is.get (c) && (c == ' ' || c == '\t' || c == ':'))
+		; // Skip whitespace and the colon.
+
+	      is.putback (c);
+	      if (c != '\n')
+		is >> value;
+	      if (is)
+		status = true;
+	      while (is.get (c) && c != '\n')
+		; // Skip to beginning of next line;
+	      break;
+	    }
+	}
+    }
+  return status;
+}
+
+// Extract one value (scalar, matrix, string, etc.) from stream IS and
+// place it in TC, returning the name of the variable.  If the value
+// is tagged as global in the file, return TRUE in GLOBAL.
+//
+// FILENAME is used for error messages.
+//
+// The data is expected to be in the following format:
+//
+// The input file must have a header followed by some data.
+//
+// All lines in the header must begin with a `#' character.
+//
+// The header must contain a list of keyword and value pairs with the
+// keyword and value separated by a colon.
+//
+// Keywords must appear in the following order:
+//
+// # name: <name>
+// # type: <type>
+// # <info>
+//
+// Where:
+//
+//  <name> : a valid identifier
+//
+//  <type> : <typename>
+//         | global <typename>
+//
+//  <typename> : scalar
+//             | complex scalar
+//             | matrix
+//             | complex matrix
+//             | string
+//             | range
+//             | string array
+//
+//  <info> : <matrix info>
+//         | <string info>
+//         | <string array info>
+//
+//  <matrix info> : # rows: <integer>
+//                : # columns: <integer>
+//
+//  <string info> : # length: <integer>
+//
+//  <string array info> : # elements: <integer>
+//                      : # length: <integer> (once before each string)
+//
+// Formatted ASCII data follows the header.
+//
+// Example:
+//
+//  # name: foo
+//  # type: matrix
+//  # rows: 2
+//  # columns: 2
+//    2  4
+//    1  3
+//
+// Example:
+//
+//  # name: foo
+//  # type: string array
+//  # elements: 5
+//  # length: 4
+//  this
+//  # length: 2
+//  is
+//  # length: 1
+//  a
+//  # length: 6
+//  string
+//  # length: 5
+//  array
+//
+// XXX FIXME XXX -- this format is fairly rigid, and doesn't allow for
+// arbitrary comments, etc.  Someone should fix that.
+
+// Ugh.  The signature of the compare method is not standard in older
+// versions of the GNU libstdc++.  Do this instead:
+
+#define SUBSTRING_COMPARE_EQ(s, pos, n, t) (s.substr (pos, n) == t)
+
+std::string
+read_ascii_data (std::istream& is, const std::string& filename, bool& global,
+		 octave_value& tc, int count)
+{
+  // Read name for this entry or break on EOF.
+
+  std::string name = extract_keyword (is, "name");
+
+  if (name.empty ())
+    {
+      if (count == 0)
+	error ("load: empty name keyword or no data found in file `%s'",
+	       filename.c_str ());
+
+      return std::string ();
+    }
+
+  if (name == CELL_ELT_TAG)
+    {
+      // This is OK -- name won't be used.
+    }
+  else if (! valid_identifier (name))
+    {
+      error ("load: bogus identifier `%s' found in file `%s'",
+	     name.c_str (), filename.c_str ());
+      return std::string ();
+    }
+
+  // Look for type keyword.
+
+  std::string tag = extract_keyword (is, "type");
+
+  if (! tag.empty ())
+    {
+      std::string typ;
+      size_t pos = tag.rfind (' ');
+
+      if (pos != NPOS)
+	{
+	  global = SUBSTRING_COMPARE_EQ (tag, 0, 6, "global");
+
+	  typ = global ? tag.substr (7) : tag;
+	}
+      else
+	typ = tag;
+
+      if (SUBSTRING_COMPARE_EQ (typ, 0, 6, "scalar"))
+	{
+	  double tmp = octave_read_double (is);
+	  if (is)
+	    tc = tmp;
+	  else
+	    error ("load: failed to load scalar constant");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 6, "matrix"))
+	{
+	  int nr = 0;
+	  int nc = 0;
+
+	  if (extract_keyword (is, "rows", nr) && nr >= 0
+	      && extract_keyword (is, "columns", nc) && nc >= 0)
+	    {
+	      if (nr > 0 && nc > 0)
+		{
+		  Matrix tmp (nr, nc);
+		  is >> tmp;
+		  if (is)
+		    tc = tmp;
+		  else
+		    error ("load: failed to load matrix constant");
+		}
+	      else if (nr == 0 || nc == 0)
+		tc = Matrix (nr, nc);
+	      else
+		panic_impossible ();
+	    }
+	  else
+	    error ("load: failed to extract number of rows and columns");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 4, "cell"))
+	{
+	  int nr = 0;
+	  int nc = 0;
+
+	  if (extract_keyword (is, "rows",    nr) && nr >= 0
+	      && extract_keyword (is, "columns", nc) && nc >= 0)
+	    {
+	      if (nr > 0 && nc > 0)
+		{
+		  Cell tmp (nr, nc);
+
+		  for (int j = 0; j < nc; j++)
+		    {
+		      for (int i = 0; i < nr; i++)
+			{
+			  octave_value t2;
+
+			  // recurse to read cell elements
+			  std::string nm
+			    = read_ascii_data (is, filename, global, t2, count);
+
+			  if (nm == CELL_ELT_TAG)
+			    {
+			      if (is)
+				tmp.elem (i, j) = t2;
+			    }
+			  else
+			    {
+			      error ("load: cell array element had unexpected name");
+			      goto cell_read_error;
+			    }
+			}
+		    }
+
+		cell_read_error:
+
+		  if (is)
+		    tc = tmp;
+		  else
+		    error ("load: failed to load cell element");
+		}
+	      else if (nr == 0 || nc == 0)
+		tc = Cell (nr, nc);
+	      else
+		panic_impossible ();
+	    }
+	  else
+	    error ("load: failed to extract number of rows and columns for cell array");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 14, "complex scalar"))
+	{
+	  Complex tmp = octave_read_complex (is);
+	  if (is)
+	    tc = tmp;
+	  else
+	    error ("load: failed to load complex scalar constant");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 14, "complex matrix"))
+	{
+	  int nr = 0;
+	  int nc = 0;
+
+	  if (extract_keyword (is, "rows", nr) && nr > 0
+	      && extract_keyword (is, "columns", nc) && nc > 0)
+	    {
+	      ComplexMatrix tmp (nr, nc);
+	      is >> tmp;
+	      if (is)
+		tc = tmp;
+	      else
+		error ("load: failed to load complex matrix constant");
+	    }
+	  else
+	    error ("load: failed to extract number of rows and columns");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 12, "string array"))
+	{
+	  int elements;
+	  if (extract_keyword (is, "elements", elements) && elements >= 0)
+	    {
+	      // XXX FIXME XXX -- need to be able to get max length
+	      // before doing anything.
+
+	      charMatrix chm (elements, 0);
+	      int max_len = 0;
+	      for (int i = 0; i < elements; i++)
+		{
+		  int len;
+		  if (extract_keyword (is, "length", len) && len >= 0)
+		    {
+		      OCTAVE_LOCAL_BUFFER (char, tmp, len+1);
+
+		      if (len > 0 && ! is.read (X_CAST (char *, tmp), len))
+			{
+			  error ("load: failed to load string constant");
+			  break;
+			}
+		      else
+			{
+			  tmp [len] = '\0';
+			  if (len > max_len)
+			    {
+			      max_len = len;
+			      chm.resize (elements, max_len, 0);
+			    }
+			  chm.insert (tmp, i, 0);
+			}
+		    }
+		  else
+		    error ("load: failed to extract string length for element %d", i+1);
+		}
+
+	      if (! error_state)
+		tc = octave_value (chm, true);
+	    }
+	  else
+	    error ("load: failed to extract number of string elements");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 6, "string"))
+	{
+	  int len;
+	  if (extract_keyword (is, "length", len) && len >= 0)
+	    {
+	      OCTAVE_LOCAL_BUFFER (char, tmp, len+1);
+
+	      if (len > 0 && ! is.read (X_CAST (char *, tmp), len))
+		{
+		  error ("load: failed to load string constant");
+		}
+	      else
+		{
+		  tmp [len] = '\0';
+
+		  if (is)
+		    tc = tmp;
+		  else
+		    error ("load: failed to load string constant");
+		}
+	    }
+	  else
+	    error ("load: failed to extract string length");
+	}
+      else if (SUBSTRING_COMPARE_EQ (typ, 0, 5, "range"))
+	{
+	  // # base, limit, range comment added by save ().
+
+	  skip_comments (is);
+	  Range tmp;
+	  is >> tmp;
+	  if (is)
+	    tc = tmp;
+	  else
+	    error ("load: failed to load range constant");
+	}
+      else
+	error ("load: unknown constant type `%s'", tag.c_str ());
+    }
+  else
+    error ("load: failed to extract keyword specifying value type");
+
+  if (error_state)
+    {
+      error ("load: reading file %s", filename.c_str ());
+      return std::string ();
+    }
+
+  return name;
+}
+
+// Save the data from TC along with the corresponding NAME, and global
+// flag MARK_AS_GLOBAL on stream OS in the plain text format described
+// above for load_ascii_data.  If NAME is empty, the name: line is not
+// generated.  PRECISION specifies the number of decimal digits to print. 
+// If STRIP_NAN_AND_INF is TRUE, rows containing NaNs are deleted,
+// and Infinite values are converted to +/-OCT_RBV (A Real Big Value,
+// but not so big that gnuplot can't handle it when trying to compute
+// axis ranges, etc.).
+//
+// Assumes ranges and strings cannot contain Inf or NaN values.
+//
+// Returns 1 for success and 0 for failure.
+
+// XXX FIXME XXX -- should probably write the help string here too.
+
+bool
+save_ascii_data (std::ostream& os, const octave_value& val_arg,
+		 const std::string& name, bool& infnan_warned,
+		 bool strip_nan_and_inf, bool mark_as_global,
+		 int precision)
+{
+  bool success = true;
+
+  if (! precision)
+    precision = Vsave_precision;
+
+  if (! name.empty ())
+    os << "# name: " << name << "\n";
+
+  long old_precision = os.precision ();
+  os.precision (precision);
+
+  octave_value val = val_arg;
+
+  if (val.is_range ())
+    {
+      Range r = val.range_value ();
+      double base = r.base ();
+      double limit = r.limit ();
+      double inc = r.inc ();
+      if (! (NINT (base) == base
+	     && NINT (limit) == limit
+	     && NINT (inc) == inc))
+	val = val.matrix_value ();
+    }	
+
+  if (val.is_string ())
+    {
+      ascii_save_type (os, "string array", mark_as_global);
+      charMatrix chm = val.char_matrix_value ();
+      int elements = chm.rows ();
+      os << "# elements: " << elements << "\n";
+      for (int i = 0; i < elements; i++)
+	{
+	  unsigned len = chm.cols ();
+	  os << "# length: " << len << "\n";
+	  std::string tstr = chm.row_as_string (i, false, true);
+	  const char *tmp = tstr.data ();
+	  if (tstr.length () > len)
+	    panic_impossible ();
+	  os.write (X_CAST (char *, tmp), len);
+	  os << "\n";
+	}
+    }
+  else if (val.is_range ())
+    {
+      ascii_save_type (os, "range", mark_as_global);
+      Range tmp = val.range_value ();
+      os << "# base, limit, increment\n";
+      octave_write_double (os, tmp.base ());
+      os << " ";
+      octave_write_double (os, tmp.limit ());
+      os << " ";
+      octave_write_double (os, tmp.inc ());
+      os << "\n";
+    }
+  else if (val.is_real_scalar ())
+    {
+      ascii_save_type (os, "scalar", mark_as_global);
+
+      double d = val.double_value ();
+
+      if (strip_nan_and_inf)
+	{
+	  if (xisnan (d))
+	    {
+	      error ("only value to plot is NaN");
+	      success = false;
+	    }
+	  else
+	    {
+	      d = xisinf (d) ? (d > 0 ? OCT_RBV : -OCT_RBV) : d;
+	      octave_write_double (os, d);
+	      os << "\n";
+	    }
+	}
+      else
+	{
+	  if (! infnan_warned && (xisnan (d) || xisinf (d)))
+	    {
+	      warning ("save: Inf or NaN values may not be reloadable");
+	      infnan_warned = true;
+	    }
+
+	  octave_write_double (os, d);
+	  os << "\n";
+	}
+    }
+  else if (val.is_real_matrix ())
+    {
+      ascii_save_type (os, "matrix", mark_as_global);
+
+      os << "# rows: " << val.rows () << "\n"
+	 << "# columns: " << val.columns () << "\n";
+
+      Matrix tmp = val.matrix_value ();
+
+      if (strip_nan_and_inf)
+	tmp = strip_infnan (tmp);
+      else if (! infnan_warned && tmp.any_element_is_inf_or_nan ())
+	{
+	  warning ("save: Inf or NaN values may not be reloadable");
+	  infnan_warned = true;
+	}
+
+      os << tmp;
+    }
+  else if (val.is_cell ())
+    {
+      ascii_save_type (os, "cell", mark_as_global);
+
+      os << "# rows: " << val.rows () << "\n"
+	 << "# columns: " << val.columns () << "\n";
+
+      Cell tmp = val.cell_value ();
+      
+      for (int j = 0; j < tmp.cols (); j++)
+	{
+	  for (int i = 0; i < tmp.rows (); i++)
+	    {
+	      octave_value o_val = tmp.elem (i, j);
+
+	      // Recurse to print sub-value.
+	      bool b = save_ascii_data (os, o_val, CELL_ELT_TAG,
+					infnan_warned, strip_nan_and_inf,
+					mark_as_global, 0);
+
+	      if (! b)
+		return os;
+	    }
+
+	  os << "\n";
+	}
+    }
+  else if (val.is_complex_scalar ())
+    {
+      ascii_save_type (os, "complex scalar", mark_as_global);
+
+      Complex c = val.complex_value ();
+
+      if (strip_nan_and_inf)
+	{
+	  if (xisnan (c))
+	    {
+	      error ("only value to plot is NaN");
+	      success = false;
+	    }
+	  else
+	    {
+	      double re = real (c);
+	      double im = imag (c);
+
+	      re = xisinf (re) ? (re > 0 ? OCT_RBV : -OCT_RBV) : re;
+	      im = xisinf (im) ? (im > 0 ? OCT_RBV : -OCT_RBV) : im;
+
+	      c = Complex (re, im);
+
+	      octave_write_complex (os, c);
+	      os << "\n";
+	    }
+	}
+      else
+	{
+	  if (! infnan_warned && (xisnan (c) || xisinf (c)))
+	    {
+	      warning ("save: Inf or NaN values may not be reloadable");
+	      infnan_warned = true;
+	    }
+
+	  octave_write_complex (os, c);
+	  os << "\n";
+	}
+    }
+  else if (val.is_complex_matrix ())
+    {
+      ascii_save_type (os, "complex matrix", mark_as_global);
+
+      os << "# rows: " << val.rows () << "\n"
+	 << "# columns: " << val.columns () << "\n";
+
+      ComplexMatrix tmp = val.complex_matrix_value ();
+
+      if (strip_nan_and_inf)
+	tmp = strip_infnan (tmp);
+      else if (! infnan_warned && tmp.any_element_is_inf_or_nan ())
+	{
+	  warning ("save: Inf or NaN values may not be reloadable");
+	  infnan_warned = true;
+	}
+
+      os << tmp;
+    }
+  else
+    gripe_wrong_type_arg ("save", val, false);
+
+  os.precision (old_precision);
+
+  return (os && success);
+}
+
+bool
+save_ascii_data_for_plotting (std::ostream& os, const octave_value& t,
+			      const std::string& name)
+{
+  bool infnan_warned = true;
+
+  return save_ascii_data (os, t, name, infnan_warned, true, false, 0);
+}
+
+// Maybe this should be a static function in tree-plot.cc?
+
+// If TC is matrix, save it on stream OS in a format useful for
+// making a 3-dimensional plot with gnuplot.  If PARAMETRIC is
+// TRUE, assume a parametric 3-dimensional plot will be generated.
+
+bool
+save_three_d (std::ostream& os, const octave_value& tc, bool parametric)
+{
+  bool fail = false;
+
+  int nr = tc.rows ();
+  int nc = tc.columns ();
+
+  if (tc.is_real_matrix ())
+    {
+      os << "# 3D data...\n"
+	 << "# type: matrix\n"
+	 << "# total rows: " << nr << "\n"
+	 << "# total columns: " << nc << "\n";
+
+      if (parametric)
+	{
+	  int extras = nc % 3;
+	  if (extras)
+	    warning ("ignoring last %d columns", extras);
+
+	  Matrix tmp = tc.matrix_value ();
+	  tmp = strip_infnan (tmp);
+	  nr = tmp.rows ();
+
+	  for (int i = 0; i < nc-extras; i += 3)
+	    {
+	      os << tmp.extract (0, i, nr-1, i+2);
+	      if (i+3 < nc-extras)
+		os << "\n";
+	    }
+	}
+      else
+	{
+	  Matrix tmp = tc.matrix_value ();
+	  tmp = strip_infnan (tmp);
+	  nr = tmp.rows ();
+
+	  for (int i = 0; i < nc; i++)
+	    {
+	      os << tmp.extract (0, i, nr-1, i);
+	      if (i+1 < nc)
+		os << "\n";
+	    }
+	}
+    }
+  else
+    {
+      ::error ("for now, I can only save real matrices in 3D format");
+      fail = true;
+    }
+
+  return (os && ! fail);
+}
+
+static int
+save_precision (void)
+{
+  double val;
+  if (builtin_real_scalar_variable ("save_precision", val)
+      && ! xisnan (val))
+    {
+      int ival = NINT (val);
+      if (ival >= 0 && ival == val)
+	{
+	  Vsave_precision = ival;
+	  return 0;
+	}
+    }
+  gripe_invalid_value_specified ("save_precision");
+  return -1;
+}
+
+void
+symbols_of_ls_oct_ascii (void)
+{
+  DEFVAR (save_precision, 15.0, save_precision,
+    "-*- texinfo -*-\n\
+@defvr {Built-in Variable} save_precision\n\
+This variable specifies the number of digits to keep when saving data in\n\
+text format.  The default value is 17.\n\
+@end defvr");
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-oct-ascii.h b/src/ls-oct-ascii.h
new file mode 100644
--- /dev/null
+++ b/src/ls-oct-ascii.h
@@ -0,0 +1,50 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_oct_ascii_h)
+#define octave_ls_oct_ascii_h 1
+
+extern std::string
+extract_keyword (std::istream& is, const char *keyword);
+
+extern std::string
+read_ascii_data (std::istream& is, const std::string& filename, bool& global,
+		 octave_value& tc, int count);
+
+extern bool
+save_ascii_data (std::ostream& os, const octave_value& val_arg,
+		 const std::string& name, bool& infnan_warned,
+		 bool strip_nan_and_inf, bool mark_as_global,
+		 int precision);
+
+extern bool
+save_ascii_data_for_plotting (std::ostream& os, const octave_value& t,
+			      const std::string& name);
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-oct-binary.cc b/src/ls-oct-binary.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-oct-binary.cc
@@ -0,0 +1,488 @@
+/*
+
+Copyright (C) 1996, 1997 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <cfloat>
+#include <cstring>
+#include <cctype>
+
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#include <memory>
+#include <string>
+
+#include "byte-swap.h"
+#include "data-conv.h"
+#include "file-ops.h"
+#include "glob-match.h"
+#include "lo-mappers.h"
+#include "lo-sstream.h"
+#include "mach-info.h"
+#include "oct-env.h"
+#include "oct-time.h"
+#include "quit.h"
+#include "str-vec.h"
+
+#include "Cell.h"
+#include "defun.h"
+#include "error.h"
+#include "gripes.h"
+#include "load-save.h"
+#include "oct-obj.h"
+#include "oct-map.h"
+#include "ov-cell.h"
+#include "pager.h"
+#include "pt-exp.h"
+#include "symtab.h"
+#include "sysdep.h"
+#include "unwind-prot.h"
+#include "utils.h"
+#include "variables.h"
+#include "version.h"
+#include "dMatrix.h"
+
+#include "ls-utils.h"
+#include "ls-oct-binary.h"
+
+// Extract one value (scalar, matrix, string, etc.) from stream IS and
+// place it in TC, returning the name of the variable.  If the value
+// is tagged as global in the file, return TRUE in GLOBAL.  If SWAP
+// is TRUE, swap bytes after reading.
+//
+// The data is expected to be in the following format:
+//
+// Header (one per file):
+// =====================
+//
+//   object               type            bytes
+//   ------               ----            -----
+//   magic number         string             10
+//
+//   float format         integer             1  
+//
+//
+// Data (one set for each item):
+// ============================
+//
+//   object               type            bytes
+//   ------               ----            -----
+//   name_length          integer             4
+//
+//   name                 string    name_length
+//
+//   doc_length           integer             4
+//
+//   doc                  string     doc_length
+//
+//   global flag          integer             1
+//
+//   data type            integer             1
+//
+//   data (one of):
+//
+//     scalar:
+//       data             real                8
+//
+//     complex scalar:
+//       data             complex            16
+//
+//     matrix:
+//       rows             integer             4
+//       columns          integer             4
+//       data             real            r*c*8
+//
+//     complex matrix:
+//       rows             integer             4
+//       columns          integer             4
+//       data             complex        r*c*16
+//
+//     string:
+//       length           int                 4
+//       data             string         length
+//
+//     range:
+//       base             real                8
+//       limit            real                8
+//       increment        real                8
+//
+//     string array
+//       elements         int                 4
+//
+//       for each element:
+//         length         int                 4
+//         data           string         length
+//
+// FILENAME is used for error messages.
+
+std::string
+read_binary_data (std::istream& is, bool swap,
+		  oct_mach_info::float_format fmt,
+		  const std::string& filename, bool& global,
+		  octave_value& tc, std::string& doc)
+{
+  std::string retval;
+
+  char tmp = 0;
+
+  FOUR_BYTE_INT name_len = 0;
+  FOUR_BYTE_INT doc_len = 0;
+
+  doc.resize (0);
+
+  // We expect to fail here, at the beginning of a record, so not
+  // being able to read another name should not result in an error.
+
+  is.read (X_CAST (char *, &name_len), 4);
+  if (! is)
+    return retval;
+  if (swap)
+    swap_4_bytes (X_CAST (char *, &name_len));
+
+  {
+    OCTAVE_LOCAL_BUFFER (char, name, name_len+1);
+    name[name_len] = '\0';
+    if (! is.read (X_CAST (char *, name), name_len))
+      goto data_read_error;
+    retval = name;
+  }
+
+  is.read (X_CAST (char *, &doc_len), 4);
+  if (! is)
+    goto data_read_error;
+  if (swap)
+    swap_4_bytes (X_CAST (char *, &doc_len));
+
+  {
+    OCTAVE_LOCAL_BUFFER (char, tdoc, doc_len+1);
+    tdoc[doc_len] = '\0';
+    if (! is.read (X_CAST (char *, tdoc), doc_len))
+      goto data_read_error;
+    doc = tdoc;
+  }
+
+  if (! is.read (X_CAST (char *, &tmp), 1))
+    goto data_read_error;
+  global = tmp ? 1 : 0;
+
+  tmp = 0;
+  if (! is.read (X_CAST (char *, &tmp), 1))
+    goto data_read_error;
+
+  switch (tmp)
+    {
+    case 1:
+      {
+	if (! is.read (X_CAST (char *, &tmp), 1))
+	  goto data_read_error;
+	double dtmp;
+	read_doubles (is, &dtmp, X_CAST (save_type, tmp), 1, swap, fmt);
+	if (error_state || ! is)
+	  goto data_read_error;
+	tc = dtmp;
+      }
+      break;
+
+    case 2:
+      {
+	FOUR_BYTE_INT nr, nc;
+	if (! is.read (X_CAST (char *, &nr), 4))
+	  goto data_read_error;
+	if (swap)
+	  swap_4_bytes (X_CAST (char *, &nr));
+	if (! is.read (X_CAST (char *, &nc), 4))
+	  goto data_read_error;
+	if (swap)
+	  swap_4_bytes (X_CAST (char *, &nc));
+	if (! is.read (X_CAST (char *, &tmp), 1))
+	  goto data_read_error;
+	Matrix m (nr, nc);
+	double *re = m.fortran_vec ();
+	int len = nr * nc;
+	read_doubles (is, re, X_CAST (save_type, tmp), len, swap, fmt);
+	if (error_state || ! is)
+	  goto data_read_error;
+	tc = m;
+      }
+      break;
+
+    case 3:
+      {
+	if (! is.read (X_CAST (char *, &tmp), 1))
+	  goto data_read_error;
+	Complex ctmp;
+	read_doubles (is, X_CAST (double *, &ctmp),
+		      X_CAST (save_type, tmp), 2, swap, fmt);
+	if (error_state || ! is)
+	  goto data_read_error;
+	tc = ctmp;
+      }
+      break;
+
+    case 4:
+      {
+	FOUR_BYTE_INT nr, nc;
+	if (! is.read (X_CAST (char *, &nr), 4))
+	  goto data_read_error;
+	if (swap)
+	  swap_4_bytes (X_CAST (char *, &nr));
+	if (! is.read (X_CAST (char *, &nc), 4))
+	  goto data_read_error;
+	if (swap)
+	  swap_4_bytes (X_CAST (char *, &nc));
+	if (! is.read (X_CAST (char *, &tmp), 1))
+	  goto data_read_error;
+	ComplexMatrix m (nr, nc);
+	Complex *im = m.fortran_vec ();
+	int len = nr * nc;
+	read_doubles (is, X_CAST (double *, im),
+		      X_CAST (save_type, tmp), 2*len, swap, fmt);
+	if (error_state || ! is)
+	  goto data_read_error;
+	tc = m;
+      }
+      break;
+
+    case 5:
+      {
+	FOUR_BYTE_INT len;
+	if (! is.read (X_CAST (char *, &len), 4))
+	  goto data_read_error;
+	if (swap)
+	  swap_4_bytes (X_CAST (char *, &len));
+	OCTAVE_LOCAL_BUFFER (char, s, len+1);
+	if (! is.read (X_CAST (char *, s), len))
+	  goto data_read_error;
+	s[len] = '\0';
+	tc = s;
+      }
+      break;
+
+    case 6:
+      {
+	if (! is.read (X_CAST (char *, &tmp), 1))
+	  goto data_read_error;
+	double bas, lim, inc;
+	if (! is.read (X_CAST (char *, &bas), 8))
+	  goto data_read_error;
+	if (swap)
+	  swap_8_bytes (X_CAST (char *, &bas));
+	if (! is.read (X_CAST (char *, &lim), 8))
+	  goto data_read_error;
+	if (swap)
+	  swap_8_bytes (X_CAST (char *, &lim));
+	if (! is.read (X_CAST (char *, &inc), 8))
+	  goto data_read_error;
+	if (swap)
+	  swap_8_bytes (X_CAST (char *, &inc));
+	Range r (bas, lim, inc);
+	tc = r;
+      }
+      break;
+
+    case 7:
+      {
+	FOUR_BYTE_INT elements;
+	if (! is.read (X_CAST (char *, &elements), 4))
+	  goto data_read_error;
+	if (swap)
+	  swap_4_bytes (X_CAST (char *, &elements));
+	charMatrix chm (elements, 0);
+	int max_len = 0;
+	for (int i = 0; i < elements; i++)
+	  {
+	    FOUR_BYTE_INT len;
+	    if (! is.read (X_CAST (char *, &len), 4))
+	      goto data_read_error;
+	    if (swap)
+	      swap_4_bytes (X_CAST (char *, &len));
+	    OCTAVE_LOCAL_BUFFER (char, btmp, len+1);
+	    if (! is.read (X_CAST (char *, btmp), len))
+	      goto data_read_error;
+	    if (len > max_len)
+	      {
+		max_len = len;
+		chm.resize (elements, max_len, 0);
+	      }
+	    btmp [len] = '\0';
+	    chm.insert (btmp, i, 0);
+	  }
+	tc = octave_value (chm, true);
+      }
+      break;
+
+    default:
+    data_read_error:
+      error ("load: trouble reading binary file `%s'", filename.c_str ());
+      break;
+    }
+
+  return retval;
+}
+
+// Save the data from TC along with the corresponding NAME, help
+// string DOC, and global flag MARK_AS_GLOBAL on stream OS in the
+// binary format described above for read_binary_data.
+
+bool
+save_binary_data (std::ostream& os, const octave_value& tc,
+		  const std::string& name, const std::string& doc,
+		  bool mark_as_global, bool save_as_floats) 
+{
+  FOUR_BYTE_INT name_len = name.length ();
+
+  os.write (X_CAST (char *, &name_len), 4);
+  os << name;
+
+  FOUR_BYTE_INT doc_len = doc.length ();
+
+  os.write (X_CAST (char *, &doc_len), 4);
+  os << doc;
+
+  char tmp;
+
+  tmp = mark_as_global;
+  os.write (X_CAST (char *, &tmp), 1);
+
+  if (tc.is_string ())
+    {
+      tmp = 7;
+      os.write (X_CAST (char *, &tmp), 1);
+      FOUR_BYTE_INT nr = tc.rows ();
+      os.write (X_CAST (char *, &nr), 4);
+      charMatrix chm = tc.char_matrix_value ();
+      for (int i = 0; i < nr; i++)
+	{
+	  FOUR_BYTE_INT len = chm.cols ();
+	  os.write (X_CAST (char *, &len), 4);
+	  std::string tstr = chm.row_as_string (i);
+	  const char *btmp = tstr.data ();
+	  os.write (X_CAST (char *, btmp), len);
+	}
+    }
+  else if (tc.is_range ())
+    {
+      tmp = 6;
+      os.write (X_CAST (char *, &tmp), 1);
+      tmp = (char) LS_DOUBLE;
+      os.write (X_CAST (char *, &tmp), 1);
+      Range r = tc.range_value ();
+      double bas = r.base ();
+      double lim = r.limit ();
+      double inc = r.inc ();
+      os.write (X_CAST (char *, &bas), 8);
+      os.write (X_CAST (char *, &lim), 8);
+      os.write (X_CAST (char *, &inc), 8);
+    }
+  else if (tc.is_real_scalar ())
+    {
+      tmp = 1;
+      os.write (X_CAST (char *, &tmp), 1);
+      tmp = (char) LS_DOUBLE;
+      os.write (X_CAST (char *, &tmp), 1);
+      double dtmp = tc.double_value ();
+      os.write (X_CAST (char *, &dtmp), 8);
+    }
+  else if (tc.is_real_matrix ())
+    {
+      tmp = 2;
+      os.write (X_CAST (char *, &tmp), 1);
+      Matrix m = tc.matrix_value ();
+      FOUR_BYTE_INT nr = m.rows ();
+      FOUR_BYTE_INT nc = m.columns ();
+      os.write (X_CAST (char *, &nr), 4);
+      os.write (X_CAST (char *, &nc), 4);
+      int len = nr * nc;
+      save_type st = LS_DOUBLE;
+      if (save_as_floats)
+	{
+	  if (m.too_large_for_float ())
+	    {
+	      warning ("save: some values too large to save as floats --");
+	      warning ("save: saving as doubles instead");
+	    }
+	  else
+	    st = LS_FLOAT;
+	}
+      else if (len > 8192) // XXX FIXME XXX -- make this configurable.
+	{
+	  double max_val, min_val;
+	  if (m.all_integers (max_val, min_val))
+	    st = get_save_type (max_val, min_val);
+	}
+      const double *mtmp = m.data ();
+      write_doubles (os, mtmp, st, len);
+    }
+  else if (tc.is_complex_scalar ())
+    {
+      tmp = 3;
+      os.write (X_CAST (char *, &tmp), 1);
+      tmp = (char) LS_DOUBLE;
+      os.write (X_CAST (char *, &tmp), 1);
+      Complex ctmp = tc.complex_value ();
+      os.write (X_CAST (char *, &ctmp), 16);
+    }
+  else if (tc.is_complex_matrix ())
+    {
+      tmp = 4;
+      os.write (X_CAST (char *, &tmp), 1);
+      ComplexMatrix m = tc.complex_matrix_value ();
+      FOUR_BYTE_INT nr = m.rows ();
+      FOUR_BYTE_INT nc = m.columns ();
+      os.write (X_CAST (char *, &nr), 4);
+      os.write (X_CAST (char *, &nc), 4);
+      int len = nr * nc;
+      save_type st = LS_DOUBLE;
+      if (save_as_floats)
+	{
+	  if (m.too_large_for_float ())
+	    {
+	      warning ("save: some values too large to save as floats --");
+	      warning ("save: saving as doubles instead");
+	    }
+	  else
+	    st = LS_FLOAT;
+	}
+      else if (len > 4096) // XXX FIXME XXX -- make this configurable.
+	{
+	  double max_val, min_val;
+	  if (m.all_integers (max_val, min_val))
+	    st = get_save_type (max_val, min_val);
+	}
+      const Complex *mtmp = m.data ();
+      write_doubles (os, X_CAST (const double *, mtmp), st, 2*len);
+    }
+  else
+    gripe_wrong_type_arg ("save", tc, false);
+
+  return os;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-oct-binary.h b/src/ls-oct-binary.h
new file mode 100644
--- /dev/null
+++ b/src/ls-oct-binary.h
@@ -0,0 +1,44 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_oct_binary_h)
+#define octave_ls_oct_binary_h 1
+
+extern bool
+save_binary_data (std::ostream& os, const octave_value& tc,
+		  const std::string& name, const std::string& doc,
+		  bool mark_as_global, bool save_as_floats);
+
+extern std::string
+read_binary_data (std::istream& is, bool swap,
+		  oct_mach_info::float_format fmt,
+		  const std::string& filename, bool& global,
+		  octave_value& tc, std::string& doc);
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
+
diff --git a/src/ls-utils.cc b/src/ls-utils.cc
new file mode 100644
--- /dev/null
+++ b/src/ls-utils.cc
@@ -0,0 +1,66 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include "data-conv.h"
+
+#include "ls-utils.h"
+
+// MAX_VAL and MIN_VAL are assumed to have integral values even though
+// they are stored in doubles.
+
+save_type
+get_save_type (double max_val, double min_val)
+{
+  save_type st = LS_DOUBLE;
+
+  // Matlab doesn't seem to load the UINT32 type correctly, so let's
+  // avoid it (and the other unsigned types, even though they may not
+  // have the same problem.  And apparently, there are problems with
+  // other smaller types as well.  If we avoid them all, then maybe we
+  // will avoid problems.  Unfortunately, we won't be able to save
+  // space...
+
+  //  if (max_val < 256 && min_val > -1)
+  //    st = LS_U_CHAR;
+  //  else if (max_val < 65536 && min_val > -1)
+  //    st = LS_U_SHORT;
+  //  else if (max_val < 4294967295UL && min_val > -1)
+  //    st = LS_U_INT;
+  //  else if (max_val < 128 && min_val >= -128)
+  //    st = LS_CHAR;
+  //  else if (max_val < 32768 && min_val >= -32768)
+  //    st = LS_SHORT;
+  //  else if (max_val <= 2147483647L && min_val >= -2147483647L)
+  //    st = LS_INT;
+
+  return st;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/
diff --git a/src/ls-utils.h b/src/ls-utils.h
new file mode 100644
--- /dev/null
+++ b/src/ls-utils.h
@@ -0,0 +1,35 @@
+/*
+
+Copyright (C) 2003 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+*/
+
+#if !defined (octave_ls_utils_h)
+#define octave_ls_utils 1
+
+extern save_type
+get_save_type (double max_val, double min_val);
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; End: ***
+*/