Mercurial > octave-dspies
changeset 18660:cce0a63afb3c
maint: Merge away extra head.
| author | Rik <rik@octave.org> |
|---|---|
| date | Fri, 18 Apr 2014 10:11:26 -0700 |
| parents | 4f43f87b7c3e (diff) 75ec138ba53b (current diff) |
| children | 1fc22871bd8b |
| files | |
| diffstat | 19 files changed, 610 insertions(+), 400 deletions(-) [+] |
line wrap: on
line diff
--- a/NEWS Tue Apr 15 14:00:30 2014 -0700 +++ b/NEWS Fri Apr 18 10:11:26 2014 -0700 @@ -9,6 +9,21 @@ has continuous 1st and 2nd derivatives while 'pchip' only has a continuous 1st derivative. + ** Integer formats used in the printf family of functions now work for + 64-bit integers and are more compatible with Matlab when printing + non-integer values. Now instead of truncating, Octave will switch + the effective format to '%g' in the following circumstances: + + * the value of an integer type (int8, uint32, etc.) value exceeds + the maximum for the format specifier. For '%d', the limit is + intmax ('int64') and for '%u' it is intmax ('uint64'). + + * round(x) != x or the value is outside the range allowed by the + integer format specifier. + + There is still one difference: Matlab switches to '%e' and Octave + is currently switching to '%g'. + ** Other new functions added in 4.2: dir_in_loadpath @@ -22,6 +37,7 @@ be removed from Octave 4.6 (or whatever version is the second major release after 4.2): + bicubic find_dir_in_path nfields
--- a/libinterp/corefcn/oct-stream.cc Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/corefcn/oct-stream.cc Fri Apr 18 10:11:26 2014 -0700 @@ -785,11 +785,14 @@ if (i < n) { + // Accept and record modifier, but don't place it in the format + // item text. All integer conversions are handled as 64-bit + // integers. + switch (s[i]) { case 'h': case 'l': case 'L': - modifier = s[i]; - *buf << s[i++]; + modifier = s[i++]; break; default: @@ -2160,7 +2163,7 @@ printf_value_cache (const octave_value_list& args, const std::string& who) : values (args), val_idx (0), elt_idx (0), - n_vals (values.length ()), n_elts (0), data (0), + n_vals (values.length ()), n_elts (0), have_data (false), curr_state (ok) { for (octave_idx_type i = 0; i < values.length (); i++) @@ -2178,7 +2181,7 @@ ~printf_value_cache (void) { } // Get the current value as a double and advance the internal pointer. - double double_value (void); + octave_value get_next_value (void); // Get the current value as an int and advance the internal pointer. int int_value (void); @@ -2197,8 +2200,8 @@ int elt_idx; int n_vals; int n_elts; - const double *data; - NDArray curr_val; + bool have_data; + octave_value curr_val; state curr_state; // Must create value cache with values! @@ -2212,29 +2215,27 @@ printf_value_cache& operator = (const printf_value_cache&); }; -double -printf_value_cache::double_value (void) +octave_value +printf_value_cache::get_next_value (void) { - double retval = 0.0; + octave_value retval; if (exhausted ()) curr_state = conversion_error; while (! exhausted ()) { - if (! data) + if (! have_data) { - octave_value tmp_val = values (val_idx); + curr_val = values (val_idx); // Force string conversion here for compatibility. - curr_val = tmp_val.array_value (true); - if (! error_state) { elt_idx = 0; - n_elts = curr_val.length (); - data = curr_val.data (); + n_elts = curr_val.numel (); + have_data = true; } else { @@ -2245,13 +2246,13 @@ if (elt_idx < n_elts) { - retval = data[elt_idx++]; + retval = curr_val.fast_elem_extract (elt_idx++); if (elt_idx >= n_elts) { elt_idx = 0; val_idx++; - data = 0; + have_data = false; } break; @@ -2259,7 +2260,7 @@ else { val_idx++; - data = 0; + have_data = false; if (n_elts == 0 && exhausted ()) curr_state = conversion_error; @@ -2276,14 +2277,19 @@ { int retval = 0; - double dval = double_value (); + octave_value val = get_next_value (); if (! error_state) { - if (D_NINT (dval) == dval) - retval = NINT (dval); - else - curr_state = conversion_error; + double dval = val.double_value (); + + if (! error_state) + { + if (D_NINT (dval) == dval) + retval = NINT (dval); + else + curr_state = conversion_error; + } } return retval; @@ -2358,37 +2364,202 @@ return retval; } -#define DO_DOUBLE_CONV_1(TYPE) \ - do \ - { \ - if (val > std::numeric_limits<TYPE>::max () \ - || val < std::numeric_limits<TYPE>::min ()) \ - { \ - std::string tfmt = fmt; \ - \ - tfmt.replace (tfmt.rfind (elt->type), 1, ".g"); \ - \ - if (elt->modifier == 'l') \ - tfmt.replace (tfmt.rfind (elt->modifier), 1, ""); \ - \ - retval += do_printf_conv (os, tfmt.c_str (), nsa, sa_1, sa_2, \ - val, who); \ - } \ - else \ - retval += do_printf_conv (os, fmt, nsa, sa_1, sa_2, \ - static_cast<TYPE> (val), who); \ - } \ - while (0) - -#define DO_DOUBLE_CONV(TQUAL) \ - do \ - { \ - if (elt->modifier == 'l') \ - DO_DOUBLE_CONV_1 (TQUAL long); \ - else \ - DO_DOUBLE_CONV_1 (TQUAL int); \ - } \ - while (0) +static bool +is_nan_or_inf (const octave_value& val) +{ + octave_value ov_isnan = val.isnan (); + octave_value ov_isinf = val.isinf (); + + return (ov_isnan.is_true () || ov_isinf.is_true ()); +} + +static bool +ok_for_signed_int_conv (const octave_value& val) +{ + uint64_t limit = std::numeric_limits<int64_t>::max (); + + if (val.is_integer_type ()) + { + if (val.is_uint64_type ()) + { + octave_uint64 ival = val.uint64_scalar_value (); + + if (ival.value () <= limit) + return true; + } + else + return true; + } + else + { + double dval = val.double_value (); + + if (dval == xround (dval) && dval <= limit) + return true; + } + + return false; +} + +static bool +ok_for_unsigned_int_conv (const octave_value& val) +{ + if (val.is_integer_type ()) + { + // Easier than dispatching here... + + octave_value ov_is_ge_zero + = do_binary_op (octave_value::op_ge, val, octave_value (0.0)); + + return ov_is_ge_zero.is_true (); + } + else + { + double dval = val.double_value (); + + uint64_t limit = std::numeric_limits<uint64_t>::max (); + + if (dval == xround (dval) && dval >= 0 && dval <= limit) + return true; + } + + return false; +} + +static std::string +switch_to_g_format (const printf_format_elt *elt) +{ + std::string tfmt = elt->text; + + tfmt.replace (tfmt.rfind (elt->type), 1, "g"); + + return tfmt; +} + +int +octave_base_stream::do_numeric_printf_conv (std::ostream& os, + const printf_format_elt *elt, + int nsa, int sa_1, int sa_2, + const octave_value& val, + const std::string& who) +{ + int retval = 0; + + const char *fmt = elt->text; + + if (is_nan_or_inf (val)) + { + double dval = val.double_value (); + + std::string tfmt = fmt; + std::string::size_type i1, i2; + + tfmt.replace ((i1 = tfmt.rfind (elt->type)), + 1, 1, 's'); + + if ((i2 = tfmt.rfind ('.')) != std::string::npos + && i2 < i1) + { + tfmt.erase (i2, i1-i2); + if (elt->prec < 0) + nsa--; + } + + const char *tval; + if (lo_ieee_isinf (dval)) + { + if (elt->flags.find ('+') != std::string::npos) + tval = (dval < 0 ? "-Inf" : "+Inf"); + else + tval = (dval < 0 ? "-Inf" : "Inf"); + } + else + { + if (elt->flags.find ('+') != std::string::npos) + tval = (lo_ieee_is_NA (dval) ? "+NA" : "+NaN"); + else + tval = (lo_ieee_is_NA (dval) ? "NA" : "NaN"); + } + + retval += do_printf_conv (os, tfmt.c_str (), nsa, sa_1, sa_2, tval, who); + } + else + { + static std::string llmod + = sizeof (long) == sizeof (int64_t) ? "l" : "ll"; + + char type = elt->type; + + switch (type) + { + case 'd': case 'i': case 'c': + if (ok_for_signed_int_conv (val)) + { + octave_int64 tval = val.int64_scalar_value (); + + // Insert "long" modifier. + std::string tfmt = fmt; + tfmt.replace (tfmt.rfind (type), 1, llmod + type); + + retval += do_printf_conv (os, tfmt.c_str (), nsa, sa_1, sa_2, + tval.value (), who); + } + else + { + std::string tfmt = switch_to_g_format (elt); + + double dval = val.double_value (); + + if (! error_state) + retval += do_printf_conv (os, tfmt.c_str (), nsa, + sa_1, sa_2, dval, who); + } + break; + + case 'o': case 'x': case 'X': case 'u': + if (ok_for_unsigned_int_conv (val)) + { + octave_uint64 tval = val.uint64_scalar_value (); + + // Insert "long" modifier. + std::string tfmt = fmt; + tfmt.replace (tfmt.rfind (type), 1, llmod + type); + + retval += do_printf_conv (os, tfmt.c_str (), nsa, sa_1, sa_2, + tval.value (), who); + } + else + { + std::string tfmt = switch_to_g_format (elt); + + double dval = val.double_value (); + + if (! error_state) + retval += do_printf_conv (os, tfmt.c_str (), nsa, + sa_1, sa_2, dval, who); + } + break; + + case 'f': case 'e': case 'E': + case 'g': case 'G': + { + double dval = val.double_value (); + + if (! error_state) + retval += do_printf_conv (os, fmt, nsa, sa_1, sa_2, dval, who); + } + break; + + default: + error ("%s: invalid format specifier", + who.c_str ()); + return -1; + break; + } + } + + return retval; +} int octave_base_stream::do_printf (printf_format_list& fmt_list, @@ -2443,8 +2614,6 @@ } } - const char *fmt = elt->text; - if (elt->type == '%') { os << "%"; @@ -2460,81 +2629,18 @@ std::string val = val_cache.string_value (); if (val_cache) - retval += do_printf_conv (os, fmt, nsa, sa_1, + retval += do_printf_conv (os, elt->text, nsa, sa_1, sa_2, val.c_str (), who); else break; } else { - double val = val_cache.double_value (); + octave_value val = val_cache.get_next_value (); if (val_cache) - { - if (lo_ieee_isnan (val) || xisinf (val)) - { - std::string tfmt = fmt; - std::string::size_type i1, i2; - - tfmt.replace ((i1 = tfmt.rfind (elt->type)), - 1, 1, 's'); - - if ((i2 = tfmt.rfind ('.')) != std::string::npos - && i2 < i1) - { - tfmt.erase (i2, i1-i2); - if (elt->prec < 0) - nsa--; - } - - const char *tval; - if (xisinf (val)) - { - if (elt->flags.find ('+') != std::string::npos) - tval = (val < 0 ? "-Inf" : "+Inf"); - else - tval = (val < 0 ? "-Inf" : "Inf"); - } - else - { - if (elt->flags.find ('+') != std::string::npos) - tval = (lo_ieee_is_NA (val) ? "+NA" : "+NaN"); - else - tval = (lo_ieee_is_NA (val) ? "NA" : "NaN"); - } - - retval += do_printf_conv (os, tfmt.c_str (), - nsa, sa_1, sa_2, - tval, who); - } - else - { - char type = elt->type; - - switch (type) - { - case 'd': case 'i': case 'c': - DO_DOUBLE_CONV (OCTAVE_EMPTY_CPP_ARG); - break; - - case 'o': case 'x': case 'X': case 'u': - DO_DOUBLE_CONV (unsigned); - break; - - case 'f': case 'e': case 'E': - case 'g': case 'G': - retval += do_printf_conv (os, fmt, nsa, - sa_1, sa_2, val, who); - break; - - default: - error ("%s: invalid format specifier", - who.c_str ()); - return -1; - break; - } - } - } + retval += do_numeric_printf_conv (os, elt, nsa, sa_1, + sa_2, val, who); else break; }
--- a/libinterp/corefcn/oct-stream.h Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/corefcn/oct-stream.h Fri Apr 18 10:11:26 2014 -0700 @@ -479,6 +479,11 @@ int flush (void); + int do_numeric_printf_conv (std::ostream& os, const printf_format_elt *elt, + int nsa, int sa_1, int sa_2, + const octave_value& val, + const std::string& who); + int do_printf (printf_format_list& fmt_list, const octave_value_list& args, const std::string& who /* = "printf" */);
--- a/libinterp/octave-value/ov-base-diag.cc Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-base-diag.cc Fri Apr 18 10:11:26 2014 -0700 @@ -530,6 +530,23 @@ template <class DMT, class MT> octave_value +octave_base_diag<DMT, MT>::fast_elem_extract (octave_idx_type n) const +{ + if (n < matrix.numel ()) + { + octave_idx_type nr = matrix.rows (); + + octave_idx_type r = n % nr; + octave_idx_type c = n / nr; + + return octave_value (matrix.elem (r, c)); + } + else + return octave_value (); +} + +template <class DMT, class MT> +octave_value octave_base_diag<DMT, MT>::to_dense (void) const { if (! dense_cache.is_defined ())
--- a/libinterp/octave-value/ov-base-diag.h Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-base-diag.h Fri Apr 18 10:11:26 2014 -0700 @@ -207,6 +207,8 @@ void print_info (std::ostream& os, const std::string& prefix) const; + octave_value fast_elem_extract (octave_idx_type n) const; + protected: DMT matrix;
--- a/libinterp/octave-value/ov-base-scalar.cc Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-base-scalar.cc Fri Apr 18 10:11:26 2014 -0700 @@ -180,6 +180,13 @@ } template <class ST> +octave_value +octave_base_scalar<ST>::fast_elem_extract (octave_idx_type n) const +{ + return (n == 0) ? octave_value (scalar) : octave_value (); +} + +template <class ST> bool octave_base_scalar<ST>::fast_elem_insert_self (void *where, builtin_type_t btyp) const
--- a/libinterp/octave-value/ov-base-scalar.h Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-base-scalar.h Fri Apr 18 10:11:26 2014 -0700 @@ -148,6 +148,8 @@ ST& scalar_ref (void) { return scalar; } + octave_value fast_elem_extract (octave_idx_type n) const; + bool fast_elem_insert_self (void *where, builtin_type_t btyp) const; protected:
--- a/libinterp/octave-value/ov-base-sparse.cc Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-base-sparse.cc Fri Apr 18 10:11:26 2014 -0700 @@ -436,6 +436,20 @@ return success; } + +template <class T> +octave_value +octave_base_sparse<T>::fast_elem_extract (octave_idx_type n) const +{ + octave_idx_type nr = matrix.rows (); + octave_idx_type nc = matrix.cols (); + + octave_idx_type i = n % nr; + octave_idx_type j = n / nr; + + return (i < nr && j < nc) ? octave_value (matrix(i,j)) : octave_value (); +} + template <class T> octave_value octave_base_sparse<T>::map (octave_base_value::unary_mapper_t umap) const
--- a/libinterp/octave-value/ov-base-sparse.h Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-base-sparse.h Fri Apr 18 10:11:26 2014 -0700 @@ -165,6 +165,8 @@ octave_idx_type *mex_get_jc (void) const { return matrix.mex_get_jc (); } + octave_value fast_elem_extract (octave_idx_type n) const; + protected: octave_value map (octave_base_value::unary_mapper_t umap) const;
--- a/libinterp/octave-value/ov-perm.cc Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-perm.cc Fri Apr 18 10:11:26 2014 -0700 @@ -448,3 +448,18 @@ return retval; } +octave_value +octave_perm_matrix::fast_elem_extract (octave_idx_type n) const +{ + if (n < matrix.numel ()) + { + octave_idx_type nr = matrix.rows (); + + octave_idx_type r = n % nr; + octave_idx_type c = n / nr; + + return octave_value (matrix.elem (r, c)); + } + else + return octave_value (); +}
--- a/libinterp/octave-value/ov-perm.h Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-perm.h Fri Apr 18 10:11:26 2014 -0700 @@ -218,6 +218,8 @@ octave_value map (unary_mapper_t umap) const { return to_dense ().map (umap); } + octave_value fast_elem_extract (octave_idx_type n) const; + protected: PermMatrix matrix;
--- a/libinterp/octave-value/ov-range.cc Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-range.cc Fri Apr 18 10:11:26 2014 -0700 @@ -677,6 +677,13 @@ return retval; } +octave_value +octave_range::fast_elem_extract (octave_idx_type n) const +{ + return (n < range.nelem ()) + ? octave_value (range.elem (n)) : octave_value (); +} + DEFUN (allow_noninteger_range_as_index, args, nargout, "-*- texinfo -*-\n\ @deftypefn {Built-in Function} {@var{val} =} allow_noninteger_range_as_index ()\n\
--- a/libinterp/octave-value/ov-range.h Tue Apr 15 14:00:30 2014 -0700 +++ b/libinterp/octave-value/ov-range.h Fri Apr 18 10:11:26 2014 -0700 @@ -290,6 +290,8 @@ return m.map (umap); } + octave_value fast_elem_extract (octave_idx_type n) const; + private: Range range;
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/scripts/deprecated/bicubic.m Fri Apr 18 10:11:26 2014 -0700 @@ -0,0 +1,260 @@ +## Copyright (C) 2005-2013 Hoxide Ma +## +## This file is part of Octave. +## +## Octave is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or (at +## your option) any later version. +## +## Octave is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with Octave; see the file COPYING. If not, see +## <http://www.gnu.org/licenses/>. + +## -*- texinfo -*- +## @deftypefn {Function File} {@var{zi} =} bicubic (@var{x}, @var{y}, @var{z}, @var{xi}, @var{yi}, @var{extrapval}) +## +## @code{bicubic} is deprecated and will be removed in Octave version 4.6. +## Use @code{interp2 (@dots{}, "spline")} for the equivalent functionality. +## +## Return a matrix @var{zi} corresponding to the bicubic +## interpolations at @var{xi} and @var{yi} of the data supplied +## as @var{x}, @var{y} and @var{z}. Points outside the grid are set +## to @var{extrapval}. +## +## See @url{http://wiki.woodpecker.org.cn/moin/Octave/Bicubic} +## for further information. +## @seealso{interp2} +## @end deftypefn + +## Bicubic interpolation method. +## Author: Hoxide Ma <hoxide_dirac@yahoo.com.cn> + +## Deprecated in version 4.2 + +function zi = bicubic (x, y, z, xi, yi, extrapval, spline_alpha) + + persistent warned = false; + if (! warned) + warned = true; + warning ("Octave:deprecated-function", + "bicubic is obsolete and will be removed from a future version of Octave, please use interp2 instead"); + endif + + if (nargin < 1 || nargin > 7) + print_usage (); + endif + + if (nargin == 7 && isscalar (spline_alpha)) + a = spline_alpha; + else + a = 0.5; + endif + + if (nargin < 6) + extrapval = NaN; + endif + + if (isa (x, "single") || isa (y, "single") || isa (z, "single") + || isa (xi, "single") || isa (yi, "single")) + myeps = eps ("single"); + else + myeps = eps (); + endif + + if (nargin <= 2) + ## bicubic (z) or bicubic (z, 2) + if (nargin == 1) + n = 1; + else + n = y; + endif + z = x; + x = []; + [rz, cz] = size (z); + s = linspace (1, cz, (cz-1) * pow2 (n) + 1); + t = linspace (1, rz, (rz-1) * pow2 (n) + 1); + elseif (nargin == 3) + if (! isvector (x) || ! isvector (y)) + error ("bicubic: XI and YI must be vector"); + endif + s = y; + t = z; + z = x; + [rz, cz] = size (z); + elseif (nargin == 5 || nargin == 6) + [rz, cz] = size (z) ; + if (isvector (x) && isvector (y)) + if (rz != length (y) || cz != length (x)) + error ("bicubic: length of X and Y must match the size of Z"); + endif + elseif (size_equal (x, y) && size_equal (x, z)) + x = x(1,:); + y = y(:,1); + else + error ("bicubic: X, Y and Z must be equal size matrices of same size"); + endif + + if (all (diff (x) < 0)) + flipx = true; + x = fliplr (x); + elseif (all (diff (x) > 0)) + flipx = false; + else + error ("bicubic:nonmonotonic", "bicubic: X values must be monotonic"); + endif + if (all (diff (y) < 0)) + flipy = true; + y = flipud (y); + elseif (all (diff (y) > 0)) + flipy = false; + else + error ("bicubic:nonmonotonic", "bicubic: Y values must be monotonic"); + endif + + ## Mark values outside the lookup table. + xfirst_ind = find (xi < x(1)); + xlast_ind = find (xi > x(cz)); + yfirst_ind = find (yi < y(1)); + ylast_ind = find (yi > y(rz)); + ## Set value outside the table preliminary to min max index. + xi(xfirst_ind) = x(1); + xi(xlast_ind) = x(cz); + yi(yfirst_ind) = y(1); + yi(ylast_ind) = y(rz); + + x = reshape (x, 1, cz); + x(cz) *= 1 + sign (x(cz)) * myeps; + if (x(cz) == 0) + x(cz) = myeps; + endif; + xi = reshape (xi, 1, length (xi)); + [m, i] = sort ([x, xi]); + o = cumsum (i <= cz); + xidx = o(find (i > cz)); + + y = reshape (y, rz, 1); + y(rz) *= 1 + sign (y(rz)) * myeps; + if (y(rz) == 0) + y(rz) = myeps; + endif; + yi = reshape (yi, length (yi), 1); + [m, i] = sort ([y; yi]); + o = cumsum (i <= rz); + yidx = o([find(i > rz)]); + + ## Set s and t used follow codes. + s = xidx + ((xi .- x(xidx)) ./ (x(xidx+1) .- x(xidx))); + t = yidx + ((yi - y(yidx)) ./ (y(yidx+1) - y(yidx))); + + if (flipx) + s = fliplr (s); + endif + if (flipy) + t = flipud (t); + endif + else + print_usage (); + endif + + if (rz < 3 || cz < 3) + error ("bicubic: Z at least a 3 by 3 matrices"); + endif + + inds = floor (s); + d = find (s == cz); + s = s - floor (s); + inds(d) = cz-1; + s(d) = 1.0; + + d = []; + indt = floor (t); + d = find (t == rz); + t = t - floor (t); + indt(d) = rz-1; + t(d) = 1.0; + d = []; + + p = zeros (size (z) + 2); + p(2:rz+1,2:cz+1) = z; + p(1,:) = (6*(1-a))*p(2,:) - 3*p(3,:) + (6*a-2)*p(4,:); + p(rz+2,:) = (6*(1-a))*p(rz+1,:) - 3*p(rz,:) + (6*a-2)*p(rz-1,:); + p(:,1) = (6*(1-a))*p(:,2) - 3*p(:,3) + (6*a-2)*p(:,4); + p(:,cz+2) = (6*(1-a))*p(:,cz+1) - 3*p(:,cz) + (6*a-2)*p(:,cz-1); + + ## Calculte the C1(t) C2(t) C3(t) C4(t) and C1(s) C2(s) C3(s) C4(s). + t2 = t.*t; + t3 = t2.*t; + + ct0 = -a .* t3 + (2 * a) .* t2 - a .* t ; # -a G0 + ct1 = (2-a) .* t3 + (-3+a) .* t2 + 1 ; # F0 - a G1 + ct2 = (a-2) .* t3 + (-2 *a + 3) .* t2 + a .* t ; # F1 + a G0 + ct3 = a .* t3 - a .* t2; # a G1 + t = []; t2 = []; t3 = []; + + s2 = s.*s; + s3 = s2.*s; + + cs0 = -a .* s3 + (2 * a) .* s2 - a .*s ; # -a G0 + cs1 = (2-a) .* s3 + (-3 + a) .* s2 + 1 ; # F0 - a G1 + cs2 = (a-2) .* s3 + (-2 *a + 3) .* s2 + a .*s ; # F1 + a G0 + cs3 = a .* s3 - a .* s2; # a G1 + s = []; s2 = []; s3 = []; + + cs0 = cs0([1,1,1,1],:); + cs1 = cs1([1,1,1,1],:); + cs2 = cs2([1,1,1,1],:); + cs3 = cs3([1,1,1,1],:); + + lent = length (ct0); + lens = columns (cs0); + zi = zeros (lent, lens); + + for i = 1:lent + it = indt(i); + int = [it, it+1, it+2, it+3]; + zi(i,:) = ([ct0(i),ct1(i),ct2(i),ct3(i)] + * (p(int,inds) .* cs0 + p(int,inds+1) .* cs1 + + p(int,inds+2) .* cs2 + p(int,inds+3) .* cs3)); + endfor + + ## Set points outside the table to extrapval. + if (! (isempty (xfirst_ind) && isempty (xlast_ind))) + zi(:, [xfirst_ind, xlast_ind]) = extrapval; + endif + if (! (isempty (yfirst_ind) && isempty (ylast_ind))) + zi([yfirst_ind; ylast_ind], :) = extrapval; + endif + +endfunction + + +%!demo +%! clf; +%! colormap ("default"); +%! A = [13,-1,12;5,4,3;1,6,2]; +%! x = [0,1,4]+10; +%! y = [-10,-9,-8]; +%! xi = linspace (min (x), max (x), 17); +%! yi = linspace (min (y), max (y), 26)'; +%! mesh (xi, yi, bicubic (x,y,A,xi,yi)); +%! [x,y] = meshgrid (x,y); +%! hold on; plot3 (x(:),y(:),A(:),"b*"); hold off; + +%!test +%! x = linspace (1, -1, 10); +%! [xx, yy] = meshgrid (x); +%! z = cos (6 * xx) + sin (6 * yy); +%! x = linspace (1, -1, 30); +%! [xx2, yy2] = meshgrid (x); +%! z1 = interp2 (xx, yy, z, xx2, yy2, "spline"); +%! z2 = interp2 (fliplr (xx), flipud (yy), fliplr (flipud(z)), +%! fliplr (xx2), flipud (yy2), "spline"); +%! z2 = fliplr (flipud (z2)); +%! assert (z1, z2, 100 * eps ()) +
--- a/scripts/deprecated/module.mk Tue Apr 15 14:00:30 2014 -0700 +++ b/scripts/deprecated/module.mk Fri Apr 18 10:11:26 2014 -0700 @@ -1,6 +1,7 @@ FCN_FILE_DIRS += deprecated deprecated_FCN_FILES = \ + deprecated/bicubic.m \ deprecated/find_dir_in_path.m \ deprecated/isstr.m \ deprecated/nfields.m
--- a/scripts/general/bicubic.m Tue Apr 15 14:00:30 2014 -0700 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,248 +0,0 @@ -## Copyright (C) 2005-2013 Hoxide Ma -## -## This file is part of Octave. -## -## Octave is free software; you can redistribute it and/or modify it -## under the terms of the GNU General Public License as published by -## the Free Software Foundation; either version 3 of the License, or (at -## your option) any later version. -## -## Octave is distributed in the hope that it will be useful, but -## WITHOUT ANY WARRANTY; without even the implied warranty of -## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -## General Public License for more details. -## -## You should have received a copy of the GNU General Public License -## along with Octave; see the file COPYING. If not, see -## <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## @deftypefn {Function File} {@var{zi} =} bicubic (@var{x}, @var{y}, @var{z}, @var{xi}, @var{yi}, @var{extrapval}) -## -## Return a matrix @var{zi} corresponding to the bicubic -## interpolations at @var{xi} and @var{yi} of the data supplied -## as @var{x}, @var{y} and @var{z}. Points outside the grid are set -## to @var{extrapval}. -## -## See @url{http://wiki.woodpecker.org.cn/moin/Octave/Bicubic} -## for further information. -## @seealso{interp2} -## @end deftypefn - -## Bicubic interpolation method. -## Author: Hoxide Ma <hoxide_dirac@yahoo.com.cn> - -function zi = bicubic (x, y, z, xi, yi, extrapval, spline_alpha) - - if (nargin < 1 || nargin > 7) - print_usage (); - endif - - if (nargin == 7 && isscalar (spline_alpha)) - a = spline_alpha; - else - a = 0.5; - endif - - if (nargin < 6) - extrapval = NaN; - endif - - if (isa (x, "single") || isa (y, "single") || isa (z, "single") - || isa (xi, "single") || isa (yi, "single")) - myeps = eps ("single"); - else - myeps = eps (); - endif - - if (nargin <= 2) - ## bicubic (z) or bicubic (z, 2) - if (nargin == 1) - n = 1; - else - n = y; - endif - z = x; - x = []; - [rz, cz] = size (z); - s = linspace (1, cz, (cz-1) * pow2 (n) + 1); - t = linspace (1, rz, (rz-1) * pow2 (n) + 1); - elseif (nargin == 3) - if (! isvector (x) || ! isvector (y)) - error ("bicubic: XI and YI must be vector"); - endif - s = y; - t = z; - z = x; - [rz, cz] = size (z); - elseif (nargin == 5 || nargin == 6) - [rz, cz] = size (z) ; - if (isvector (x) && isvector (y)) - if (rz != length (y) || cz != length (x)) - error ("bicubic: length of X and Y must match the size of Z"); - endif - elseif (size_equal (x, y) && size_equal (x, z)) - x = x(1,:); - y = y(:,1); - else - error ("bicubic: X, Y and Z must be equal size matrices of same size"); - endif - - if (all (diff (x) < 0)) - flipx = true; - x = fliplr (x); - elseif (all (diff (x) > 0)) - flipx = false; - else - error ("bicubic:nonmonotonic", "bicubic: X values must be monotonic"); - endif - if (all (diff (y) < 0)) - flipy = true; - y = flipud (y); - elseif (all (diff (y) > 0)) - flipy = false; - else - error ("bicubic:nonmonotonic", "bicubic: Y values must be monotonic"); - endif - - ## Mark values outside the lookup table. - xfirst_ind = find (xi < x(1)); - xlast_ind = find (xi > x(cz)); - yfirst_ind = find (yi < y(1)); - ylast_ind = find (yi > y(rz)); - ## Set value outside the table preliminary to min max index. - xi(xfirst_ind) = x(1); - xi(xlast_ind) = x(cz); - yi(yfirst_ind) = y(1); - yi(ylast_ind) = y(rz); - - x = reshape (x, 1, cz); - x(cz) *= 1 + sign (x(cz)) * myeps; - if (x(cz) == 0) - x(cz) = myeps; - endif; - xi = reshape (xi, 1, length (xi)); - [m, i] = sort ([x, xi]); - o = cumsum (i <= cz); - xidx = o(find (i > cz)); - - y = reshape (y, rz, 1); - y(rz) *= 1 + sign (y(rz)) * myeps; - if (y(rz) == 0) - y(rz) = myeps; - endif; - yi = reshape (yi, length (yi), 1); - [m, i] = sort ([y; yi]); - o = cumsum (i <= rz); - yidx = o([find(i > rz)]); - - ## Set s and t used follow codes. - s = xidx + ((xi .- x(xidx)) ./ (x(xidx+1) .- x(xidx))); - t = yidx + ((yi - y(yidx)) ./ (y(yidx+1) - y(yidx))); - - if (flipx) - s = fliplr (s); - endif - if (flipy) - t = flipud (t); - endif - else - print_usage (); - endif - - if (rz < 3 || cz < 3) - error ("bicubic: Z at least a 3 by 3 matrices"); - endif - - inds = floor (s); - d = find (s == cz); - s = s - floor (s); - inds(d) = cz-1; - s(d) = 1.0; - - d = []; - indt = floor (t); - d = find (t == rz); - t = t - floor (t); - indt(d) = rz-1; - t(d) = 1.0; - d = []; - - p = zeros (size (z) + 2); - p(2:rz+1,2:cz+1) = z; - p(1,:) = (6*(1-a))*p(2,:) - 3*p(3,:) + (6*a-2)*p(4,:); - p(rz+2,:) = (6*(1-a))*p(rz+1,:) - 3*p(rz,:) + (6*a-2)*p(rz-1,:); - p(:,1) = (6*(1-a))*p(:,2) - 3*p(:,3) + (6*a-2)*p(:,4); - p(:,cz+2) = (6*(1-a))*p(:,cz+1) - 3*p(:,cz) + (6*a-2)*p(:,cz-1); - - ## Calculte the C1(t) C2(t) C3(t) C4(t) and C1(s) C2(s) C3(s) C4(s). - t2 = t.*t; - t3 = t2.*t; - - ct0 = -a .* t3 + (2 * a) .* t2 - a .* t ; # -a G0 - ct1 = (2-a) .* t3 + (-3+a) .* t2 + 1 ; # F0 - a G1 - ct2 = (a-2) .* t3 + (-2 *a + 3) .* t2 + a .* t ; # F1 + a G0 - ct3 = a .* t3 - a .* t2; # a G1 - t = []; t2 = []; t3 = []; - - s2 = s.*s; - s3 = s2.*s; - - cs0 = -a .* s3 + (2 * a) .* s2 - a .*s ; # -a G0 - cs1 = (2-a) .* s3 + (-3 + a) .* s2 + 1 ; # F0 - a G1 - cs2 = (a-2) .* s3 + (-2 *a + 3) .* s2 + a .*s ; # F1 + a G0 - cs3 = a .* s3 - a .* s2; # a G1 - s = []; s2 = []; s3 = []; - - cs0 = cs0([1,1,1,1],:); - cs1 = cs1([1,1,1,1],:); - cs2 = cs2([1,1,1,1],:); - cs3 = cs3([1,1,1,1],:); - - lent = length (ct0); - lens = columns (cs0); - zi = zeros (lent, lens); - - for i = 1:lent - it = indt(i); - int = [it, it+1, it+2, it+3]; - zi(i,:) = ([ct0(i),ct1(i),ct2(i),ct3(i)] - * (p(int,inds) .* cs0 + p(int,inds+1) .* cs1 - + p(int,inds+2) .* cs2 + p(int,inds+3) .* cs3)); - endfor - - ## Set points outside the table to extrapval. - if (! (isempty (xfirst_ind) && isempty (xlast_ind))) - zi(:, [xfirst_ind, xlast_ind]) = extrapval; - endif - if (! (isempty (yfirst_ind) && isempty (ylast_ind))) - zi([yfirst_ind; ylast_ind], :) = extrapval; - endif - -endfunction - - -%!demo -%! clf; -%! colormap ("default"); -%! A = [13,-1,12;5,4,3;1,6,2]; -%! x = [0,1,4]+10; -%! y = [-10,-9,-8]; -%! xi = linspace (min (x), max (x), 17); -%! yi = linspace (min (y), max (y), 26)'; -%! mesh (xi, yi, bicubic (x,y,A,xi,yi)); -%! [x,y] = meshgrid (x,y); -%! hold on; plot3 (x(:),y(:),A(:),"b*"); hold off; - -%!test -%! x = linspace (1, -1, 10); -%! [xx, yy] = meshgrid (x); -%! z = cos (6 * xx) + sin (6 * yy); -%! x = linspace (1, -1, 30); -%! [xx2, yy2] = meshgrid (x); -%! z1 = interp2 (xx, yy, z, xx2, yy2, "spline"); -%! z2 = interp2 (fliplr (xx), flipud (yy), fliplr (flipud(z)), -%! fliplr (xx2), flipud (yy2), "spline"); -%! z2 = fliplr (flipud (z2)); -%! assert (z1, z2, 100 * eps ()) -
--- a/scripts/general/module.mk Tue Apr 15 14:00:30 2014 -0700 +++ b/scripts/general/module.mk Fri Apr 18 10:11:26 2014 -0700 @@ -7,7 +7,6 @@ general_FCN_FILES = \ general/accumarray.m \ general/accumdim.m \ - general/bicubic.m \ general/bincoeff.m \ general/bitcmp.m \ general/bitget.m \
--- a/scripts/plot/draw/stemleaf.m Tue Apr 15 14:00:30 2014 -0700 +++ b/scripts/plot/draw/stemleaf.m Fri Apr 18 10:11:26 2014 -0700 @@ -185,7 +185,7 @@ hl = xs(hlidx); # lower hinge hu = xs(huidx); # upper hinge h_spread = hu - hl; # h_spread: difference between hinges - step = 1.5*h_spread; # step: 1.5 * h_spread + step = fix(1.5*h_spread); # step: 1.5 * h_spread i_fence_l = hl - step; # inner fences: outside hinges + step o_fence_l = hl - 2*step; # outer fences: outside hinges + 2*step i_fence_h = hu + step; @@ -349,13 +349,13 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " #138|___________________" , -%! " M 69| 52 |" , +%! " #138|___________________" , +%! " M 69| 52 |" , %! " H 35| 30 116| 86" , %! " 1 | -28 146|" , %! " _______" , %! " ______| 129|_______" , -%! " f| -99 245|" , +%! " f| -99 245|" , %! " | 0 0| out" , %! " F| -228 374|" , %! " | 0 0| far" , @@ -389,13 +389,13 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 14|___________________" , -%! " M 7| 22 |" , +%! " # 14|___________________" , +%! " M 7| 22 |" , %! " H 4| 12 42| 30" , %! " 1 | 5 52|" , %! " _______" , %! " ______| 45|_______" , -%! " f| -33 87|" , +%! " f| -33 87|" , %! " | 0 0| out" , %! " F| -78 132|" , %! " | 0 0| far" , @@ -418,13 +418,13 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 14|___________________" , -%! " M 7| -28 |" , +%! " # 14|___________________" , +%! " M 7| -28 |" , %! " H 4| -42 -12| 30" , %! " 1 | -52 -5|" , %! " _______" , %! " ______| 45|_______" , -%! " f| -87 33|" , +%! " f| -87 33|" , %! " | 0 0| out" , %! " F| -132 78|" , %! " | 0 0| far" , @@ -446,15 +446,15 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 15|___________________" , -%! " M 8| 22 |" , +%! " # 15|___________________" , +%! " M 8| 22 |" , %! " H 4| 11 42| 31" , %! " 1 | 0 52|" , %! " _______" , %! " ______| 46|_______" , -%! " f| -35 88|" , +%! " f| -35 88|" , %! " | 0 0| out" , -%! " F| -82 135|" , +%! " F| -81 134|" , %! " | 0 0| far" , %! " " , %! " 0 | 20" , @@ -475,15 +475,15 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 15|___________________" , -%! " M 8| -22 |" , +%! " # 15|___________________" , +%! " M 8| -22 |" , %! " H 4| -42 -11| 31" , %! " 1 | -52 0|" , %! " _______" , %! " ______| 46|_______" , -%! " f| -88 35|" , +%! " f| -88 35|" , %! " | 0 0| out" , -%! " F| -135 82|" , +%! " F| -134 81|" , %! " | 0 0| far" , %! " " , %! " -5 | 2" , @@ -503,15 +503,15 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 5|___________________" , -%! " M 3| 0 |" , +%! " # 5|___________________" , +%! " M 3| 0 |" , %! " H 2| -7 0| 7" , %! " 1 | -9 0|" , %! " _______" , %! " ______| 10|_______" , -%! " f| -17 10|" , +%! " f| -17 10|" , %! " | 0 0| out" , -%! " F| -28 21|" , +%! " F| -27 20|" , %! " | 0 0| far" , %! " " , %! " -0 | 9700" , @@ -527,15 +527,15 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 4|___________________" , -%! " M 2| -7 |" , +%! " # 4|___________________" , +%! " M 2| -7 |" , %! " H 1| -9 0| 9" , %! " 1 | -9 0|" , %! " _______" , %! " ______| 13|_______" , -%! " f| -22 13|" , +%! " f| -22 13|" , %! " | 0 0| out" , -%! " F| -36 27|" , +%! " F| -35 26|" , %! " | 0 0| far" , %! " " , %! " -0 | 970" , @@ -552,7 +552,7 @@ %! " " , %! " Fenced Letter Display" , %! " " , -%! " # 17|___________________" , +%! " # 17|___________________" , %! " M 9| 895 |" , %! " H 5| 795 1499| 704" , %! " 1 | 150 1995|" ,
--- a/scripts/time/datestr.m Tue Apr 15 14:00:30 2014 -0700 +++ b/scripts/time/datestr.m Fri Apr 18 10:11:26 2014 -0700 @@ -252,7 +252,8 @@ df = regexprep (df, '[Ss][Ss]', "%S"); - df = strrep (df, "FFF", sprintf ("%03d", 1000 * (v(i,6) - fix (v(i,6))))); + df = strrep (df, "FFF", sprintf ("%03d", + round (1000 * (v(i,6) - fix (v(i,6)))))); df = strrep (df, 'QQ', sprintf ("Q%d", fix ((v(i,2) + 2) / 3)));