Mercurial > hg > octave-nkf

annotate scripts/general/interp2.m @ 9399:a5f6b5800f86

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
fix bug in recursive lookup
author Jaroslav Hajek <highegg@gmail.com>
date Fri, 26 Jun 2009 07:13:28 +0200
parents 1bf0ce0930be
children 4f4873f6f873
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
8920
eb63fbe60fab update copyright notices
John W. Eaton <jwe@octave.org>
parents: 8828
diff changeset
1 ## Copyright (C) 2000, 2006, 2007, 2008 Kai Habel
eb63fbe60fab update copyright notices
John W. Eaton <jwe@octave.org>
parents: 8828
diff changeset
2 ## Copyright (C) 2009 Jaroslav Hajek
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
3 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
4 ## This file is part of Octave.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
5 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
6 ## Octave is free software; you can redistribute it and/or modify it
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
7 ## under the terms of the GNU General Public License as published by
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
8 ## the Free Software Foundation; either version 3 of the License, or (at
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
9 ## your option) any later version.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
10 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
11 ## Octave is distributed in the hope that it will be useful, but
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
12 ## WITHOUT ANY WARRANTY; without even the implied warranty of
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
14 ## General Public License for more details.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
15 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
16 ## You should have received a copy of the GNU General Public License
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
17 ## along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
18 ## <http://www.gnu.org/licenses/>.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
19
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
20 ## -*- texinfo -*-
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
21 ## @deftypefn {Function File} {@var{zi} =} interp2 (@var{x}, @var{y}, @var{z}, @var{xi}, @var{yi})
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
22 ## @deftypefnx {Function File} {@var{zi} =} interp2 (@var{Z}, @var{xi}, @var{yi})
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
23 ## @deftypefnx {Function File} {@var{zi} =} interp2 (@var{Z}, @var{n})
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
24 ## @deftypefnx {Function File} {@var{zi} =} interp2 (@dots{}, @var{method})
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
25 ## @deftypefnx {Function File} {@var{zi} =} interp2 (@dots{}, @var{method}, @var{extrapval})
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
26 ##
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
27 ## Two-dimensional interpolation. @var{x}, @var{y} and @var{z} describe a
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
28 ## surface function. If @var{x} and @var{y} are vectors their length
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
29 ## must correspondent to the size of @var{z}. @var{x} and @var{y} must be
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
30 ## monotonic. If they are matrices they must have the @code{meshgrid}
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
31 ## format.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
32 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
33 ## @table @code
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
34 ## @item interp2 (@var{x}, @var{y}, @var{Z}, @var{xi}, @var{yi}, @dots{})
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
35 ## Returns a matrix corresponding to the points described by the
8478
5aafd61e4c3c [docs] @var{XI}, @var{YI}. => @var{xi}, @var{yi}.
Brian Gough <bjg@gnu.org>
parents: 8347
diff changeset
36 ## matrices @var{xi}, @var{yi}.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
37 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
38 ## If the last argument is a string, the interpolation method can
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
39 ## be specified. The method can be 'linear', 'nearest' or 'cubic'.
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
40 ## If it is omitted 'linear' interpolation is assumed.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
41 ##
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
42 ## @item interp2 (@var{z}, @var{xi}, @var{yi})
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
43 ## Assumes @code{@var{x} = 1:rows (@var{z})} and @code{@var{y} =
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
44 ## 1:columns (@var{z})}
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
45 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
46 ## @item interp2 (@var{z}, @var{n})
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
47 ## Interleaves the matrix @var{z} n-times. If @var{n} is omitted a value
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
48 ## of @code{@var{n} = 1} is assumed.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
49 ## @end table
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
50 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
51 ## The variable @var{method} defines the method to use for the
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
52 ## interpolation. It can take one of the following values
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
53 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
54 ## @table @asis
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
55 ## @item 'nearest'
6218
c1b66b74937d [project @ 2006-12-27 17:43:50 by jwe]
jwe
parents: 6157
diff changeset
56 ## Return the nearest neighbor.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
57 ## @item 'linear'
6218
c1b66b74937d [project @ 2006-12-27 17:43:50 by jwe]
jwe
parents: 6157
diff changeset
58 ## Linear interpolation from nearest neighbors.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
59 ## @item 'pchip'
6218
c1b66b74937d [project @ 2006-12-27 17:43:50 by jwe]
jwe
parents: 6157
diff changeset
60 ## Piece-wise cubic hermite interpolating polynomial (not implemented yet).
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
61 ## @item 'cubic'
6218
c1b66b74937d [project @ 2006-12-27 17:43:50 by jwe]
jwe
parents: 6157
diff changeset
62 ## Cubic interpolation from four nearest neighbors.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
63 ## @item 'spline'
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
64 ## Cubic spline interpolation--smooth first and second derivatives
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
65 ## throughout the curve.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
66 ## @end table
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
67 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
68 ## If a scalar value @var{extrapval} is defined as the final value, then
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
69 ## values outside the mesh as set to this value. Note that in this case
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
70 ## @var{method} must be defined as well. If @var{extrapval} is not
6742
ebf96cc00ee9 [project @ 2007-06-18 16:27:52 by jwe]
jwe
parents: 6702
diff changeset
71 ## defined then NA is assumed.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
72 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
73 ## @seealso{interp1}
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
74 ## @end deftypefn
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
75
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
76 ## Author: Kai Habel <kai.habel@gmx.de>
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
77 ## 2005-03-02 Thomas Weber <weber@num.uni-sb.de>
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
78 ## * Add test cases
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
79 ## 2005-03-02 Paul Kienzle <pkienzle@users.sf.net>
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
80 ## * Simplify
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
81 ## 2005-04-23 Dmitri A. Sergatskov <dasergatskov@gmail.com>
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
82 ## * Modified demo and test for new gnuplot interface
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
83 ## 2005-09-07 Hoxide <hoxide_dirac@yahoo.com.cn>
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
84 ## * Add bicubic interpolation method
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
85 ## * Fix the eat line bug when the last element of XI or YI is
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
86 ## negative or zero.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
87 ## 2005-11-26 Pierre Baldensperger <balden@libertysurf.fr>
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
88 ## * Rather big modification (XI,YI no longer need to be
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
89 ## "meshgridded") to be consistent with the help message
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
90 ## above and for compatibility.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
91
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
92 function ZI = interp2 (varargin)
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
93 Z = X = Y = XI = YI = n = [];
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
94 method = "linear";
6742
ebf96cc00ee9 [project @ 2007-06-18 16:27:52 by jwe]
jwe
parents: 6702
diff changeset
95 extrapval = NA;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
96
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
97 switch (nargin)
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
98 case 1
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
99 Z = varargin{1};
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
100 case 2
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
101 if (ischar (varargin{2}))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
102 [Z, method] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
103 else
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
104 [Z, n] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
105 endif
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
106 case 3
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
107 if (ischar (varargin{3}))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
108 [Z, n, method] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
109 else
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
110 [Z, XI, YI] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
111 endif
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
112 case 4
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
113 if (ischar (varargin{4}))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
114 [Z, XI, YI, method] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
115 else
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
116 [Z, n, method, extrapval] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
117 endif
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
118 case 5
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
119 if (ischar (varargin{4}))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
120 [Z, XI, YI, method, extrapval] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
121 else
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
122 [X, Y, Z, XI, YI] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
123 endif
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
124 case 6
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
125 [X, Y, Z, XI, YI, method] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
126 case 7
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
127 [X, Y, Z, XI, YI, method, extrapval] = deal (varargin{:});
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
128 otherwise
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
129 print_usage ();
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
130 endswitch
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
131
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
132 ## Type checking.
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
133 if (!ismatrix (Z))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
134 error ("interp2 expected matrix Z");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
135 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
136 if (!isempty (n) && !isscalar (n))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
137 error ("interp2 expected scalar n");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
138 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
139 if (!ischar (method))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
140 error ("interp2 expected string 'method'");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
141 endif
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
142 if (ischar (extrapval) || strcmp (extrapval, "extrap"))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
143 extrapval = [];
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
144 elseif (!isscalar (extrapval))
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
145 error ("interp2 expected n extrapval");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
146 endif
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
147
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
148 ## Define X, Y, XI, YI if needed
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
149 [zr, zc] = size (Z);
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
150 if (isempty (X))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
151 X = 1:zc;
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
152 Y = 1:zr;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
153 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
154 if (! isnumeric (X) || ! isnumeric (Y))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
155 error ("interp2 expected numeric X, Y");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
156 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
157 if (! isempty (n))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
158 p = 2^n;
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
159 XI = (p:p*zc)/p;
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
160 YI = (p:p*zr)'/p;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
161 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
162 if (! isnumeric (XI) || ! isnumeric (YI))
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
163 error ("interp2 expected numeric XI, YI");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
164 endif
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
165
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
166
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
167 if (strcmp (method, "linear") || strcmp (method, "nearest") ...
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
168 || strcmp (method, "pchip"))
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
169
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
170 ## If X and Y vectors produce a grid from them
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
171 if (isvector (X) && isvector (Y))
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
172 X = X(:); Y = Y(:);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
173 elseif (size_equal (X, Y))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
174 X = X(1,:)'; Y = Y(:,1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
175 else
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
176 error ("X and Y must be matrices of same size");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
177 endif
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
178 if (columns (Z) != length (X) || rows (Z) != length (Y))
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
179 error ("X and Y size must match Z dimensions");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
180 endif
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
181
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
182 ## If Xi and Yi are vectors of different orientation build a grid
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
183 if ((rows (XI) == 1 && columns (YI) == 1)
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
184 || (columns (XI) == 1 && rows (YI) == 1))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
185 [XI, YI] = meshgrid (XI, YI);
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
186 elseif (! size_equal (XI, YI))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
187 error ("XI and YI must be matrices of same size");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
188 endif
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
189
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
190 ## if XI, YI are vectors, X and Y should share their orientation.
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
191 if (rows (XI) == 1)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
192 if (rows (X) != 1)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
193 X = X.';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
194 endif
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
195 if (rows (Y) != 1)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
196 Y = Y.';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
197 endif
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
198 elseif (columns (XI) == 1)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
199 if (columns (X) != 1)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
200 X = X.';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
201 endif
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
202 if (columns (Y) != 1)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
203 Y = Y.';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
204 endif
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
205 endif
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
206
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
207 xidx = lookup (X, XI, "lr");
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
208 yidx = lookup (Y, YI, "lr");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
209
9399
a5f6b5800f86 fix bug in recursive lookup
Jaroslav Hajek <highegg@gmail.com>
parents: 9051
diff changeset
210 if (min (xidx) <= 0)
a5f6b5800f86 fix bug in recursive lookup
Jaroslav Hajek <highegg@gmail.com>
parents: 9051
diff changeset
211 save bad_data X XI
a5f6b5800f86 fix bug in recursive lookup
Jaroslav Hajek <highegg@gmail.com>
parents: 9051
diff changeset
212 endif
a5f6b5800f86 fix bug in recursive lookup
Jaroslav Hajek <highegg@gmail.com>
parents: 9051
diff changeset
213
a5f6b5800f86 fix bug in recursive lookup
Jaroslav Hajek <highegg@gmail.com>
parents: 9051
diff changeset
214
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
215 if (strcmp (method, "linear"))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
216 ## each quad satisfies the equation z(x,y)=a+b*x+c*y+d*xy
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
217 ##
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
218 ## a-b
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
219 ## | |
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
220 ## c-d
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
221 a = Z(1:(zr - 1), 1:(zc - 1));
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
222 b = Z(1:(zr - 1), 2:zc) - a;
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
223 c = Z(2:zr, 1:(zc - 1)) - a;
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
224 d = Z(2:zr, 2:zc) - a - b - c;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
225
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
226 idx = sub2ind (size (a), yidx, xidx);
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
227
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
228 ## scale XI, YI values to a 1-spaced grid
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
229 Xsc = (XI - X(xidx)) ./ (X(xidx + 1) - X(xidx));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
230 Ysc = (YI - Y(yidx)) ./ (Y(yidx + 1) - Y(yidx));
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
231
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
232 ## apply plane equation
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
233 ZI = a(idx) + b(idx).*Xsc + c(idx).*Ysc + d(idx).*Xsc.*Ysc;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
234
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
235 elseif (strcmp (method, "nearest"))
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
236 ii = (XI - X(xidx) > X(xidx + 1) - XI);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
237 jj = (YI - Y(yidx) > Y(yidx + 1) - YI);
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
238 idx = sub2ind (size (Z), yidx+jj, xidx+ii);
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
239 ZI = Z(idx);
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
240
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
241 elseif (strcmp (method, "pchip"))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
242
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
243 if (length (X) < 2 || length (Y) < 2)
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
244 error ("interp2: pchip2 requires at least 2 points in each dimension")
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
245 endif
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
246
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
247 ## first order derivatives
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
248 DX = __pchip_deriv__ (X, Z, 2);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
249 DY = __pchip_deriv__ (Y, Z, 1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
250 ## Compute mixed derivatives row-wise and column-wise, use the average.
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
251 DXY = (__pchip_deriv__ (X, DY, 2) + __pchip_deriv__ (Y, DX, 1))/2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
252
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
253 ## do the bicubic interpolation
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
254 hx = diff (X); hx = hx(xidx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
255 hy = diff (Y); hy = hy(yidx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
256
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
257 tx = (XI - X(xidx)) ./ hx;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
258 ty = (YI - Y(yidx)) ./ hy;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
259
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
260 ## construct the cubic hermite base functions in x, y
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
261
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
262 ## formulas:
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
263 ## b{1,1} = ( 2*t.^3 - 3*t.^2 + 1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
264 ## b{2,1} = h.*( t.^3 - 2*t.^2 + t );
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
265 ## b{1,2} = (-2*t.^3 + 3*t.^2 );
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
266 ## b{2,2} = h.*( t.^3 - t.^2 );
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
267
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
268 ## optimized equivalents of the above:
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
269 t1 = tx.^2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
270 t2 = tx.*t1 - t1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
271 xb{2,2} = hx.*t2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
272 t1 = t2 - t1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
273 xb{2,1} = hx.*(t1 + tx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
274 t2 += t1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
275 xb{1,2} = -t2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
276 xb{1,1} = t2 + 1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
277
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
278 t1 = ty.^2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
279 t2 = ty.*t1 - t1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
280 yb{2,2} = hy.*t2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
281 t1 = t2 - t1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
282 yb{2,1} = hy.*(t1 + ty);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
283 t2 += t1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
284 yb{1,2} = -t2;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
285 yb{1,1} = t2 + 1;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
286
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
287 ZI = zeros (size (XI));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
288 for i = 1:2
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
289 for j = 1:2
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
290 zidx = sub2ind (size (Z), yidx+(j-1), xidx+(i-1));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
291 ZI += xb{1,i} .* yb{1,j} .* Z(zidx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
292 ZI += xb{2,i} .* yb{1,j} .* DX(zidx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
293 ZI += xb{1,i} .* yb{2,j} .* DY(zidx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
294 ZI += xb{2,i} .* yb{2,j} .* DXY(zidx);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
295 endfor
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
296 endfor
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
297
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
298 endif
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
299
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
300 if (! isempty (extrapval))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
301 ## set points outside the table to 'extrapval'
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
302 if (X (1) < X (end))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
303 if (Y (1) < Y (end))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
304 ZI (XI < X(1,1) | XI > X(end) | YI < Y(1,1) | YI > Y(end)) = ...
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
305 extrapval;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
306 else
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
307 ZI (XI < X(1) | XI > X(end) | YI < Y(end) | YI > Y(1)) = ...
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
308 extrapval;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
309 endif
6979
2883ea1c5c18 [project @ 2007-10-08 20:23:48 by dbateman]
dbateman
parents: 6742
diff changeset
310 else
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
311 if (Y (1) < Y (end))
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
312 ZI (XI < X(end) | XI > X(1) | YI < Y(1) | YI > Y(end)) = ...
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
313 extrapval;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
314 else
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
315 ZI (XI < X(1,end) | XI > X(1) | YI < Y(end) | YI > Y(1)) = ...
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
316 extrapval;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
317 endif
6979
2883ea1c5c18 [project @ 2007-10-08 20:23:48 by dbateman]
dbateman
parents: 6742
diff changeset
318 endif
2883ea1c5c18 [project @ 2007-10-08 20:23:48 by dbateman]
dbateman
parents: 6742
diff changeset
319 endif
2883ea1c5c18 [project @ 2007-10-08 20:23:48 by dbateman]
dbateman
parents: 6742
diff changeset
320
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
321 else
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
322
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
323 ## If X and Y vectors produce a grid from them
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
324 if (isvector (X) && isvector (Y))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
325 X = X(:).';
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
326 Y = Y(:);
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
327 if (!isequal ([length(X), length(Y)], size(Z)))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
328 error ("X and Y size must match Z dimensions");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
329 endif
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
330 elseif (!size_equal (X, Y))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
331 error ("X and Y must be matrices of same size");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
332 if (! size_equal (X, Z))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
333 error ("X and Y size must match Z dimensions");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
334 endif
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
335 endif
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
336
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
337 ## If Xi and Yi are vectors of different orientation build a grid
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
338 if ((rows (XI) == 1 && columns (YI) == 1)
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
339 || (columns (XI) == 1 && rows (YI) == 1))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
340 ## Do nothing
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
341 elseif (! size_equal (XI, YI))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
342 error ("XI and YI must be matrices of same size");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
343 endif
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
344
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
345 ## FIXME bicubic/__splinen__ don't handle arbitrary XI, YI
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
346 if (strcmp (method, "cubic"))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
347 ZI = bicubic (X, Y, Z, XI(1,:), YI(:,1), extrapval);
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
348
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
349 elseif (strcmp (method, "spline"))
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
350 ZI = __splinen__ ({Y(:,1).', X(1,:)}, Z, {YI(:,1), XI(1,:)}, extrapval,
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
351 "spline");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
352 else
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
353 error ("interpolation method not recognized");
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
354 endif
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
355
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
356 endif
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
357 endfunction
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
358
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
359 %!demo
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
360 %! A=[13,-1,12;5,4,3;1,6,2];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
361 %! x=[0,1,4]; y=[10,11,12];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
362 %! xi=linspace(min(x),max(x),17);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
363 %! yi=linspace(min(y),max(y),26)';
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
364 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'linear'));
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
365 %! [x,y] = meshgrid(x,y);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
366 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
367
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
368 %!demo
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
369 %! [x,y,A] = peaks(10);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
370 %! x = x(1,:)'; y = y(:,1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
371 %! xi=linspace(min(x),max(x),41);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
372 %! yi=linspace(min(y),max(y),41)';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
373 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'linear'));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
374 %! [x,y] = meshgrid(x,y);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
375 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
376
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
377 %!demo
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
378 %! A=[13,-1,12;5,4,3;1,6,2];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
379 %! x=[0,1,4]; y=[10,11,12];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
380 %! xi=linspace(min(x),max(x),17);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
381 %! yi=linspace(min(y),max(y),26)';
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
382 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'nearest'));
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
383 %! [x,y] = meshgrid(x,y);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
384 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
385
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
386 %!demo
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
387 %! [x,y,A] = peaks(10);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
388 %! x = x(1,:)'; y = y(:,1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
389 %! xi=linspace(min(x),max(x),41);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
390 %! yi=linspace(min(y),max(y),41)';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
391 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'nearest'));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
392 %! [x,y] = meshgrid(x,y);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
393 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
394
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
395 %!demo
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
396 %! A=[13,-1,12;5,4,3;1,6,2];
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
397 %! x=[0,1,2]; y=[10,11,12];
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
398 %! xi=linspace(min(x),max(x),17);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
399 %! yi=linspace(min(y),max(y),26)';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
400 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'pchip'));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
401 %! [x,y] = meshgrid(x,y);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
402 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
403
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
404 %!demo
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
405 %! [x,y,A] = peaks(10);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
406 %! x = x(1,:)'; y = y(:,1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
407 %! xi=linspace(min(x),max(x),41);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
408 %! yi=linspace(min(y),max(y),41)';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
409 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'pchip'));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
410 %! [x,y] = meshgrid(x,y);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
411 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
412
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
413 %!demo
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
414 %! A=[13,-1,12;5,4,3;1,6,2];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
415 %! x=[0,1,2]; y=[10,11,12];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
416 %! xi=linspace(min(x),max(x),17);
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
417 %! yi=linspace(min(y),max(y),26)';
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
418 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'cubic'));
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
419 %! [x,y] = meshgrid(x,y);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
420 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
421
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
422 %!demo
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
423 %! [x,y,A] = peaks(10);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
424 %! x = x(1,:)'; y = y(:,1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
425 %! xi=linspace(min(x),max(x),41);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
426 %! yi=linspace(min(y),max(y),41)';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
427 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'cubic'));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
428 %! [x,y] = meshgrid(x,y);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
429 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
430
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
431 %!demo
6702
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
432 %! A=[13,-1,12;5,4,3;1,6,2];
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
433 %! x=[0,1,2]; y=[10,11,12];
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
434 %! xi=linspace(min(x),max(x),17);
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
435 %! yi=linspace(min(y),max(y),26)';
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
436 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'spline'));
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
437 %! [x,y] = meshgrid(x,y);
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
438 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
b2391d403ed2 [project @ 2007-06-12 21:39:26 by dbateman]
dbateman
parents: 6448
diff changeset
439
8712
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
440 %!demo
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
441 %! [x,y,A] = peaks(10);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
442 %! x = x(1,:)'; y = y(:,1);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
443 %! xi=linspace(min(x),max(x),41);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
444 %! yi=linspace(min(y),max(y),41)';
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
445 %! mesh(xi,yi,interp2(x,y,A,xi,yi,'spline'));
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
446 %! [x,y] = meshgrid(x,y);
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
447 %! hold on; plot3(x(:),y(:),A(:),"b*"); hold off;
010e15c7be9a support pchip method in interp2
Jaroslav Hajek <highegg@gmail.com>
parents: 8479
diff changeset
448
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
449 %!test % simple test
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
450 %! x = [1,2,3];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
451 %! y = [4,5,6,7];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
452 %! [X, Y] = meshgrid(x,y);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
453 %! Orig = X.^2 + Y.^3;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
454 %! xi = [1.2,2, 1.5];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
455 %! yi = [6.2, 4.0, 5.0]';
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
456 %!
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
457 %! Expected = ...
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
458 %! [243, 245.4, 243.9;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
459 %! 65.6, 68, 66.5;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
460 %! 126.6, 129, 127.5];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
461 %! Result = interp2(x,y,Orig, xi, yi);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
462 %!
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
463 %! assert(Result, Expected, 1000*eps);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
464
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
465 %!test % 2^n form
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
466 %! x = [1,2,3];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
467 %! y = [4,5,6,7];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
468 %! [X, Y] = meshgrid(x,y);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
469 %! Orig = X.^2 + Y.^3;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
470 %! xi = [1:0.25:3]; yi = [4:0.25:7]';
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
471 %! Expected = interp2(x,y,Orig, xi, yi);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
472 %! Result = interp2(Orig,2);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
473 %!
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
474 %! assert(Result, Expected, 10*eps);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
475
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
476 %!test % matrix slice
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
477 %! A = eye(4);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
478 %! assert(interp2(A,[1:4],[1:4]),[1,1,1,1]);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
479
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
480 %!test % non-gridded XI,YI
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
481 %! A = eye(4);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
482 %! assert(interp2(A,[1,2;3,4],[1,3;2,4]),[1,0;0,1]);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
483
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
484 %!test % for values outside of boundaries
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
485 %! x = [1,2,3];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
486 %! y = [4,5,6,7];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
487 %! [X, Y] = meshgrid(x,y);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
488 %! Orig = X.^2 + Y.^3;
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
489 %! xi = [0,4];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
490 %! yi = [3,8]';
6742
ebf96cc00ee9 [project @ 2007-06-18 16:27:52 by jwe]
jwe
parents: 6702
diff changeset
491 %! assert(interp2(x,y,Orig, xi, yi),[NA,NA;NA,NA]);
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
492 %! assert(interp2(x,y,Orig, xi, yi,'linear', 0),[0,0;0,0]);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
493
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
494 %!test % for values at boundaries
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
495 %! A=[1,2;3,4];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
496 %! x=[0,1];
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
497 %! y=[2,3]';
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
498 %! assert(interp2(x,y,A,x,y,'linear'), A);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
499 %! assert(interp2(x,y,A,x,y,'nearest'), A);
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
500