Mercurial > hg > octave-nkf
view src/data.cc @ 2475:1d7925d6bede
[project @ 1996-11-07 04:36:00 by jwe]
author | jwe |
---|---|
date | Thu, 07 Nov 1996 04:46:54 +0000 |
parents | 449f35baba49 |
children | e93b4c32457b |
line wrap: on
line source
/* Copyright (C) 1996 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. */ /* The function builtin_pwd adapted from a similar function from GNU Bash, the Bourne Again SHell, copyright (C) 1987, 1989, 1991 Free Software Foundation, Inc. */ #ifdef HAVE_CONFIG_H #include <config.h> #endif #include <cfloat> #include <cmath> #include <string> #include "lo-ieee.h" #include "str-vec.h" #include "defun.h" #include "error.h" #include "gripes.h" #include "help.h" #include "oct-map.h" #include "ov.h" #include "variables.h" #include "oct-obj.h" #include "utils.h" #ifndef MIN #define MIN(a,b) ((a) < (b) ? (a) : (b)) #endif #ifndef ABS #define ABS(x) (((x) < 0) ? (-x) : (x)) #endif // Should expressions like ones (-1, 5) result in an empty matrix or // an error? A positive value means yes. A negative value means // yes, but print a warning message. Zero means it should be // considered an error. static int Vtreat_neg_dim_as_zero; DEFUN (all, args, , "all (X): are all elements of X nonzero?") { octave_value_list retval; int nargin = args.length (); if (nargin == 1 && args(0).is_defined ()) retval = args(0).all (); else print_usage ("all"); return retval; } DEFUN (any, args, , "any (X): are any elements of X nonzero?") { octave_value_list retval; int nargin = args.length (); if (nargin == 1 && args(0).is_defined ()) retval = args(0).any (); else print_usage ("any"); return retval; } // These mapping functions may also be useful in other places, eh? typedef double (*d_dd_fcn) (double, double); static Matrix map (d_dd_fcn f, double x, const Matrix& y) { int nr = y.rows (); int nc = y.columns (); Matrix retval (nr, nc); for (int j = 0; j < nc; j++) for (int i = 0; i < nr; i++) retval (i, j) = f (x, y (i, j)); return retval; } static Matrix map (d_dd_fcn f, const Matrix& x, double y) { int nr = x.rows (); int nc = x.columns (); Matrix retval (nr, nc); for (int j = 0; j < nc; j++) for (int i = 0; i < nr; i++) retval (i, j) = f (x (i, j), y); return retval; } static Matrix map (d_dd_fcn f, const Matrix& x, const Matrix& y) { int x_nr = x.rows (); int x_nc = x.columns (); int y_nr = y.rows (); int y_nc = y.columns (); assert (x_nr == y_nr && x_nc == y_nc); Matrix retval (x_nr, x_nc); for (int j = 0; j < x_nc; j++) for (int i = 0; i < x_nr; i++) retval (i, j) = f (x (i, j), y (i, j)); return retval; } DEFUN (atan2, args, , "atan2 (Y, X): atan (Y / X) in range -pi to pi") { octave_value_list retval; int nargin = args.length (); if (nargin == 2 && args(0).is_defined () && args(1).is_defined ()) { octave_value arg_y = args(0); octave_value arg_x = args(1); int y_nr = arg_y.rows (); int y_nc = arg_y.columns (); int x_nr = arg_x.rows (); int x_nc = arg_x.columns (); int arg_y_empty = empty_arg ("atan2", y_nr, y_nc); int arg_x_empty = empty_arg ("atan2", x_nr, x_nc); if (arg_y_empty > 0 && arg_x_empty > 0) return Matrix (); else if (arg_y_empty || arg_x_empty) return retval; int y_is_scalar = (y_nr == 1 && y_nc == 1); int x_is_scalar = (x_nr == 1 && x_nc == 1); if (y_is_scalar && x_is_scalar) { double y = arg_y.double_value (); if (! error_state) { double x = arg_x.double_value (); if (! error_state) retval = atan2 (y, x); } } else if (y_is_scalar) { double y = arg_y.double_value (); if (! error_state) { Matrix x = arg_x.matrix_value (); if (! error_state) retval = map (atan2, y, x); } } else if (x_is_scalar) { Matrix y = arg_y.matrix_value (); if (! error_state) { double x = arg_x.double_value (); if (! error_state) retval = map (atan2, y, x); } } else if (y_nr == x_nr && y_nc == x_nc) { Matrix y = arg_y.matrix_value (); if (! error_state) { Matrix x = arg_x.matrix_value (); if (! error_state) retval = map (atan2, y, x); } } else error ("atan2: nonconformant matrices"); } else print_usage ("atan2"); return retval; } DEFUN (cumprod, args, , "cumprod (X): cumulative products") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { octave_value arg = args(0); if (arg.is_real_type ()) { Matrix tmp = arg.matrix_value (); if (! error_state) retval(0) = tmp.cumprod (); } else if (arg.is_complex_type ()) { ComplexMatrix tmp = arg.complex_matrix_value (); if (! error_state) retval(0) = tmp.cumprod (); } else { gripe_wrong_type_arg ("cumprod", arg); return retval; } } else print_usage ("cumprod"); return retval; } DEFUN (cumsum, args, , "cumsum (X): cumulative sums") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { octave_value arg = args(0); if (arg.is_real_type ()) { Matrix tmp = arg.matrix_value (); if (! error_state) retval(0) = tmp.cumsum (); } else if (arg.is_complex_type ()) { ComplexMatrix tmp = arg.complex_matrix_value (); if (! error_state) retval(0) = tmp.cumsum (); } else { gripe_wrong_type_arg ("cumsum", arg); return retval; } } else print_usage ("cumsum"); return retval; } static octave_value make_diag (const Matrix& v, int k) { int nr = v.rows (); int nc = v.columns (); assert (nc == 1 || nr == 1); octave_value retval; int roff = 0; int coff = 0; if (k > 0) { roff = 0; coff = k; } else if (k < 0) { roff = -k; coff = 0; } if (nr == 1) { int n = nc + ABS (k); Matrix m (n, n, 0.0); for (int i = 0; i < nc; i++) m (i+roff, i+coff) = v (0, i); retval = octave_value (m); } else { int n = nr + ABS (k); Matrix m (n, n, 0.0); for (int i = 0; i < nr; i++) m (i+roff, i+coff) = v (i, 0); retval = octave_value (m); } return retval; } static octave_value make_diag (const ComplexMatrix& v, int k) { int nr = v.rows (); int nc = v.columns (); assert (nc == 1 || nr == 1); octave_value retval; int roff = 0; int coff = 0; if (k > 0) { roff = 0; coff = k; } else if (k < 0) { roff = -k; coff = 0; } if (nr == 1) { int n = nc + ABS (k); ComplexMatrix m (n, n, 0.0); for (int i = 0; i < nc; i++) m (i+roff, i+coff) = v (0, i); retval = octave_value (m); } else { int n = nr + ABS (k); ComplexMatrix m (n, n, 0.0); for (int i = 0; i < nr; i++) m (i+roff, i+coff) = v (i, 0); retval = octave_value (m); } return retval; } static octave_value make_diag (const octave_value& arg) { octave_value retval; if (arg.is_real_type ()) { Matrix m = arg.matrix_value (); if (! error_state) { int nr = m.rows (); int nc = m.columns (); if (nr == 0 || nc == 0) retval = Matrix (); else if (nr == 1 || nc == 1) retval = make_diag (m, 0); else { ColumnVector v = m.diag (); if (v.capacity () > 0) retval = v; } } else gripe_wrong_type_arg ("diag", arg); } else if (arg.is_complex_type ()) { ComplexMatrix cm = arg.complex_matrix_value (); if (! error_state) { int nr = cm.rows (); int nc = cm.columns (); if (nr == 0 || nc == 0) retval = Matrix (); else if (nr == 1 || nc == 1) retval = make_diag (cm, 0); else { ComplexColumnVector v = cm.diag (); if (v.capacity () > 0) retval = v; } } else gripe_wrong_type_arg ("diag", arg); } else gripe_wrong_type_arg ("diag", arg); return retval; } static octave_value make_diag (const octave_value& a, const octave_value& b) { octave_value retval; double tmp = b.double_value (); if (error_state) { error ("diag: invalid second argument"); return retval; } int k = NINT (tmp); int n = ABS (k) + 1; if (a.is_real_type ()) { if (a.is_scalar_type ()) { double d = a.double_value (); if (k == 0) retval = d; else if (k > 0) { Matrix m (n, n, 0.0); m (0, k) = d; retval = m; } else if (k < 0) { Matrix m (n, n, 0.0); m (-k, 0) = d; retval = m; } } else if (a.is_matrix_type ()) { Matrix m = a.matrix_value (); int nr = m.rows (); int nc = m.columns (); if (nr == 0 || nc == 0) retval = Matrix (); else if (nr == 1 || nc == 1) retval = make_diag (m, k); else { ColumnVector d = m.diag (k); retval = d; } } else gripe_wrong_type_arg ("diag", a); } else if (a.is_complex_type ()) { if (a.is_scalar_type ()) { Complex c = a.complex_value (); if (k == 0) retval = c; else if (k > 0) { ComplexMatrix m (n, n, 0.0); m (0, k) = c; retval = m; } else if (k < 0) { ComplexMatrix m (n, n, 0.0); m (-k, 0) = c; retval = m; } } else if (a.is_matrix_type ()) { ComplexMatrix cm = a.complex_matrix_value (); int nr = cm.rows (); int nc = cm.columns (); if (nr == 0 || nc == 0) retval = Matrix (); else if (nr == 1 || nc == 1) retval = make_diag (cm, k); else { ComplexColumnVector d = cm.diag (k); retval = d; } } else gripe_wrong_type_arg ("diag", a); } else gripe_wrong_type_arg ("diag", a); return retval; } DEFUN (diag, args, , "diag (X [,k]): form/extract diagonals") { octave_value_list retval; int nargin = args.length (); if (nargin == 1 && args(0).is_defined ()) retval = make_diag (args(0)); else if (nargin == 2 && args(0).is_defined () && args(1).is_defined ()) retval = make_diag (args(0), args(1)); else print_usage ("diag"); return retval; } DEFUN (prod, args, , "prod (X): products") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { octave_value arg = args(0); if (arg.is_real_type ()) { Matrix tmp = arg.matrix_value (); if (! error_state) retval(0) = tmp.prod (); } else if (arg.is_complex_type ()) { ComplexMatrix tmp = arg.complex_matrix_value (); if (! error_state) retval(0) = tmp.prod (); } else { gripe_wrong_type_arg ("prod", arg); return retval; } } else print_usage ("prod"); return retval; } DEFUN (size, args, nargout, "[m, n] = size (x): return rows and columns of X\n\ \n\ d = size (x): return number of rows and columns of x as a row vector\n\ \n\ m = size (x, 1): return number of rows in x\n\ m = size (x, 2): return number of columns in x") { octave_value_list retval; int nargin = args.length (); if (nargin == 1 && nargout < 3) { int nr = args(0).rows (); int nc = args(0).columns (); if (nargout == 0 || nargout == 1) { Matrix m (1, 2); m (0, 0) = nr; m (0, 1) = nc; retval = m; } else if (nargout == 2) { retval(1) = (double) nc; retval(0) = (double) nr; } } else if (nargin == 2 && nargout < 2) { int nd = NINT (args(1).double_value ()); if (error_state) error ("size: expecting scalar as second argument"); else { if (nd == 1) retval(0) = (double) (args(0).rows ()); else if (nd == 2) retval(0) = (double) (args(0).columns ()); else error ("size: invalid second argument -- expecting 1 or 2"); } } else print_usage ("size"); return retval; } DEFUN (sum, args, , "sum (X): sum of elements") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { octave_value arg = args(0); if (arg.is_real_type ()) { Matrix tmp = arg.matrix_value (); if (! error_state) retval(0) = tmp.sum (); } else if (arg.is_complex_type ()) { ComplexMatrix tmp = arg.complex_matrix_value (); if (! error_state) retval(0) = tmp.sum (); } else { gripe_wrong_type_arg ("sum", arg); return retval; } } else print_usage ("sum"); return retval; } DEFUN (sumsq, args, , "sumsq (X): sum of squares of elements") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { octave_value arg = args(0); if (arg.is_real_type ()) { Matrix tmp = arg.matrix_value (); if (! error_state) retval(0) = tmp.sumsq (); } else if (arg.is_complex_type ()) { ComplexMatrix tmp = arg.complex_matrix_value (); if (! error_state) retval(0) = tmp.sumsq (); } else { gripe_wrong_type_arg ("sumsq", arg); return retval; } } else print_usage ("sumsq"); return retval; } DEFUN (is_struct, args, , "is_struct (x): return nonzero if x is a structure") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { octave_value arg = args(0); if (arg.is_map ()) retval = 1.0; else retval = 0.0; } else print_usage ("is_struct"); return retval; } DEFUN (struct_elements, args, , "struct_elements (S)\n\ \n\ Return a list of the names of the elements of the structure S.") { octave_value_list retval; int nargin = args.length (); if (nargin == 1) { if (args (0).is_map ()) { Octave_map m = args(0).map_value (); retval(0) = m.make_name_list (); } else gripe_wrong_type_arg ("struct_elements", args (0)); } else print_usage ("struct_elements"); return retval; } DEFUN (struct_contains, args, , "struct_contains (S, NAME)\n\ \n\ Return nonzero if S is a structure with element NAME.\n\ S must be a structure and NAME must be a string.") { octave_value_list retval; int nargin = args.length (); if (nargin == 2) { retval = 0.0; if (args(0).is_map () && args(1).is_string ()) { string s = args(1).string_value (); octave_value tmp = args(0).struct_elt_val (s, true); retval = (double) tmp.is_defined (); } else print_usage ("struct_contains"); } else print_usage ("struct_contains"); return retval; } static void check_dimensions (int& nr, int& nc, const char *warnfor) { if (nr < 0 || nc < 0) { if (Vtreat_neg_dim_as_zero) { nr = (nr < 0) ? 0 : nr; nc = (nc < 0) ? 0 : nc; if (Vtreat_neg_dim_as_zero < 0) warning ("%s: converting negative dimension to zero", warnfor); } else error ("%s: can't create a matrix with negative dimensions", warnfor); } } static void get_dimensions (const octave_value& a, const char *warn_for, int& nr, int& nc) { if (a.is_scalar_type ()) { double tmp = a.double_value (); nr = nc = NINT (tmp); } else { nr = a.rows (); nc = a.columns (); if ((nr == 1 && nc == 2) || (nr == 2 && nc == 1)) { ColumnVector v = a.vector_value (); if (error_state) return; nr = NINT (v (0)); nc = NINT (v (1)); } else warning ("%s (A): use %s (size (A)) instead", warn_for, warn_for); } check_dimensions (nr, nc, warn_for); // May set error_state. } static void get_dimensions (const octave_value& a, const octave_value& b, const char *warn_for, int& nr, int& nc) { nr = NINT (a.double_value ()); nc = NINT (b.double_value ()); if (error_state) error ("%s: expecting two scalar arguments", warn_for); else check_dimensions (nr, nc, warn_for); // May set error_state. } static octave_value fill_matrix (const octave_value& a, double val, const char *warn_for) { int nr, nc; get_dimensions (a, warn_for, nr, nc); if (error_state) return octave_value (); Matrix m (nr, nc, val); return m; } static octave_value fill_matrix (const octave_value& a, const octave_value& b, double val, const char *warn_for) { int nr, nc; get_dimensions (a, b, warn_for, nr, nc); // May set error_state. if (error_state) return octave_value (); Matrix m (nr, nc, val); return m; } DEFUN (ones, args, , "ones (N), ones (N, M), ones (X): create a matrix of all ones") { octave_value_list retval; int nargin = args.length (); switch (nargin) { case 0: retval = 1.0; break; case 1: retval = fill_matrix (args(0), 1.0, "ones"); break; case 2: retval = fill_matrix (args(0), args(1), 1.0, "ones"); break; default: print_usage ("ones"); break; } return retval; } DEFUN (zeros, args, , "zeros (N), zeros (N, M), zeros (X): create a matrix of all zeros") { octave_value_list retval; int nargin = args.length (); switch (nargin) { case 0: retval = 0.0; break; case 1: retval = fill_matrix (args(0), 0.0, "zeros"); break; case 2: retval = fill_matrix (args(0), args(1), 0.0, "zeros"); break; default: print_usage ("zeros"); break; } return retval; } static octave_value identity_matrix (const octave_value& a) { int nr, nc; get_dimensions (a, "eye", nr, nc); // May set error_state. if (error_state) return octave_value (); Matrix m (nr, nc, 0.0); if (nr > 0 && nc > 0) { int n = MIN (nr, nc); for (int i = 0; i < n; i++) m (i, i) = 1.0; } return m; } static octave_value identity_matrix (const octave_value& a, const octave_value& b) { int nr, nc; get_dimensions (a, b, "eye", nr, nc); // May set error_state. if (error_state) return octave_value (); Matrix m (nr, nc, 0.0); if (nr > 0 && nc > 0) { int n = MIN (nr, nc); for (int i = 0; i < n; i++) m (i, i) = 1.0; } return m; } DEFUN (eye, args, , "eye (N), eye (N, M), eye (X): create an identity matrix") { octave_value_list retval; int nargin = args.length (); switch (nargin) { case 0: retval = 1.0; break; case 1: retval = identity_matrix (args(0)); break; case 2: retval = identity_matrix (args(0), args(1)); break; default: print_usage ("eye"); break; } return retval; } DEFUN (linspace, args, , "usage: linspace (x1, x2, n)\n\ \n\ Return a vector of n equally spaced points between x1 and x2\n\ inclusive.\n\ \n\ If the final argument is omitted, n = 100 is assumed.\n\ \n\ All three arguments must be scalars.\n\ \n\ See also: logspace") { octave_value_list retval; int nargin = args.length (); int npoints = 100; if (nargin != 2 && nargin != 3) { print_usage ("linspace"); return retval; } if (nargin == 3) { double n = args(2).double_value (); if (! error_state) npoints = NINT (n); } if (! error_state) { if (npoints > 1) { octave_value arg_1 = args(0); octave_value arg_2 = args(1); if (arg_1.is_complex_type () || arg_2.is_complex_type ()) { Complex x1 = arg_1.complex_value (); Complex x2 = arg_2.complex_value (); if (! error_state) { ComplexRowVector rv = linspace (x1, x2, npoints); if (! error_state) retval (0) = octave_value (rv, 0); } } else { double x1 = arg_1.double_value (); double x2 = arg_2.double_value (); if (! error_state) { RowVector rv = linspace (x1, x2, npoints); if (! error_state) retval (0) = octave_value (rv, 0); } } } else error ("linspace: npoints must be greater than 2"); } return retval; } static int treat_neg_dim_as_zero (void) { Vtreat_neg_dim_as_zero = check_preference ("treat_neg_dim_as_zero"); return 0; } void symbols_of_data (void) { DEFCONST (I, Complex (0.0, 1.0), 0, 0, "sqrt (-1)"); DEFCONST (Inf, octave_Inf, 0, 0, "infinity"); DEFCONST (J, Complex (0.0, 1.0), 0, 0, "sqrt (-1)"); DEFCONST (NaN, octave_NaN, 0, 0, "not a number"); #if defined (M_E) double e_val = M_E; #else double e_val = exp (1.0); #endif DEFCONST (e, e_val, 0, 0, "exp (1)"); DEFCONST (eps, DBL_EPSILON, 0, 0, "machine precision"); DEFCONST (i, Complex (0.0, 1.0), 1, 0, "sqrt (-1)"); DEFCONST (inf, octave_Inf, 0, 0, "infinity"); DEFCONST (j, Complex (0.0, 1.0), 1, 0, "sqrt (-1)"); DEFCONST (nan, octave_NaN, 0, 0, "not a number"); #if defined (M_PI) double pi_val = M_PI; #else double pi_val = 4.0 * atan (1.0); #endif DEFCONST (pi, pi_val, 0, 0, "ratio of the circumference of a circle to its diameter"); DEFCONST (realmax, DBL_MAX, 0, 0, "realmax (): return largest representable floating point number"); DEFCONST (realmin, DBL_MIN, 0, 0, "realmin (): return smallest representable floating point number"); DEFVAR (treat_neg_dim_as_zero, 0.0, 0, treat_neg_dim_as_zero, "convert negative dimensions to zero"); } /* ;;; Local Variables: *** ;;; mode: C++ *** ;;; End: *** */