Mercurial > octave
annotate scripts/geometry/griddatan.m @ 28210:bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
* griddata.m: Add v4 method algorithm, documentation notes, and BIST tests.
* griddatan.m: Updated method error messages and added BIST tests.
* NEWS: Announce support for "v4" method in Matlab compatibility section.
author | Nicholas R. Jankowski <jankowskin@asme.org> |
---|---|
date | Thu, 09 Apr 2020 16:21:05 -0400 |
parents | 2e6dc7e2b191 |
children | 289882040316 |
rev | line source |
---|---|
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
1 ######################################################################## |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
2 ## |
27919
1891570abac8
update Octave Project Developers copyright for the new year
John W. Eaton <jwe@octave.org>
parents:
27918
diff
changeset
|
3 ## Copyright (C) 2007-2020 The Octave Project Developers |
27918
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27845
diff
changeset
|
4 ## |
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
5 ## See the file COPYRIGHT.md in the top-level directory of this |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
6 ## distribution or <https://octave.org/copyright/>. |
6823 | 7 ## |
8 ## This file is part of Octave. | |
9 ## | |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
10 ## Octave is free software: you can redistribute it and/or modify it |
6823 | 11 ## under the terms of the GNU General Public License as published by |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
12 ## the Free Software Foundation, either version 3 of the License, or |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
13 ## (at your option) any later version. |
6823 | 14 ## |
15 ## Octave is distributed in the hope that it will be useful, but | |
16 ## WITHOUT ANY WARRANTY; without even the implied warranty of | |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
17 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
18 ## GNU General Public License for more details. |
6823 | 19 ## |
20 ## You should have received a copy of the GNU General Public License | |
7016 | 21 ## along with Octave; see the file COPYING. If not, see |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
22 ## <https://www.gnu.org/licenses/>. |
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
23 ## |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
24 ######################################################################## |
6823 | 25 |
26 ## -*- texinfo -*- | |
20852
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20158
diff
changeset
|
27 ## @deftypefn {} {@var{yi} =} griddatan (@var{x}, @var{y}, @var{xi}) |
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20158
diff
changeset
|
28 ## @deftypefnx {} {@var{yi} =} griddatan (@var{x}, @var{y}, @var{xi}, @var{method}) |
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20158
diff
changeset
|
29 ## @deftypefnx {} {@var{yi} =} griddatan (@var{x}, @var{y}, @var{xi}, @var{method}, @var{options}) |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
30 ## |
6823 | 31 ## Generate a regular mesh from irregular data using interpolation. |
20158
7503499a252b
doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
32 ## |
6823 | 33 ## The function is defined by @code{@var{y} = f (@var{x})}. |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
34 ## The interpolation points are all @var{xi}. |
6823 | 35 ## |
17281
bc924baa2c4e
doc: Add new @qcode macro for code samples which are quoted.
Rik <rik@octave.org>
parents:
14370
diff
changeset
|
36 ## The interpolation method can be @qcode{"nearest"} or @qcode{"linear"}. |
bc924baa2c4e
doc: Add new @qcode macro for code samples which are quoted.
Rik <rik@octave.org>
parents:
14370
diff
changeset
|
37 ## If method is omitted it defaults to @qcode{"linear"}. |
19593
446c46af4b42
strip trailing whitespace from most source files
John W. Eaton <jwe@octave.org>
parents:
17744
diff
changeset
|
38 ## |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
39 ## The optional argument @var{options} is passed directly to Qhull when |
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
40 ## computing the Delaunay triangulation used for interpolation. See |
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
41 ## @code{delaunayn} for information on the defaults and how to pass different |
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
42 ## values. |
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
43 ## @seealso{griddata, griddata3, delaunayn} |
6823 | 44 ## @end deftypefn |
45 | |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
46 function yi = griddatan (x, y, xi, method = "linear", varargin) |
6826 | 47 |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
48 if (nargin < 3) |
6826 | 49 print_usage (); |
6823 | 50 endif |
51 | |
6826 | 52 if (ischar (method)) |
53 method = tolower (method); | |
6823 | 54 endif |
55 | |
6826 | 56 [m, n] = size (x); |
57 [mi, ni] = size (xi); | |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
58 |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
59 if (n != ni || rows (y) != m || columns (y) != 1) |
6823 | 60 error ("griddatan: dimensional mismatch"); |
61 endif | |
62 | |
63 ## triangulate data | |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
64 tri = delaunayn (x, varargin{:}); |
6823 | 65 |
10548
479536c5bb10
Replace lowercase nan with NaN for visual cue in scripts
Rik <code@nomad.inbox5.com>
parents:
8920
diff
changeset
|
66 yi = NaN (mi, 1); |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
67 |
6826 | 68 if (strcmp (method, "nearest")) |
6823 | 69 ## search index of nearest point |
6826 | 70 idx = dsearchn (x, tri, xi); |
19833
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
71 valid = ! isnan (idx); |
6823 | 72 yi(valid) = y(idx(valid)); |
73 | |
6826 | 74 elseif (strcmp (method, "linear")) |
6823 | 75 ## search for every point the enclosing triangle |
6826 | 76 [tri_list, bary_list] = tsearchn (x, tri, xi); |
6823 | 77 |
78 ## only keep the points within triangles. | |
19833
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
79 valid = ! isnan (tri_list); |
27244
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
80 tri_list = tri_list(valid); |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
81 bary_list = bary_list(valid, :); |
6826 | 82 nr_t = rows (tri_list); |
6823 | 83 |
84 ## Use barycentric coordinate of point to calculate yi | |
27244
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
85 if (isscalar (tri_list)) |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
86 ## Special case required by orientation rules for vector/vector index. |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
87 yi(valid) = sum (y(tri(tri_list,:)).' .* bary_list, 2); |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
88 else |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
89 yi(valid) = sum (y(tri(tri_list,:)) .* bary_list, 2); |
27239
bee80e27dcc5
griddatan.m: Fix 4-D interpolation (bug #45542, bug #56515).
Juho Iipponen <juho.iipponen@helsinki.fi>
parents:
26376
diff
changeset
|
90 endif |
6823 | 91 |
28210
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
92 elseif (any (strcmp (method, {"cubic", "v4"}))) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
93 error ('griddata: "%s" METHOD only valid for 2-D inputs using "griddata"', method); |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
94 |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
95 elseif (strcmp (method, "natural")) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
96 ## FIXME: implement missing interpolation method 'natural' for 3-D inputs. |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
97 error ('griddatan: "natural" interpolation METHOD not yet implemented'); |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
98 |
6823 | 99 else |
28210
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
100 error ('griddatan: unknown interpolation METHOD: "%s"', method); |
6823 | 101 endif |
102 | |
103 endfunction | |
104 | |
105 | |
8153
ec0a13863eb7
Only run tests that depend on HDF5 and QHull if Octave was actually
Soren Hauberg <hauberg@gmail.com>
parents:
7017
diff
changeset
|
106 %!testif HAVE_QHULL |
14363
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
107 %! [xx,yy] = meshgrid (linspace (-1,1,32)); |
6823 | 108 %! xi = [xx(:), yy(:)]; |
14363
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
109 %! x = 2*rand (100,2) - 1; |
6823 | 110 %! x = [x;1,1;1,-1;-1,-1;-1,1]; |
14363
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
111 %! y = sin (2 * sum (x.^2,2)); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
112 %! zz = griddatan (x,y,xi,"linear"); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
113 %! zz2 = griddata (x(:,1),x(:,2),y,xi(:,1),xi(:,2),"linear"); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
114 %! assert (zz, zz2, 1e-10); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
115 |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
116 %!testif HAVE_QHULL |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
117 %! [xx,yy] = meshgrid (linspace (-1,1,32)); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
118 %! xi = [xx(:), yy(:)]; |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
119 %! x = 2*rand (100,2) - 1; |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
120 %! x = [x;1,1;1,-1;-1,-1;-1,1]; |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
121 %! y = sin (2*sum (x.^2,2)); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
122 %! zz = griddatan (x,y,xi,"nearest"); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
123 %! zz2 = griddata (x(:,1),x(:,2),y,xi(:,1),xi(:,2),"nearest"); |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
124 %! assert (zz, zz2, 1e-10); |
27244
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
125 |
28011
b31a118729ed
update bug status for tests
John W. Eaton <jwe@octave.org>
parents:
27978
diff
changeset
|
126 %!testif HAVE_QHULL <*56515> |
27244
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
127 %! x = [ 0, 0; 1, 1; 0, 1; 1, 0 ]; |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
128 %! y = [ 1; 2; 3; 4 ]; |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
129 %! xi = [ .5, .5 ]; |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
130 %! yi = griddatan (x, y, xi); |
28210
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
131 |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
132 ## Test input validation |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
133 %!error griddatan () |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
134 %!error griddatan (1) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
135 %!error griddatan (1,2) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
136 %!error griddatan (1,2,3) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
137 %!error <OPTIONS argument must be a string> griddatan (1,2,3,4,5) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
138 %!error <unknown interpolation METHOD> griddatan (1,2,3,4) |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
139 %!#error <"v4" METHOD only valid for 2-D inputs> |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
140 %! griddatan (ones(2,2,2), 2*ones(2,2,2), 3, "v4") |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
141 %!error <"cubic" METHOD only valid for> griddatan (1, 2, 3, "cubic") |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
142 %!error <"natural" .* not yet implemented> griddatan (1, 2, 3, "natural") |
bb929d5a34cb
griddata.m: Added support for "v4" biharmonic spline interpolation method (bug #33539).
Nicholas R. Jankowski <jankowskin@asme.org>
parents:
28027
diff
changeset
|
143 %!error <unknown interpolation METHOD: "foobar"> griddatan (1, 2, 3, "foobar") |