Mercurial > octave
annotate scripts/geometry/griddatan.m @ 27923:bd51beb6205e
update formatting of copyright notices
* Use <https://octave.org/copyright/> instead of
<https://octave.org/COPYRIGHT.html/>.
* For consistency with other comments in the Octave sources, use
C++-style comments for copyright blocks in C and C++ files.
* Use delimiters above and below copyright blocks that are appropriate
for the language used in the file.
* Eliminate extra spacing inside copyright blocks.
* lex.ll (looks_like_copyright): Also allow newlines and carriage
returns before the word "Copyright".
* scripts/mk-doc.pl (gethelp): Also skip empty comment lines.
* bp-table.cc, type.m: Adjust tests.
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Wed, 08 Jan 2020 11:59:41 -0500 |
parents | 1891570abac8 |
children | a4268efb7334 |
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 | |
7017 | 46 ## Author: David Bateman <dbateman@free.fr> |
47 | |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
48 function yi = griddatan (x, y, xi, method = "linear", varargin) |
6826 | 49 |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
50 if (nargin < 3) |
6826 | 51 print_usage (); |
6823 | 52 endif |
53 | |
6826 | 54 if (ischar (method)) |
55 method = tolower (method); | |
6823 | 56 endif |
57 | |
6826 | 58 [m, n] = size (x); |
59 [mi, ni] = size (xi); | |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
60 |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
61 if (n != ni || rows (y) != m || columns (y) != 1) |
6823 | 62 error ("griddatan: dimensional mismatch"); |
63 endif | |
64 | |
27845
25f914321e7b
Use different Qhull options for griddataX to get better precision (bug #50494).
Rik <rik@octave.org>
parents:
27244
diff
changeset
|
65 ## Bug #50494 (loss of precision unless non-default Qhull options used) |
25f914321e7b
Use different Qhull options for griddataX to get better precision (bug #50494).
Rik <rik@octave.org>
parents:
27244
diff
changeset
|
66 if (isempty (varargin) && n <= 3) |
25f914321e7b
Use different Qhull options for griddataX to get better precision (bug #50494).
Rik <rik@octave.org>
parents:
27244
diff
changeset
|
67 varargin{1} = {"Qt", "Qbb", "Qc"}; |
25f914321e7b
Use different Qhull options for griddataX to get better precision (bug #50494).
Rik <rik@octave.org>
parents:
27244
diff
changeset
|
68 endif |
25f914321e7b
Use different Qhull options for griddataX to get better precision (bug #50494).
Rik <rik@octave.org>
parents:
27244
diff
changeset
|
69 |
6823 | 70 ## triangulate data |
14370
fbdee844550c
griddatan.m: Restore ability to pass options to underlying Qhull
Rik <octave@nomad.inbox5.com>
parents:
14363
diff
changeset
|
71 tri = delaunayn (x, varargin{:}); |
6823 | 72 |
10548
479536c5bb10
Replace lowercase nan with NaN for visual cue in scripts
Rik <code@nomad.inbox5.com>
parents:
8920
diff
changeset
|
73 yi = NaN (mi, 1); |
11587
c792872f8942
all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
74 |
6826 | 75 if (strcmp (method, "nearest")) |
6823 | 76 ## search index of nearest point |
6826 | 77 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
|
78 valid = ! isnan (idx); |
6823 | 79 yi(valid) = y(idx(valid)); |
80 | |
6826 | 81 elseif (strcmp (method, "linear")) |
6823 | 82 ## search for every point the enclosing triangle |
6826 | 83 [tri_list, bary_list] = tsearchn (x, tri, xi); |
6823 | 84 |
85 ## 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
|
86 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
|
87 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
|
88 bary_list = bary_list(valid, :); |
6826 | 89 nr_t = rows (tri_list); |
6823 | 90 |
91 ## 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
|
92 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
|
93 ## 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
|
94 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
|
95 else |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
96 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
|
97 endif |
6823 | 98 |
99 else | |
11472
1740012184f9
Use uppercase for variable names in error() strings to match Info documentation. Only m-files done.
Rik <octave@nomad.inbox5.com>
parents:
10548
diff
changeset
|
100 error ("griddatan: unknown interpolation 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 |
3f354ef16400
griddatan.m: Allow a single query point as input (bug #45542, bug #56515).
Rik <rik@octave.org>
parents:
27239
diff
changeset
|
126 %!testif HAVE_QHULL <56515> |
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); |