--- a/liboctave/Array.cc
+++ b/liboctave/Array.cc
@@ -51,11 +51,19 @@
   : rep (a.rep), dimensions (dv), 
     slice_data (a.slice_data), slice_len (a.slice_len)
 {
+  if (dv.numel () != a.numel ())
+    {
+      std::string dimensions_str = dimensions.str ();
+      std::string new_dims_str = dv.str ();
+
+      (*current_liboctave_error_handler)
+        ("reshape: can't reshape %s array to %s array",
+         dimensions_str.c_str (), new_dims_str.c_str ());
+    }
+
+  // This goes here because if an exception is thrown by the above,
+  // destructor will be never called.
   rep->count++;
-
-  if (a.numel () < dv.numel ())
-    (*current_liboctave_error_handler)
-      ("Array::Array (const Array&, const dim_vector&): dimension mismatch");
 }
 
 template <class T>
@@ -414,30 +422,51 @@
   return foo;
 }
 
+
+template <class T>
+Array<T>
+Array<T>::column (octave_idx_type k) const
+{
+  octave_idx_type r = dimensions(0);
+#ifdef BOUNDS_CHECKING
+  if (k < 0 || k * r >= numel ())
+    range_error ("column", k);
+#endif
+
+  return Array<T> (*this, dim_vector (r, 1), k*r, k*r + r);
+}
+
 template <class T>
 Array<T>
-Array<T>::reshape (const dim_vector& new_dims) const
+Array<T>::page (octave_idx_type k) const
+{
+  octave_idx_type r = dimensions(0), c = dimensions (1), p = r*c;
+#ifdef BOUNDS_CHECKING
+  if (k < 0 || k * p >= numel ())
+    range_error ("page", k);
+#endif
+
+  return Array<T> (*this, dim_vector (r, c), k*p, k*p + p);
+}
+
+template <class T>
+Array<T>
+Array<T>::linearize (void) const
 {
-  Array<T> retval;
-
-  if (dimensions != new_dims)
-    {
-      if (dimensions.numel () == new_dims.numel ())
-	retval = Array<T> (*this, new_dims);
-      else
-	{
-	  std::string dimensions_str = dimensions.str ();
-	  std::string new_dims_str = new_dims.str ();
-
-	  (*current_liboctave_error_handler)
-	    ("reshape: can't reshape %s array to %s array",
-	     dimensions_str.c_str (), new_dims_str.c_str ());
-	}
-    }
-  else
-    retval = *this;
-
-  return retval;
+  octave_idx_type n = numel ();
+  return Array<T> (*this, dim_vector (n, 1), 0, n);
+}
+
+template <class T>
+Array<T>
+Array<T>::linear_slice (octave_idx_type lo, octave_idx_type up) const
+{
+#ifdef BOUNDS_CHECKING
+  if (lo < 0 || up > numel ())
+    range_error ("linear_slice", lo, up);
+#endif
+  if (up < lo) up = lo;
+  return Array<T> (*this, dim_vector (up - lo, 1), lo, up);
 }
 
 // Helper class for multi-d dimension permuting (generalized transpose).

--- a/liboctave/Array.h
+++ b/liboctave/Array.h
@@ -151,7 +151,7 @@
     {
       rep->count++;
       slice_data = a.slice_data + l;
-      slice_len = std::min (u, a.slice_len) - l;
+      slice_len = u - l;
     }
 
 private:
@@ -239,6 +239,7 @@
       fill (val);
     }
 
+  // Reshape constructor.
   Array (const Array<T>& a, const dim_vector& dv);
 
   ~Array (void)
@@ -455,7 +456,23 @@
   T operator () (const Array<octave_idx_type>& ra_idx) const { return elem (ra_idx); }
 #endif
 
-  Array<T> reshape (const dim_vector& new_dims) const;
+  // Fast extractors. All of these produce shallow copies.
+  // Warning: none of these do check bounds, unless BOUNDS_CHECKING is on!
+
+  // Extract column: A(:,k+1).
+  Array<T> column (octave_idx_type k) const;
+  // Extract page: A(:,:,k+1).
+  Array<T> page (octave_idx_type k) const;
+
+  // Give this array as a column vector: A(:).
+  Array<T> linearize (void) const;
+
+  // Extract a slice from this array as a column vector: A(:)(lo+1:up).
+  // Must be 0 <= lo && up <= numel. May be up < lo.
+  Array<T> linear_slice (octave_idx_type lo, octave_idx_type up) const;
+
+  Array<T> reshape (const dim_vector& new_dims) const
+    { return Array<T> (*this, new_dims); }
 
   Array<T> permute (const Array<octave_idx_type>& vec, bool inv = false) const;
   Array<T> ipermute (const Array<octave_idx_type>& vec) const

--- a/liboctave/ArrayN.h
+++ b/liboctave/ArrayN.h
@@ -84,6 +84,16 @@
       return *this;
     }
 
+  ArrayN<T> column (octave_idx_type k) const
+    { return Array<T>::column (k); }
+  ArrayN<T> page (octave_idx_type k) const
+    { return Array<T>::page (k); }
+
+  ArrayN<T> linearize (void) const
+    { return Array<T>::linearize (); }
+  ArrayN<T> linear_slice (octave_idx_type lo, octave_idx_type up) const
+    { return Array<T>::linear_slice (lo, up); }
+
   ArrayN<T> reshape (const dim_vector& new_dims) const
     { return Array<T>::reshape (new_dims); }
 

--- a/liboctave/ChangeLog
+++ b/liboctave/ChangeLog
@@ -1,3 +1,10 @@
+2009-10-16  Jaroslav Hajek  <highegg@gmail.com>
+
+	* Array.cc (Array<T>::column, Array<T>::page, Array<T>::linearize,
+	Array<T>::linear_slice): New methods.
+	* Array.h: Declare them
+	* ArrayN.h: Forward them.
+
 2009-10-14  Jaroslav Hajek  <highegg@gmail.com>
 
 	* oct-sort.cc (octave_sort<T>::nth_element): New overloaded method.

--- a/src/ChangeLog
+++ b/src/ChangeLog
@@ -1,3 +1,7 @@
+2009-10-16  Jaroslav Hajek  <highegg@gmail.com>
+
+	* DLD-FUNCTIONS/cellfun.cc (Fnum2cell): Use Array<T>::column here.
+
 2009-10-15  Jaroslav Hajek  <highegg@gmail.com>
 
 	* DLD-FUNCTIONS/lookup.cc (Flookup): Update docstring.

--- a/src/DLD-FUNCTIONS/cellfun.cc
+++ b/src/DLD-FUNCTIONS/cellfun.cc
@@ -993,8 +993,7 @@
       Cell retval (celldv);
       for (octave_idx_type i = 0; i < nelc; i++)
         {
-          NDA tmp (parray.index (idx_vector::colon, idx_vector (i)));
-          retval.xelem (i) = tmp.reshape (arraydv);
+          retval.xelem (i) = NDA (parray.column (i).reshape (arraydv));
         }
 
       return retval;