summaryrefslogtreecommitdiffstats
path: root/scilab
diff options
context:
space:
mode:
authorSamuel GOUGEON <sgougeon@free.fr>2021-01-18 07:25:45 +0100
committerClément David <clement.david@esi-group.com>2021-03-17 09:02:50 +0100
commitaa892b9fc8def5f8452312b893d4de72e1b1e2cb (patch)
tree6a84044cf686b062348125a10e730cf172662b6b /scilab
parentc55427777028825aae0e05210c01fd859daf1e21 (diff)
downloadscilab-aa892b9fc8def5f8452312b893d4de72e1b1e2cb.zip
scilab-aa892b9fc8def5f8452312b893d4de72e1b1e2cb.tar.gz
* Bug 16629 fixed: interp1() fixed + complex + extended extrapolations
http://bugzilla.scilab.org/16629 interp1 page fixed & updated (PDF): http://bugzilla.scilab.org/attachment.cgi?id=5217 Change-Id: Iec59c6c2ffa312d6d282c52828b40966bc8871d9
Diffstat (limited to 'scilab')
-rw-r--r--scilab/CHANGES.md7
-rw-r--r--scilab/modules/helptools/etc/images_md5.txt2
-rw-r--r--scilab/modules/helptools/images/cell_view.pngbin0 -> 6206 bytes
-rw-r--r--scilab/modules/helptools/images/interp1_1.pngbin14171 -> 16355 bytes
-rw-r--r--scilab/modules/helptools/images/library_view.pngbin0 -> 10345 bytes
-rw-r--r--scilab/modules/helptools/images/list_view.pngbin0 -> 4769 bytes
-rw-r--r--scilab/modules/helptools/images/structure_view.pngbin0 -> 7300 bytes
-rw-r--r--scilab/modules/interpolation/help/en_US/interp1.xml357
-rw-r--r--scilab/modules/interpolation/help/ja_JP/interp1.xml200
-rw-r--r--scilab/modules/interpolation/help/pt_BR/interp1.xml182
-rw-r--r--scilab/modules/interpolation/help/ru_RU/interp1.xml356
-rw-r--r--scilab/modules/interpolation/macros/interp1.sci704
-rw-r--r--scilab/modules/interpolation/tests/unit_tests/interp1.tst654
13 files changed, 1587 insertions, 875 deletions
diff --git a/scilab/CHANGES.md b/scilab/CHANGES.md
index 3c27968..136da4a 100644
--- a/scilab/CHANGES.md
+++ b/scilab/CHANGES.md
@@ -200,6 +200,9 @@ Feature changes and additions
200* `unique` is enabled for any 2D sparse arrays, in simple, 'c' and 'r' modes. 200* `unique` is enabled for any 2D sparse arrays, in simple, 'c' and 'r' modes.
201* `%chars` constant added, to easily access to some selected sets of unicode symbols. 201* `%chars` constant added, to easily access to some selected sets of unicode symbols.
202* Lists are displayed in a more compact and comprehensive way. 202* Lists are displayed in a more compact and comprehensive way.
203* `interp1` is upgraded:
204 - complex numbers `y` now supported: the real and imaginary parts are interpolated separately.
205 - extrapolation option extended: `edgevalue` mode added for all interpolations; `periodic` mode added for linear and spline interpolations. `linear` mode added for the spline interpolations.
203 206
204Help pages: 207Help pages:
205----------- 208-----------
@@ -354,10 +357,14 @@ Bug Fixes
354* [#16553](https://bugzilla.scilab.org/16553): `unique(["" ""])` returned `["" ""]`. 357* [#16553](https://bugzilla.scilab.org/16553): `unique(["" ""])` returned `["" ""]`.
355* [#16557](https://bugzilla.scilab.org/16557): `macr2tree` + `tree2code` translated `e={2}` into `"e=1"` and `e={2,"ab"}` into `"e=[2,"ab"]"`. 358* [#16557](https://bugzilla.scilab.org/16557): `macr2tree` + `tree2code` translated `e={2}` into `"e=1"` and `e={2,"ab"}` into `"e=[2,"ab"]"`.
356* [#16559](https://bugzilla.scilab.org/16553): `isempty(A)` was true for sparse matrix of dimension 2^16 or larger. 359* [#16559](https://bugzilla.scilab.org/16553): `isempty(A)` was true for sparse matrix of dimension 2^16 or larger.
360<<<<<<< HEAD
357* [#16622](https://bugzilla.scilab.org/16622): `inv` could no longer be overloaded for hypermatrices of decimal or complex numbers. 361* [#16622](https://bugzilla.scilab.org/16622): `inv` could no longer be overloaded for hypermatrices of decimal or complex numbers.
358* [#16623](https://bugzilla.scilab.org/16623): `rand(2,2,2)^2` yielded a wrong result instead of trying to call the `%s_p_s` overload for input hypermatrices. 362* [#16623](https://bugzilla.scilab.org/16623): `rand(2,2,2)^2` yielded a wrong result instead of trying to call the `%s_p_s` overload for input hypermatrices.
359* [#16644](https://bugzilla.scilab.org/16644): `input("message:")` yielded a wrong error message about `mprintf` in case of non-interpretable input. 363* [#16644](https://bugzilla.scilab.org/16644): `input("message:")` yielded a wrong error message about `mprintf` in case of non-interpretable input.
360* [#16654](https://bugzilla.scilab.org/16654): `interp` was leaking memory. 364* [#16654](https://bugzilla.scilab.org/16654): `interp` was leaking memory.
365=======
366* [#16629](https://bugzilla.scilab.org/16629): `interp1`'s documentation did not tell the spline edges conditions ; extrapolation modes were poorly explained ; the description of the result's size was completely wrong ; x as an option was not documented. A wrong extrapolation value could silently return a wrong result. There was some dead code like `if varargin(5)==%nan`. A bugged error message yielded its own error. When x is implicit, the argument index in error messages could be wrong. `periodic` and `edgevalue` extrapolation modes were not available. `linear` extrapolation was not available for splines. When `xp` is an hypermatrix with `size(xp,1)==1`, the size of the result was irregular/wrong.
367>>>>>>> ae66d0ee461... * Bug 16629 fixed: interp1() fixed + complex + extended extrapolations
361 368
362 369
363### Bugs fixed in 6.1.0: 370### Bugs fixed in 6.1.0:
diff --git a/scilab/modules/helptools/etc/images_md5.txt b/scilab/modules/helptools/etc/images_md5.txt
index 459a3ae..674e87c 100644
--- a/scilab/modules/helptools/etc/images_md5.txt
+++ b/scilab/modules/helptools/etc/images_md5.txt
@@ -802,7 +802,7 @@ hsv2rgb_1.png=7d20c259e94301d9763fbddb7bff4784
802hsvcolormap_1.png=11918d88bcc793200af0b9f1b58b0554 802hsvcolormap_1.png=11918d88bcc793200af0b9f1b58b0554
803iir_1.png=e675c2755f68ddc4611c849895b63012 803iir_1.png=e675c2755f68ddc4611c849895b63012
804intdec_1.png=da3896e7d2e8468dd33edf57ccbe4480 804intdec_1.png=da3896e7d2e8468dd33edf57ccbe4480
805interp1_1.png=0e9a3f4319b2818ce4921b9bc3008d80 805interp1_1.png=e76a8e0a6070fbe078060844707c49fe
806interp2d_1.png=f4af61bc3faf493d778169ec7decc7ae 806interp2d_1.png=f4af61bc3faf493d778169ec7decc7ae
807interp2d_2.png=a62363849a3a9a3571dc1942fcebfbd3 807interp2d_2.png=a62363849a3a9a3571dc1942fcebfbd3
808interp2d_3.png=a9f2a0b0e59ac8cedae3ad22e0b2cc46 808interp2d_3.png=a9f2a0b0e59ac8cedae3ad22e0b2cc46
diff --git a/scilab/modules/helptools/images/cell_view.png b/scilab/modules/helptools/images/cell_view.png
new file mode 100644
index 0000000..2905e6c
--- /dev/null
+++ b/scilab/modules/helptools/images/cell_view.png
Binary files differ
diff --git a/scilab/modules/helptools/images/interp1_1.png b/scilab/modules/helptools/images/interp1_1.png
index 7002a47..2b0659b 100644
--- a/scilab/modules/helptools/images/interp1_1.png
+++ b/scilab/modules/helptools/images/interp1_1.png
Binary files differ
diff --git a/scilab/modules/helptools/images/library_view.png b/scilab/modules/helptools/images/library_view.png
new file mode 100644
index 0000000..d3df3ad
--- /dev/null
+++ b/scilab/modules/helptools/images/library_view.png
Binary files differ
diff --git a/scilab/modules/helptools/images/list_view.png b/scilab/modules/helptools/images/list_view.png
new file mode 100644
index 0000000..0a9d969
--- /dev/null
+++ b/scilab/modules/helptools/images/list_view.png
Binary files differ
diff --git a/scilab/modules/helptools/images/structure_view.png b/scilab/modules/helptools/images/structure_view.png
new file mode 100644
index 0000000..106a0a9
--- /dev/null
+++ b/scilab/modules/helptools/images/structure_view.png
Binary files differ
diff --git a/scilab/modules/interpolation/help/en_US/interp1.xml b/scilab/modules/interpolation/help/en_US/interp1.xml
index 1408549..2b6b370 100644
--- a/scilab/modules/interpolation/help/en_US/interp1.xml
+++ b/scilab/modules/interpolation/help/en_US/interp1.xml
@@ -2,8 +2,8 @@
2<!-- 2<!--
3 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab 3 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
4 * Copyright (C) 2008 - INRIA - Farid BELAHCENE 4 * Copyright (C) 2008 - INRIA - Farid BELAHCENE
5 *
6 * Copyright (C) 2012 - 2016 - Scilab Enterprises 5 * Copyright (C) 2012 - 2016 - Scilab Enterprises
6 * Copyright (C) 2021 - Samuel GOUGEON - Le Mans Université
7 * 7 *
8 * This file is hereby licensed under the terms of the GNU GPL v2.0, 8 * This file is hereby licensed under the terms of the GNU GPL v2.0,
9 * pursuant to article 5.3.4 of the CeCILL v.2.1. 9 * pursuant to article 5.3.4 of the CeCILL v.2.1.
@@ -19,49 +19,182 @@
19 xmlns:scilab="http://www.scilab.org" xml:id="interp1" xml:lang="en"> 19 xmlns:scilab="http://www.scilab.org" xml:id="interp1" xml:lang="en">
20 <refnamediv> 20 <refnamediv>
21 <refname>interp1</refname> 21 <refname>interp1</refname>
22 <refpurpose>one_dimension interpolation function</refpurpose> 22 <refpurpose>1D interpolation in nearest, linear or spline mode</refpurpose>
23 </refnamediv> 23 </refnamediv>
24 <refsynopsisdiv> 24 <refsynopsisdiv>
25 <title>Syntax</title> 25 <title>Syntax</title>
26 <synopsis> 26 <synopsis>
27 yp = interp1(y, xp)
27 yp = interp1(x, y, xp) 28 yp = interp1(x, y, xp)
28 yp = interp1(x, y, xp, method) 29 yp = interp1(.., xp, method)
29 yp = interp1(x, y, xp, method, extrapolation) 30 yp = interp1(.., xp, method, extrapolation)
30 </synopsis> 31 </synopsis>
31 </refsynopsisdiv> 32 </refsynopsisdiv>
32 <refsection> 33 <refsection>
33 <title>Arguments</title> 34 <title>Arguments</title>
34 <variablelist> 35 <variablelist>
35 <varlistentry> 36 <varlistentry>
36 <term>xp</term> 37 <term>x</term>
37 <listitem> 38 <listitem>
38 <para>reals scalar, vector or matrix (or hypermatrix)</para> 39 vector of at least 2 real numbers: Abscissas of known interpolation nodes,
40 without duplicates. By default,
41 <itemizedlist>
42 <listitem>
43 if <varname>y</varname> is a vector: <literal>x=1:length(y)</literal>.
44 </listitem>
45 <listitem>
46 if <varname>y</varname> is a matrix or an hypermatrix:
47 <literal>x=1:size(y,1)</literal>.
48 </listitem>
49 </itemizedlist>
50 <para/>
39 </listitem> 51 </listitem>
40 </varlistentry> 52 </varlistentry>
41 <varlistentry> 53 <varlistentry>
42 <term>x</term> 54 <term>y</term>
43 <listitem> 55 <listitem>
44 <para>reals vector</para> 56 vector, matrix or hypermatrix of real or complex numbers: values at known
57 interpolation nodes, at the corresponding <varname>x</varname> abscissas.
58 <itemizedlist>
59 <listitem>
60 if <varname>y</varname> is a vector, <varname>x</varname> and
61 <varname>y</varname> must have the same length.
62 </listitem>
63 <listitem>
64 if <varname>y</varname> is a matrix or an hypermatrix, we must have
65 <literal>length(x)==size(y,1)</literal>. Each column of
66 <varname>y</varname> is then interpolated versus the same
67 <varname>x</varname> abscissas, for the given <varname>xp</varname>.
68 </listitem>
69 </itemizedlist>
70 <para/>
45 </listitem> 71 </listitem>
46 </varlistentry> 72 </varlistentry>
47 <varlistentry> 73 <varlistentry>
48 <term>method</term> 74 <term>xp</term>
49 <listitem> 75 <listitem>
50 <para>(optional) string defining the interpolation method</para> 76 scalar, vector, matrix or hypermatrix or decimal numbers: abscissas of
77 points whose values <varname>yp</varname> must be computed according
78 to data of interpolation nodes.
79 <para/>
51 </listitem> 80 </listitem>
52 </varlistentry> 81 </varlistentry>
53 <varlistentry> 82 <varlistentry>
54 <term>extrapolation</term> 83 <term>yp</term>
55 <listitem> 84 <listitem>
56 <para>(optional) string, or real value defining the yp(j) components 85 vector, matrix, or hypermatrix of numbers: interpolated <varname>y</varname>
57 for xp(j) values outside [x1,xn] interval. 86 values at the given <varname>xp</varname>.
58 </para> 87 <itemizedlist>
88 <listitem>
89 if <varname>y</varname> is a vector: <varname>yp</varname> has the
90 size of <varname>xp</varname>.
91 <para/>
92 </listitem>
93 <listitem>
94 if <varname>y</varname> is a matrix or an hypermatrix:
95 <itemizedlist>
96 <listitem>
97 if <varname>xp</varname> is a scalar or a vector:
98 <literal>size(yp)</literal> is <literal>[length(xp) size(y)(2:$)]</literal>
99 </listitem>
100 <listitem>
101 if <varname>xp</varname> is a matrix or an hypermatrix:
102 <literal>size(yp)</literal> is <literal>[size(xp) size(y)(2:$)]</literal>
103 </listitem>
104 </itemizedlist>
105 </listitem>
106 </itemizedlist>
107 <para/>
59 </listitem> 108 </listitem>
60 </varlistentry> 109 </varlistentry>
61 <varlistentry> 110 <varlistentry>
62 <term>yp</term> 111 <term>method</term>
63 <listitem> 112 <listitem>
64 <para>vector, or matrix (or hypermatrix)</para> 113 string defining the interpolation method. Possible values and processing are:
114 <para/>
115 <table>
116 <tr>
117 <td><emphasis role="bold">"linear"</emphasis>:</td>
118 <td>
119 linear interpolation between consecutive nodes, used by default.
120 </td>
121 </tr>
122 <tr>
123 <td><emphasis role="bold">"spline"</emphasis>:</td>
124 <td>
125 interpolation by cubic splines
126 </td>
127 </tr>
128 <tr>
129 <td><emphasis role="bold">"nearest"</emphasis>:</td>
130 <td>
131 <para>for each value xp(j), yp(j) takes the value or y(i)
132 corresponding to x(i) the nearest neighbor of xp(j)
133 </para>
134 </td>
135 </tr>
136 </table>
137 <para/>
138 </listitem>
139 </varlistentry>
140 <varlistentry>
141 <term>extrapolation</term>
142 <listitem>
143 string or number defining the yp(j) components for xp(j) values outside the
144 [x(1)=min(x),x($)=max(x)] interval. We suppose here-below that <varname>x</varname>
145 and <varname>y</varname> have already been sorted accordingly.
146 <para/>
147 <table>
148 <tr>
149 <td><emphasis role="bold">"extrap"</emphasis>:</td>
150 <td>
151 <literal>interp1(x,y,xp, method, "extrap")</literal> is equivalent
152 to <literal>interp1(x,y,xp, method, method)</literal>.
153 </td>
154 </tr>
155 <tr>
156 <td><emphasis role="bold">"linear"</emphasis>:</td>
157 <td>
158 Can be used with the <literal>"spline"</literal> (and obviously
159 <literal>"linear"</literal>) interpolation methods.
160 </td>
161 </tr>
162 <tr>
163 <td><emphasis role="bold">"periodic"</emphasis>:</td>
164 <td>
165 This extrapolation type can be used with the
166 <literal>"linear"</literal> or <literal>"spline"</literal>
167 interpolation methods. Then: if <varname>y</varname> is a
168 vector, <literal>y(1)==y($)</literal> is required ; otherwise
169 <literal>y(1,:)==y($,:)</literal> is required.
170 </td>
171 </tr>
172 <tr>
173 <td><emphasis role="bold">"edgevalue"</emphasis>:</td>
174 <td>
175 Then <literal>yp(i)=y(1)</literal> for every
176 <literal>xp(i)&lt;x(1)</literal>, and
177 <literal>yp(i)=y($)</literal> for every <literal>xp(i)>x($)</literal>.
178 </td>
179 </tr>
180 <tr>
181 <td><emphasis role="bold">padding</emphasis>:</td>
182 <td>
183 <varname>padding</varname> is a decimal or complex number
184 used to set <literal>yp(i)=padding</literal> for every
185 <literal>xp(i) ∉ [min(x),max(x)]</literal>. Example:
186 <literal>yi=interp1(x,y,xp,method, 0)</literal>.
187 </td>
188 </tr>
189 <tr>
190 <td><emphasis role="bold">(none)</emphasis>:</td>
191 <td>
192 By default, the extrapolation is performed by splines when
193 splines are used for the interpolation, and by padding
194 with %nan when the interpolation is linear or by "nearest" node.
195 </td>
196 </tr>
197 </table>
65 </listitem> 198 </listitem>
66 </varlistentry> 199 </varlistentry>
67 </variablelist> 200 </variablelist>
@@ -69,109 +202,87 @@
69 <refsection> 202 <refsection>
70 <title>Description</title> 203 <title>Description</title>
71 <para> 204 <para>
72 Given <literal>(x,y,xp)</literal>, this function performs the yp 205 Given <literal>(x,y,xp)</literal>, this function computes the yp
73 components corresponding to xp by the interpolation(linear by default) 206 components corresponding to xp by the interpolation between known
74 defined by x and y. 207 data provided by (x,y) nodes.
75 </para> 208 </para>
76 <para> 209 <para>
77 If <literal>yp</literal> is a vector then the length of xp is equal 210 <varname>x</varname> is priorly sorted in ascending order, and
78 to the length of <literal>yp,</literal> if <literal>yp</literal> is a 211 <varname>y</varname> values or per column are then sorted accordingly.
79 matrix then <literal>xp</literal> have the same length than the length of
80 each columns of yp, if <literal>yp</literal> is a hypermatrix then the
81 length of <literal>xp</literal> have the same length than the first
82 dimension of <literal>yp. </literal>
83 </para>
84 <para>If size(y)=[C,M1,M2,M3,...,Mj] and size(xp)=[N1,N2,N3,...,Nk] then
85 size(yp)=[N1,N2,..,Nk,M1,M2,...Mj] and length of x must be equal to
86 size(y,1)
87 </para> 212 </para>
88 <para> 213 <para>
89 The <literal>method</literal> parameter sets the evaluation rule for 214 <emphasis role="bold">Interpolation of complex values</emphasis>: When
90 interpolation 215 <varname>y</varname> is complex, its real and imaginary parts are interpolated
216 separately, and then added to build the complex <varname>yp</varname>.
91 </para> 217 </para>
92 <variablelist>
93 <varlistentry>
94 <term>"linear"</term>
95 <listitem>
96 <para>
97 the interpolation is defined by linear method (see <link linkend="interpln">interpln</link>)
98 </para>
99 </listitem>
100 </varlistentry>
101 <varlistentry>
102 <term>"spline"</term>
103 <listitem>
104 <para>the interpolation is defined by cubic spline interpolation (
105 see <link linkend="splin">splin</link> , <link linkend="interp">interp</link>)
106 </para>
107 </listitem>
108 </varlistentry>
109 <varlistentry>
110 <term>"nearest"</term>
111 <listitem>
112 <para>for each value xp(j), yp(j) takes the value or y(i)
113 corresponding to x(i) the nearest neighbor of xp(j)
114 </para>
115 </listitem>
116 </varlistentry>
117 </variablelist>
118 <para> 218 <para>
119 The <literal>extrapolation</literal> parameter sets the evaluation 219 <emphasis role="bold">interp1(x,y,xp,"nearest")</emphasis>: For any xp
120 rule for extrapolation, i.e for <literal>xp(i) </literal>not in [x1,xn] 220 at the middle of an <literal>[x(i),x(i+1)]</literal> interval, the upper
121 interval 221 bound <literal>x(i+1)</literal> is considered as the nearest
222 <varname>x</varname> value, and <literal>yp=y(i+1)</literal> is assigned.
122 </para> 223 </para>
123 <variablelist> 224 <refsect3>
124 <varlistentry> 225 <title>linear interpolations</title>
125 <term>"extrap"</term> 226 They are performed through the <literal>linear_interpn(..)</literal> function,
126 <listitem> 227 with the corresponding <literal>"edgevalue"→"C0"</literal>,
127 <para>the extrapolation is performed by the defined method. 228 <literal>"linear"→"natural"</literal>, <literal>"periodic"→"periodic"</literal>
128 yp=interp1(x,y,xp,method,"extrap") 229 extrapolation option.
129 </para> 230 </refsect3>
130 </listitem> 231 <refsect3>
131 </varlistentry> 232 <title>spline interpolations</title>
132 <varlistentry> 233 <para>
133 <term>real value</term> 234 <emphasis role="bold">interp1(..,xp,"spline")</emphasis> or
134 <listitem> 235 <emphasis role="bold">interp1(..,xp,"spline","spline")</emphasis> or
135 <para>you can choose a real value for extrapolation, in this way 236 <emphasis role="bold">interp1(..,xp,"spline","extrap")</emphasis>
136 yp(i) takes this value for xp(i) not in [x1,xn] interval, for 237 use <literal>not_a_knot</literal> edges conditions. Extrapolation is performed
137 example 0 (but also nan or inf). yi=interp1(x,y,xp,method, 0) 238 by using both spline polynomials computed at the <literal>(x,y)</literal> edges.
138 </para> 239 </para>
139 </listitem> 240 <para>
140 </varlistentry> 241 <emphasis role="bold">interp1(..,xp,"spline","edgevalue")</emphasis>
141 <varlistentry> 242 uses <literal>not_a_knot</literal> edges conditions and then calls
142 <term>by default</term> 243 <literal>interp(..,"C0")</literal> to perform the actual interpolation
143 <listitem> 244 and extrapolation.
144 <para>the extrapolation is performed by the defined method (for 245 </para>
145 spline method), and by nan for linear and nearest methods. 246 <para>
146 yp=interp1(x,y,xp,method) 247 <emphasis role="bold">interp1(..,xp,"spline","periodic")</emphasis>
147 </para> 248 calls both <literal>splin(..)</literal> and then <literal>interp(..)</literal>
148 </listitem> 249 with their <literal>"periodic"</literal> option.
149 </varlistentry> 250 </para>
150 </variablelist> 251 <para>
252 <emphasis role="bold">interp1(..,xp,"spline","linear")</emphasis>
253 calls <literal>splin(..,"natural")</literal> for linear edges conditions,
254 and then feeds <literal>interp(..,"linear")</literal>.
255 </para>
256 </refsect3>
151 </refsection> 257 </refsection>
152 <refsection> 258 <refsection>
153 <title>Examples</title> 259 <title>Examples</title>
154 <programlisting role="example"><![CDATA[ 260 <programlisting role="example"><![CDATA[
155x = linspace(0, 3, 20); 261x = linspace(0, 10, 11)';
156y = x.^2; 262y = sin(x);
157xx = linspace(0, 3, 100); 263xx = linspace(0,10,1000)';
158yy1 = interp1(x, y, xx, 'linear'); 264yy2 = interp1(x, y, xx, 'linear');
159yy2 = interp1(x, y, xx, 'spline'); 265yy1 = interp1(x, y, xx, 'nearest');
160yy3 = interp1(x, y, xx, 'nearest'); 266yy3 = interp1(x, y, xx, 'spline');
161plot(xx', [yy1; yy2; yy3]', x, y, '*') 267clf
162xtitle('interpolation of square function') 268h = plot(xx, [yy1 yy2 yy3], x, y, '.')
163legend(['linear','spline','nearest'], "in_upper_left"); 269h(1).mark_size = 8;
270title "Interpolation of a poorly sampled sin() function" fontsize 3
271legend(['nearest','linear','spline','nodes'], "in_lower_left");
164 ]]></programlisting> 272 ]]></programlisting>
165 <scilab:image> 273 <scilab:image>
166 x=linspace(0,3,20); 274 x = linspace(0, 10, 11)';
167 y=x.^2; 275 y = sin(x);
168 xx=linspace(0,3,100); 276 xx = linspace(0,10,1000)';
169 yy1=interp1(x,y,xx,'linear'); 277 yy2 = interp1(x, y, xx, 'linear');
170 yy2=interp1(x,y,xx,'spline'); 278 yy1 = interp1(x, y, xx, 'nearest');
171 yy3=interp1(x,y,xx,'nearest'); 279 yy3 = interp1(x, y, xx, 'spline');
172 plot(xx',[yy1;yy2;yy3]',x,y,'*') 280 clf
173 xtitle('interpolation of square function') 281 h = plot(xx, [yy1 yy2 yy3], x, y, '.')
174 legend(['linear','spline','nearest'],a=2) 282 h(1).mark_size = 8;
283 title "Interpolation of a poorly sampled sin() function" fontsize 3
284 legend(['nearest','linear','spline','nodes'], "in_lower_left");
285 gcf().axes_size = [600 400];
175 </scilab:image> 286 </scilab:image>
176 </refsection> 287 </refsection>
177 <refsection role="see also"> 288 <refsection role="see also">
@@ -181,10 +292,10 @@ legend(['linear','spline','nearest'], "in_upper_left");
181 <link linkend="interp">interp</link> 292 <link linkend="interp">interp</link>
182 </member> 293 </member>
183 <member> 294 <member>
184 <link linkend="interpln">interpln</link> 295 <link linkend="splin">splin</link>
185 </member> 296 </member>
186 <member> 297 <member>
187 <link linkend="splin">splin</link> 298 <link linkend="linear_interpn">linear_interpn</link>
188 </member> 299 </member>
189 </simplelist> 300 </simplelist>
190 </refsection> 301 </refsection>
@@ -192,8 +303,34 @@ legend(['linear','spline','nearest'], "in_upper_left");
192 <title>History</title> 303 <title>History</title>
193 <revhistory> 304 <revhistory>
194 <revision> 305 <revision>
195 <revnumber>5.4.0</revnumber> 306 <revnumber>6.1.1</revnumber>
196 <revremark>previously, imaginary part of input arguments were implicitly ignored.</revremark> 307 <revremark>
308 <itemizedlist>
309 <listitem>
310 For complex <varname>y</varname> values, <literal>imag(y)</literal> is
311 no longer ignored: <literal>real(y)</literal> and <literal>imag(y)</literal>
312 parts are separately interpolated.
313 </listitem>
314 <listitem>
315 <literal>"periodic"</literal> extrapolation added for the linear and
316 spline interpolations.
317 </listitem>
318 <listitem>
319 <literal>"edgevalue"</literal> extrapolation added for all nearest,
320 linear and spline interpolations.
321 </listitem>
322 <listitem>
323 <literal>"linear"</literal> extrapolation added for the spline
324 interpolation.
325 </listitem>
326 <listitem>
327 When <varname>xp</varname> is an hypermatrix and
328 <literal>size(xp,1)==1</literal>, <literal>size(yp)</literal>
329 is now always <literal>[size(xp) size(y)(2,$)</literal> instead of
330 <literal>[size(xp,2:$), size(y)(2,$)</literal>.
331 </listitem>
332 </itemizedlist>
333 </revremark>
197 </revision> 334 </revision>
198 </revhistory> 335 </revhistory>
199 </refsection> 336 </refsection>
diff --git a/scilab/modules/interpolation/help/ja_JP/interp1.xml b/scilab/modules/interpolation/help/ja_JP/interp1.xml
deleted file mode 100644
index 547ad0c..0000000
--- a/scilab/modules/interpolation/help/ja_JP/interp1.xml
+++ /dev/null
@@ -1,200 +0,0 @@
1<?xml version="1.0" encoding="UTF-8"?>
2<!--
3 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
4 * Copyright (C) 2008 - INRIA - Farid BELAHCENE
5 *
6 * Copyright (C) 2012 - 2016 - Scilab Enterprises
7 *
8 * This file is hereby licensed under the terms of the GNU GPL v2.0,
9 * pursuant to article 5.3.4 of the CeCILL v.2.1.
10 * This file was originally licensed under the terms of the CeCILL v2.1,
11 * and continues to be available under such terms.
12 * For more information, see the COPYING file which you should have received
13 * along with this program.
14 *
15 -->
16<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink"
17 xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml"
18 xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook"
19 xmlns:scilab="http://www.scilab.org" xml:id="interp1" xml:lang="ja">
20 <refnamediv>
21 <refname>interp1</refname>
22 <refpurpose>一次元補間関数</refpurpose>
23 </refnamediv>
24 <refsynopsisdiv>
25 <title>呼び出し手順</title>
26 <synopsis>
27 yp = interp1(x, y, xp)
28 yp = interp1(x, y, xp, method)
29 yp = interp1(x, y, xp, method, extrapolation)
30 </synopsis>
31 </refsynopsisdiv>
32 <refsection>
33 <title>パラメータ</title>
34 <variablelist>
35 <varlistentry>
36 <term>xp</term>
37 <listitem>
38 <para>実数スカラー, ベクトルまたは行列 (またはハイパー行列)</para>
39 </listitem>
40 </varlistentry>
41 <varlistentry>
42 <term>x</term>
43 <listitem>
44 <para>実数ベクトル</para>
45 </listitem>
46 </varlistentry>
47 <varlistentry>
48 <term>method</term>
49 <listitem>
50 <para>(オプションの) 補間方法を定義する文字列</para>
51 </listitem>
52 </varlistentry>
53 <varlistentry>
54 <term>extrapolation</term>
55 <listitem>
56 <para>(オプションの) 文字列, または領域 [x1,xn] の外側の xp(j)の値について
57 yp(j) の要素を定義する実数値.
58 </para>
59 </listitem>
60 </varlistentry>
61 <varlistentry>
62 <term>yp</term>
63 <listitem>
64 <para>ベクトルまたは行列 (またはハイパー行列)</para>
65 </listitem>
66 </varlistentry>
67 </variablelist>
68 </refsection>
69 <refsection>
70 <title>説明</title>
71 <para>
72 <literal>(x,y,xp)</literal>を指定すると, この関数は
73 xp に対応する yp の要素を x および y により定義された
74 (デフォルトでは線形)補間により計算します.
75 </para>
76 <para>
77 <literal>yp</literal> がベクトルの場合,
78 xp の長さは <literal>yp</literal> の長さに等しくなります.
79 <literal>yp</literal> が行列の場合,
80 <literal>xp</literal> は yp の各列の長さと同じ長さとなります.
81 <literal>yp</literal> がハイパー行列の場合,
82 <literal>xp</literal> の長さは<literal>yp</literal> の最初の次元と同じ
83 になります.
84 </para>
85 <para>size(y)=[C,M1,M2,M3,...,Mj] かつ size(xp)=[N1,N2,N3,...,Nk] の場合,
86 size(yp)=[N1,N2,..,Nk,M1,M2,...Mj] となり, x の長さは size(y,1) に等しくなります.
87 </para>
88 <para>
89 <literal>method</literal> パラメータは補間の計算手法を設定します.
90 </para>
91 <variablelist>
92 <varlistentry>
93 <term>"linear"</term>
94 <listitem>
95 <para>
96 補間は線形手法により定義されます ( <link linkend="interpln">interpln</link>参照)
97 </para>
98 </listitem>
99 </varlistentry>
100 <varlistentry>
101 <term>"spline"</term>
102 <listitem>
103 <para>補間は3次スプライン補間により定義されます (
104 <link linkend="splin">splin</link> , <link linkend="interp">interp</link>参照)
105 </para>
106 </listitem>
107 </varlistentry>
108 <varlistentry>
109 <term>"nearest"</term>
110 <listitem>
111 <para>各 xp(j), yp(j) の値は,
112 xp(j)の最近傍にある x(i) に対応する y(i) の値を選択します.
113 </para>
114 </listitem>
115 </varlistentry>
116 </variablelist>
117 <para>
118 <literal>extrapolation</literal> パラメータは捕外用,すなわち
119 <literal>xp(i) </literal> が領域[x1,xn]内にない時
120 の計算手法を設定します.
121 </para>
122 <variablelist>
123 <varlistentry>
124 <term>"extrap"</term>
125 <listitem>
126 <para>捕外が定義された手法により実行されます.
127 yp=interp1(x,y,xp,method,"extrap")
128 </para>
129 </listitem>
130 </varlistentry>
131 <varlistentry>
132 <term>実数値</term>
133 <listitem>
134 <para>捕外用の実数値を選択できます.
135 この場合, 区間 [x1,xn] の中にない xp(i) について yp(i) は
136 この値とします.
137 yi=interp1(x,y,xp,method, 0)
138 </para>
139 </listitem>
140 </varlistentry>
141 <varlistentry>
142 <term>デフォルトの動作</term>
143 <listitem>
144 <para>捕外は(スプライン手法の場合)定義された手法により行われます
145 線形または最近傍手法の場合は nan となります.
146 yp=interp1(x,y,xp,method)
147 </para>
148 </listitem>
149 </varlistentry>
150 </variablelist>
151 </refsection>
152 <refsection>
153 <title>例</title>
154 <programlisting role="example"><![CDATA[
155x = linspace(0, 3, 20);
156y = x^2;
157xx = linspace(0, 3, 100);
158yy1 = interp1(x, y, xx, 'linear');
159yy2 = interp1(x, y, xx, 'spline');
160yy3 = interp1(x, y, xx, 'nearest');
161plot(xx, [yy1; yy2; yy3], x, y, '*')
162xtitle('interpolation of square function')
163legend(['linear','spline','nearest'], "in_upper_left");
164 ]]></programlisting>
165 <scilab:image>
166 x=linspace(0,3,20);
167 y=x.^2;
168 xx=linspace(0,3,100);
169 yy1=interp1(x,y,xx,'linear');
170 yy2=interp1(x,y,xx,'spline');
171 yy3=interp1(x,y,xx,'nearest');
172 plot(xx',[yy1;yy2;yy3]',x,y,'*')
173 xtitle('interpolation of square function')
174 legend(['linear','spline','nearest'],a=2)
175 </scilab:image>
176 </refsection>
177 <refsection role="see also">
178 <title>参照</title>
179 <simplelist type="inline">
180 <member>
181 <link linkend="interp">interp</link>
182 </member>
183 <member>
184 <link linkend="interpln">interpln</link>
185 </member>
186 <member>
187 <link linkend="splin">splin</link>
188 </member>
189 </simplelist>
190 </refsection>
191 <refsection>
192 <title>履歴</title>
193 <revhistory>
194 <revision>
195 <revnumber>5.4.0</revnumber>
196 <revremark>以前では, 入力引数の虚部は暗黙のうちに無視されていました.</revremark>
197 </revision>
198 </revhistory>
199 </refsection>
200</refentry>
diff --git a/scilab/modules/interpolation/help/pt_BR/interp1.xml b/scilab/modules/interpolation/help/pt_BR/interp1.xml
deleted file mode 100644
index 2bae0d1..0000000
--- a/scilab/modules/interpolation/help/pt_BR/interp1.xml
+++ /dev/null
@@ -1,182 +0,0 @@
1<?xml version="1.0" encoding="UTF-8"?>
2<!--
3 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
4 * Copyright (C) 2008 - INRIA - Farid BELAHCENE
5 *
6 * Copyright (C) 2012 - 2016 - Scilab Enterprises
7 *
8 * This file is hereby licensed under the terms of the GNU GPL v2.0,
9 * pursuant to article 5.3.4 of the CeCILL v.2.1.
10 * This file was originally licensed under the terms of the CeCILL v2.1,
11 * and continues to be available under such terms.
12 * For more information, see the COPYING file which you should have received
13 * along with this program.
14 *
15 -->
16<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink"
17 xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml"
18 xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook"
19 xmlns:scilab="http://www.scilab.org" xml:id="interp1" xml:lang="pt">
20 <refnamediv>
21 <refname>interp1</refname>
22 <refpurpose>função de interpolação 1d</refpurpose>
23 </refnamediv>
24 <refsynopsisdiv>
25 <title>Seqüência de Chamamento</title>
26 <synopsis>
27 yp = interp1(x, y, xp)
28 yp = interp1(x, y, xp, method)
29 yp = interp1(x, y, xp, method, extrapolation)
30 </synopsis>
31 </refsynopsisdiv>
32 <refsection>
33 <title>Parâmetros</title>
34 <variablelist>
35 <varlistentry>
36 <term>xp</term>
37 <listitem>
38 <para>escalar real, vetor ou matriz (ou hipermatriz) de reais
39 </para>
40 </listitem>
41 </varlistentry>
42 <varlistentry>
43 <term>x</term>
44 <listitem>
45 <para>vetor de reais</para>
46 </listitem>
47 </varlistentry>
48 <varlistentry>
49 <term>method</term>
50 <listitem>
51 <para>(opcional) string definindo o método de interpolação </para>
52 </listitem>
53 </varlistentry>
54 <varlistentry>
55 <term>extrapolation</term>
56 <listitem>
57 <para>(opcional) string ou valor real defindo os componentes yp(j)
58 para os valores xp(j) fora do intervalo [x1,xn].
59 </para>
60 </listitem>
61 </varlistentry>
62 <varlistentry>
63 <term>yp</term>
64 <listitem>
65 <para>vetor ou matriz (ou hipermatriz) </para>
66 </listitem>
67 </varlistentry>
68 </variablelist>
69 </refsection>
70 <refsection>
71 <title>Descrição</title>
72 <para>
73 Dados <literal>(x,y,xp)</literal>, esta função faz corresponder os
74 yp componentes aos xp por interpolação (linear por padrão) definida por x
75 e y.
76 </para>
77 <para>
78 Se <literal>yp</literal> é um vetor, então o comprimento de xp é
79 igual ao comprimento de <literal>yp</literal>. Se <literal>yp</literal> é
80 uma matriz, então <literal>xp</literal> tem o mesmo comprimento que cada
81 uma das colunas de yp. Se <literal>yp</literal> é uma hipermatriz, então o
82 comprimento de <literal>xp</literal> é o mesmo da primeira dimensão de
83 <literal>yp. </literal>
84 </para>
85 <para>Se size(y)=[C,M1,M2,M3,...,Mj] e size(xp)= [N1,N2,N3,...,Nk] então
86 size(yp)= [N1,N2,..,Nk,M1,M2,...Mj] e o comprimento de x deve ser igual a
87 size(y,1)
88 </para>
89 <para>
90 O parâmetro <literal>method</literal> ajusta a regra avaliação para
91 interpolação.
92 </para>
93 <variablelist>
94 <varlistentry>
95 <term>"linear"</term>
96 <listitem>
97 <para>
98 a interpolação é definida pelo método linear (ver <link linkend="interpln">interpln</link>)
99 </para>
100 </listitem>
101 </varlistentry>
102 <varlistentry>
103 <term>"spline"</term>
104 <listitem>
105 <para>
106 definição de interpolação por spline cúbico (ver <link linkend="splin">splin</link> , <link linkend="interp">interp</link>)
107 </para>
108 </listitem>
109 </varlistentry>
110 <varlistentry>
111 <term>"nearest"</term>
112 <listitem>
113 <para>para cada valor xp(j), yp(j) toma o valor ou y(i)
114 correspondente ao x(i), o vízinho mais próximo de xp(j)
115 </para>
116 <para/>
117 </listitem>
118 </varlistentry>
119 </variablelist>
120 <para>
121 O parâmetro <literal>extrapolation</literal> ajusta a regra de
122 avaliação para extrapolação, i.e para <literal>xp(i) </literal>fora do
123 intervalo [x1,xn]
124 </para>
125 <variablelist>
126 <varlistentry>
127 <term>"extrap"</term>
128 <listitem>
129 <para>a extrapolação é realizada pelo método definido.
130 yp=interp1(x,y,xp,method,"extrap")
131 </para>
132 </listitem>
133 </varlistentry>
134 <varlistentry>
135 <term>valor real</term>
136 <listitem>
137 <para>você pode escolher um valor real para extrapolação. Deste
138 modo, yp(i) toma este valor para xp(i) fora do intervalo [x1,xn],
139 por exemplo 0 (mas também nan ou inf). yi=interp1(x,y,xp,method, 0)
140 </para>
141 </listitem>
142 </varlistentry>
143 <varlistentry>
144 <term>por padrão</term>
145 <listitem>
146 <para>a extrapolação é realizada pelo método definido (para método
147 spline), e por nan para os métodos "linear" e "nearest".
148 yp=interp1(x,y,xp,method)
149 </para>
150 </listitem>
151 </varlistentry>
152 </variablelist>
153 </refsection>
154 <refsection>
155 <title>Exemplos</title>
156 <programlisting role="example"><![CDATA[
157x = linspace(0, 3, 20);
158y = x^2;
159xx = linspace(0, 3, 100);
160yy1 = interp1(x, y, xx, 'linear');
161yy2 = interp1(x, y, xx, 'spline');
162yy3 = interp1(x, y, xx, 'nearest');
163plot(xx, [yy1; yy2; yy3], x, y, '*')
164xtitle('interpolation of square function')
165legend(['linear','spline','nearest'], "in_upper_left");
166 ]]></programlisting>
167 </refsection>
168 <refsection role="see also">
169 <title>Ver Também</title>
170 <simplelist type="inline">
171 <member>
172 <link linkend="interp">interp</link>
173 </member>
174 <member>
175 <link linkend="interpln">interpln</link>
176 </member>
177 <member>
178 <link linkend="splin">splin</link>
179 </member>
180 </simplelist>
181 </refsection>
182</refentry>
diff --git a/scilab/modules/interpolation/help/ru_RU/interp1.xml b/scilab/modules/interpolation/help/ru_RU/interp1.xml
new file mode 100644
index 0000000..d0c573f
--- /dev/null
+++ b/scilab/modules/interpolation/help/ru_RU/interp1.xml
@@ -0,0 +1,356 @@
1<?xml version="1.0" encoding="UTF-8"?>
2<!--
3 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
4 * Copyright (C) 2008 - INRIA - Farid BELAHCENE
5 * Copyright (C) 2012 - 2016 - Scilab Enterprises
6 * Copyright (C) 2021 - Samuel GOUGEON - Le Mans Université
7 *
8 * This file is hereby licensed under the terms of the GNU GPL v2.0,
9 * pursuant to article 5.3.4 of the CeCILL v.2.1.
10 * This file was originally licensed under the terms of the CeCILL v2.1,
11 * and continues to be available under such terms.
12 * For more information, see the COPYING file which you should have received
13 * along with this program.
14 *
15 -->
16<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink"
17 xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml"
18 xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook"
19 xmlns:scilab="http://www.scilab.org" xml:id="interp1" xml:lang="ru">
20 <refnamediv>
21 <refname>interp1</refname>
22 <refpurpose>одномерная интерполяция в режиме nearest, linear или spline</refpurpose>
23 </refnamediv>
24 <refsynopsisdiv>
25 <title>Синтаксис</title>
26 <synopsis>
27 yp = interp1(y, xp)
28 yp = interp1(x, y, xp)
29 yp = interp1(.., xp, method)
30 yp = interp1(.., xp, method, extrapolation)
31 </synopsis>
32 </refsynopsisdiv>
33 <refsection>
34 <title>Аргументы</title>
35 <variablelist>
36 <varlistentry>
37 <term>x</term>
38 <listitem>
39 вектор по меньшей мере из двух вещественных чисел: ось абсцисс известных узлов
40 интерполяции, без повторов. По умолчанию
41 <itemizedlist>
42 <listitem>
43 если <varname>y</varname> является вектором: <literal>x=1:length(y)</literal>.
44 </listitem>
45 <listitem>
46 если <varname>y</varname> является матрицей или гиперматрицей:
47 <literal>x=1:size(y,1)</literal>.
48 </listitem>
49 </itemizedlist>
50 <para/>
51 </listitem>
52 </varlistentry>
53 <varlistentry>
54 <term>y</term>
55 <listitem>
56 вектор, матрица или гиперматрица вещественных или комплексных
57 чисел: значения в узлах интерполяции в соответствующих значениях
58 оси абсцисс <varname>x</varname>.
59 <itemizedlist>
60 <listitem>
61 если <varname>y</varname> является вектором, то <varname>x</varname>
62 и <varname>y</varname> должны быть одной длины.
63 </listitem>
64 <listitem>
65 если <varname>y</varname> является матрицей или гиперматрицей, то мы
66 должны иметь <literal>length(x)==size(y,1)</literal>. Каждый столбец
67 <varname>y</varname> тогда интерполирован в зависимости от тех же
68 значений оси абсцисс <varname>x</varname>, для указанного
69 <varname>xp</varname>.
70 </listitem>
71 </itemizedlist>
72 <para/>
73 </listitem>
74 </varlistentry>
75 <varlistentry>
76 <term>xp</term>
77 <listitem>
78 скаляр, вектор, матрица или гиперматрица десятичных чисел:
79 ось абсцисс точек, чьи значения <varname>yp</varname> должны быть
80 вычислены в соответствии с данными узлов интерполяции.
81 <para/>
82 </listitem>
83 </varlistentry>
84 <varlistentry>
85 <term>yp</term>
86 <listitem>
87 вектор, матрица или гиперматрица чисел: интерполированные значения
88 <varname>y</varname> в указанном <varname>xp</varname>.
89 <itemizedlist>
90 <listitem>
91 если <varname>y</varname> является вектором:
92 <varname>yp</varname> имеет размер <varname>xp</varname>.
93 <para/>
94 </listitem>
95 <listitem>
96 если <varname>y</varname> является матрицей или гиперматрицей:
97 <itemizedlist>
98 <listitem>
99 если <varname>xp</varname> является скаляром или вектором:
100 <literal>size(yp)</literal> равен
101 <literal>[length(xp) size(y)(2:$)]</literal>
102 </listitem>
103 <listitem>
104 если <varname>xp</varname> является матрицей или
105 гиперматрицей: <literal>size(yp)</literal> равен
106 <literal>[size(xp) size(y)(2:$)]</literal>
107 </listitem>
108 </itemizedlist>
109 </listitem>
110 </itemizedlist>
111 <para/>
112 </listitem>
113 </varlistentry>
114 <varlistentry>
115 <term>method</term>
116 <listitem>
117 строковое значение, определяющее метод интерполяции. Возможными
118 значениями и обработкой являются:
119 <para/>
120 <table>
121 <tr>
122 <td><emphasis role="bold">"linear"</emphasis>:</td>
123 <td>
124 линейная интерполяция между последовательными узлами, используется
125 по умолчанию.
126 </td>
127 </tr>
128 <tr>
129 <td><emphasis role="bold">"spline"</emphasis>:</td>
130 <td>
131 интерполяция кубическими сплайнами
132 </td>
133 </tr>
134 <tr>
135 <td><emphasis role="bold">"nearest"</emphasis>:</td>
136 <td>
137 <para>
138 для каждого значения <literal>xp(j)</literal> <literal>yp(j)</literal>
139 принимает значение <literal>y(i)</literal>, соответствующее
140 <literal>x(i)</literal>, ближайшего соседа
141 <literal>xp(j)</literal>
142 </para>
143 </td>
144 </tr>
145 </table>
146 <para/>
147 </listitem>
148 </varlistentry>
149 <varlistentry>
150 <term>extrapolation</term>
151 <listitem>
152 строковое значение или число, определяющее элементы <literal>yp(j)</literal>
153 для значений <literal>xp(j)</literal> за пределами интервала
154 <literal>[x(1)=min(x),x($)=max(x)]</literal>. Мы полагаем здесь и далее,
155 что <varname>x</varname> и <varname>y</varname> уже соответственно
156 отсортированы.
157 <para/>
158 <table>
159 <tr>
160 <td><emphasis role="bold">"extrap"</emphasis>:</td>
161 <td>
162 <literal>interp1(x,y,xp, method, "extrap")</literal> эквивалентен
163 <literal>interp1(x,y,xp, method, method)</literal>.
164 </td>
165 </tr>
166 <tr>
167 <td><emphasis role="bold">"linear"</emphasis>:</td>
168 <td>
169 Может использоваться с методами интерполяции
170 <literal>"spline"</literal> (и, очевидно,
171 <literal>"linear"</literal>).
172 </td>
173 </tr>
174 <tr>
175 <td><emphasis role="bold">"periodic"</emphasis>:</td>
176 <td>
177 Этот тип экстраполяции может использоваться с методами
178 интерполяции <literal>"linear"</literal> или
179 <literal>"spline"</literal>. Тогда: если <varname>y</varname>
180 является вектором, то требуется, чтобы
181 <literal>y(1)==y($)</literal>; в противном случае требуется,
182 чтобы <literal>y(1,:)==y($,:)</literal>.
183 </td>
184 </tr>
185 <tr>
186 <td><emphasis role="bold">"edgevalue"</emphasis>:</td>
187 <td>
188 Тогда <literal>yp(i)=y(1)</literal> для каждого
189 <literal>xp(i)&lt;x(1)</literal>, и
190 <literal>yp(i)=y($)</literal> для каждого
191 <literal>xp(i)>x($)</literal>.
192 </td>
193 </tr>
194 <tr>
195 <td><emphasis role="bold">padding</emphasis>:</td>
196 <td>
197 <varname>padding</varname> является десятичным или
198 комплексным числом, используемым для установки
199 <literal>yp(i)=padding</literal> для каждого
200 <literal>xp(i) ∉ [min(x),max(x)]</literal>. Например:
201 <literal>yi=interp1(x,y,xp,method, 0)</literal>.
202 </td>
203 </tr>
204 <tr>
205 <td><emphasis role="bold">(none)</emphasis>:</td>
206 <td>
207 По умолчанию экстраполяция выполняется сплайнами, когда
208 сплайны используются для интерполяции и дополнением
209 значениями <literal>%nan</literal>, если интерполяция
210 линейна или по "ближайшему" узлу.
211 </td>
212 </tr>
213 </table>
214 </listitem>
215 </varlistentry>
216 </variablelist>
217 </refsection>
218 <refsection>
219 <title>Описание</title>
220 <para>
221 Указывая <literal>(x,y,xp)</literal>, данная функция вычисляет элементы
222 <varname>yp</varname>, соответствующие <varname>xp</varname> с помощью интерполяции
223 между известными данными, указанными в узлах <literal>(x,y)</literal>.
224 </para>
225 <para>
226 <varname>x</varname> предварительно отсортирован в порядке возрастания, а значения
227 <varname>y</varname> либо по столбцам тогда сортируются соответственно.
228 </para>
229 <para>
230 <emphasis role="bold">Интерполяция комплексных значений</emphasis>:
231 Если <varname>y</varname> является комплексным, то его вещественная и мнимая
232 части интерполируются отдельно, и потом суммируются для формирования
233 комплексного <varname>yp</varname>.
234 </para>
235 <para>
236 <emphasis role="bold">interp1(x,y,xp,"nearest")</emphasis>:
237 Для любого <varname>xp</varname> в середине интервала
238 <literal>[x(i),x(i+1)]</literal> верхняя граница <literal>x(i+1)</literal>
239 рассматривается как ближайшее значение <varname>x</varname>, и присваивается
240 <literal>yp=y(i+1)</literal>.
241 </para>
242 <refsect3>
243 <title>линейные интерполяции</title>
244 Они выполняются через функцию <literal>linear_interpn(..)</literal>,
245 с соответствующей опцией интерполяции <literal>"edgevalue"→"C0"</literal>,
246 <literal>"linear"→"natural"</literal>, <literal>"periodic"→"periodic"</literal>.
247 </refsect3>
248 <refsect3>
249 <title>интерполяции сплайнами</title>
250 <para>
251 <emphasis role="bold">interp1(..,xp,"spline")</emphasis> или
252 <emphasis role="bold">interp1(..,xp,"spline","spline")</emphasis> или
253 <emphasis role="bold">interp1(..,xp,"spline","extrap")</emphasis>
254 используют условия граней <literal>не_узел</literal>.
255 Экстраполяция выполняется с помощью обоих сплайновых полиномов,
256 вычисленных на гранях <literal>(x,y)</literal>.
257 </para>
258 <para>
259 <emphasis role="bold">interp1(..,xp,"spline","edgevalue")</emphasis>
260 использует условия граней <literal>не_узел</literal>, а затем вызывает
261 <literal>interp(..,"C0")</literal>, чтобы выполнить фактическую
262 интерполяцию и экстраполяцию.
263 </para>
264 <para>
265 <emphasis role="bold">interp1(..,xp,"spline","periodic")</emphasis>
266 вызывает оба <literal>splin(..)</literal>, а затем <literal>interp(..)</literal>
267 с их опцией <literal>"periodic"</literal>.
268 </para>
269 <para>
270 <emphasis role="bold">interp1(..,xp,"spline","linear")</emphasis>
271 вызывает <literal>splin(..,"natural")</literal> для условий линейных граней, а затем
272 передаётся в <literal>interp(..,"linear")</literal>.
273 </para>
274 </refsect3>
275 </refsection>
276 <refsection>
277 <title>Примеры</title>
278 <programlisting role="example"><![CDATA[
279x = linspace(0, 10, 11)';
280y = sin(x);
281xx = linspace(0,10,1000)';
282yy2 = interp1(x, y, xx, 'linear');
283yy1 = interp1(x, y, xx, 'nearest');
284yy3 = interp1(x, y, xx, 'spline');
285clf
286h = plot(xx, [yy1 yy2 yy3], x, y, '.')
287h(1).mark_size = 8;
288title "Interpolation of a poorly sampled sin() function" fontsize 3
289legend(['nearest','linear','spline','nodes'], "in_lower_left");
290 ]]></programlisting>
291 <scilab:image>
292 x = linspace(0, 10, 11)';
293 y = sin(x);
294 xx = linspace(0,10,1000)';
295 yy2 = interp1(x, y, xx, 'linear');
296 yy1 = interp1(x, y, xx, 'nearest');
297 yy3 = interp1(x, y, xx, 'spline');
298 clf
299 h=plot(xx, [yy1 yy2 yy3],x, y, '.')
300 h(1).mark_size = 8;
301 title "Interpolation of a poorly sampled sin() function" fontsize 3
302 legend(['nearest','linear','spline','nodes'], "in_lower_left");
303 gcf().axes_size = [600 400];
304 </scilab:image>
305 </refsection>
306 <refsection role="see also">
307 <title>Смотрите также</title>
308 <simplelist type="inline">
309 <member>
310 <link linkend="interp">interp</link>
311 </member>
312 <member>
313 <link linkend="splin">splin</link>
314 </member>
315 <member>
316 <link linkend="linear_interpn">linear_interpn</link>
317 </member>
318 </simplelist>
319 </refsection>
320 <refsection>
321 <title>История</title>
322 <revhistory>
323 <revision>
324 <revnumber>6.1.1</revnumber>
325 <revremark>
326 <itemizedlist>
327 <listitem>
328 Для комплексных значений <varname>y</varname>,
329 <literal>imag(y)</literal> более не игнорируется:
330 <literal>real(y)</literal> и <literal>imag(y)</literal>
331 части интерполируются раздельно.
332 </listitem>
333 <listitem>
334 Экстраполяция <literal>"periodic"</literal> добавлена для
335 линейной и сплайновой интерполяций.
336 </listitem>
337 <listitem>
338 Экстраполяция <literal>"edgevalue"</literal> добавлена для
339 всех интерполяций ближайшего соседа, линейной и сплайновой.
340 </listitem>
341 <listitem>
342 Экстраполяция <literal>"linear"</literal> добавлена для
343 сплайновой интерполяции.
344 </listitem>
345 <listitem>
346 Когда <varname>xp</varname> является гиперматрицей и
347 <literal>size(xp,1)==1</literal>, <literal>size(yp)</literal>
348 теперь всегда <literal>[size(xp) size(y)(2,$)</literal> вместо
349 <literal>[size(xp,2:$), size(y)(2,$)</literal>.
350 </listitem>
351 </itemizedlist>
352 </revremark>
353 </revision>
354 </revhistory>
355 </refsection>
356</refentry>
diff --git a/scilab/modules/interpolation/macros/interp1.sci b/scilab/modules/interpolation/macros/interp1.sci
index 56dac19..1709911 100644
--- a/scilab/modules/interpolation/macros/interp1.sci
+++ b/scilab/modules/interpolation/macros/interp1.sci
@@ -1,7 +1,7 @@
1// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab 1// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
2// Copyright (C) INRIA - Farid BELAHCENE 2// Copyright (C) INRIA - Farid BELAHCENE
3//
4// Copyright (C) 2012 - 2016 - Scilab Enterprises 3// Copyright (C) 2012 - 2016 - Scilab Enterprises
4// Copyright (C) 2021 - Samuel GOUGEON - Le Mans Université
5// 5//
6// This file is hereby licensed under the terms of the GNU GPL v2.0, 6// This file is hereby licensed under the terms of the GNU GPL v2.0,
7// pursuant to article 5.3.4 of the CeCILL v.2.1. 7// pursuant to article 5.3.4 of the CeCILL v.2.1.
@@ -10,433 +10,373 @@
10// For more information, see the COPYING file which you should have received 10// For more information, see the COPYING file which you should have received
11// along with this program. 11// along with this program.
12 12
13function yi=interp1(varargin) 13function yi = interp1(varargin)
14 14 // yi = interp1(x, y, xi)
15 // yi=interp1(x,y,xi[,method[,interpolation) 15 // yi = interp1(x, y, xi, method)
16 // This function performs the yi values corresponding to xi by interpolation defined by x and y. 16 // yi = interp1(x, y, xi, method, extrapolation)
17 // Inputs : 17 // yi = interp1(y, xi)
18 // x , y : given data, x is a reals vector, y is a vector, matrix, or hypermatrix of reals 18 // yi = interp1(y, xi, method)
19 // if y is a vector, the length of x must be equal to the length of y, 19 // yi = interp1(y, xi, method, extrapolation)
20 // else the size of the first dimension of y must be equal to length of x.
21 // xi : a vector, matrix, or hypermatrix of reals
22 // Output
23 // yi : reals vector, matrix or hypermatrix, the values corresponding to xi by interpolation defined by x and y
24 // if size(y)=[C,N1,N2,N3,....] and size(xi)=[M1,M2,M3,M4] then size(xi)=[M1,M2,M3,M4,N1,N2,N3,N4,..], and length of x must be equal to C.
25 // Several kind of intepolations may be computed by selecting the appropriate method parameter:
26 // The methods are:
27 // linear : this is the default method (using the interp Scilab function)
28 // spline : this is the cubic spline interpolation (using interpln and splin Scilab functions)
29 // nearest : yi take the values corresponding to the nearest neighbor of xi
30 // 20 //
31 // Several kind of extrapolations may be computed : 21 // x: real vector
32 // 'extrap' : the extrapolation points is performed by the defined method 22 // y: real or complex vector of length(x) or array of size(y,1)==length(x)
33 // real value : you can choose a real value for extrapolation, in this way yp(i) takes this value for xp(i) not in [x1,xn] interval, for example 0 (but also nan or inf). 23 // xi: real scalar, vector, matrix or hypermatrix of any size.
34 // by default the extrapolation is performed by the defined method (for spline method), and by nan for linear and nearest method. 24 // method = "linear" (default) | "nearest" | "spline"
35 // F.B 25 // extrapolation = "extrap" (= method) | real or complex padding scalar |
26 // * if method = "linear": "linear" | "periodic" | "edgevalue"
27 // * if method = "spline": "linear" | "periodic" | "edgevalue" | "spline"
28 // * if method = "nearest": "nearest" = "edgevalue"
36 29
37 rhs=size(varargin) 30 // ========================
31 // CHECKING INPUT ARGUMENTS
32 // ========================
33 rhs = size(varargin)
38 // 2 <= Number of input arguments <= 5 34 // 2 <= Number of input arguments <= 5
39 if rhs<2 | rhs>5 then 35 if rhs<2 | rhs>5 then
40 error(msprintf(gettext("%s: Wrong number of input arguments: Must be between %d and %d.\n"),"interp1",2,5)); 36 msg = _("%s: Wrong number of input arguments: Must be between %d and %d.\n")
37 error(msprintf(msg, "interp1", 2, 5));
41 end 38 end
42 39
43 //if yi=interp1(x,y,xi,..) not change 40 // Checking x
44 //if yi=interp1(y,xi,...) replace input argument by yi=interp1(x,y,xi,..),whith x=1:size(y,1) by default 41 // ----------
45 if rhs==2 | (rhs>2 & type(varargin(3))==10) then 42 iy = 2 // index of y argin, for error messages
46 if isvector(varargin(1)) then 43 if rhs >= 3 & type(varargin(1))<>0 & type(varargin(3))<>10 then
47 X=1:size(varargin(1),"*") 44 // x is explicit: interp1(x, y, xi,..)
48 elseif size(size(varargin(1)),"*")==2 then 45 x = varargin(1)
49 if (size(varargin(1),1)>1 & size(varargin(1),2)>1) then 46 if type(x)<>1 | ~isreal(x) then
50 X=1:size(varargin(1),1) 47 msg = _("%s: Argument #%d: Real numbers expected.\n")
51 end 48 error(msprintf(msg,"interp1",1))
49 elseif ~isvector(x)
50 msg = _("%s: Argument #%d: Vector expected.\n")
51 error(msprintf(msg, "interp1", 1))
52 elseif or(isnan(x)) then
53 msg = _("%s: Argument #%d: Nan value forbidden.\n")
54 error(msprintf(msg, "interp1", 1))
55 end
56 x = matrix(x,1,-1)
57 y = varargin(2)
58 else // Default x:
59 if type(varargin(1))<>0 // interp1(y, xi,..)
60 iy = 1
61 y = varargin(1)
62 else // interp1(, y, xi,..)
63 y = varargin(2)
64 end
65 if isvector(y) then
66 x = 1:size(y,"*")
67 y = y(:)
52 else 68 else
53 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector or matrix expected.\n"),"interp1",2)); 69 x = 1:size(y,1)
54 end 70 end
55 for i=rhs:-1:1 71 end
56 varargin(i+1)=varargin(i) 72 if length(x) <= 1
73 if iy==2
74 msg = _("%s: Argument #%d: At least 2 elements expected.\n")
75 error(msprintf(msg, "interp1", 1))
76 else
77 msg = _("%s: Argument #%d: vector or at least 2 rows expected.\n")
78 error(msprintf(msg, "interp1", 2))
57 end 79 end
58 varargin(1)=X
59 end 80 end
81 // Sorting x in increasing order:
82 [?, k] = gsort(x,"g","i")
83 x = x(k) // should anyway be a row
60 84
61 //............................ 85 // Checking y
62 //ininialisation of xi 86 // ----------
63 //............................ 87 if type(y) <> 1 then
64 //xi components must be reals 88 msg = _("%s: Argument #%d: Decimal or complex numbers expected.\n")
65 xi=varargin(3) 89 error(msprintf(msg, "interp1", iy));
66 if type(xi)<>1 then
67 error(msprintf(gettext("%s: Wrong type for input argument #%d: Array of reals expected.\n"),"interp1",3));
68 end
69 //delete the dimension of xi equal to one after the second dimension
70 //or the first dimension
71 xisize=size(xi);
72 while size(xisize,"*")>=2 & xisize($)==1
73 xisize=xisize(1:$-1);
74 end 90 end
75 xisizetemp=xisize 91 if isvector(y) then
76 if size(xisize,"*")>=2 then 92 y = y(:)
77 if xisize(1)==1 then 93 if length(y) <> length(x)
78 xisize=xisize(2:$); 94 msg = _("%s: Arguments #%d and #%d: Same numbers of elements expected.\n")
95 error(msprintf(msg, "interp1", 1, 2))
96 end
97 else
98 if size(y,1) <> length(x) then
99 msg = _("%s: Arguments #%d and #%d: Same numbers of rows expected.\n")
100 error(msprintf(msg, "interp1", 1, 2))
79 end 101 end
80 end 102 end
103 // Sorting y according to x:
104 cmd = "y = y(k," + strcat(repmat(":",[1 ndims(y)-1]),",")+")"
105 execstr(cmd)
81 106
82 //------------------------- 107 // Checking xi
83 //Initialisation of x, y 108 // -----------
84 //------------------------- 109 if type(varargin(iy+1))<>1 | ~isreal(varargin(iy+1)) then
85 x=varargin(1); 110 msg = _("%s: Argument #%d: Decimal numbers expected.\n")
86 y=varargin(2); 111 error(msprintf(msg, "interp1", iy+1))
87 //x must be real vector
88 if type(x)<>1 then
89 error(msprintf(gettext("%s: Wrong type for input argument #%d: Array of reals expected.\n"),"interp1",1));
90 end 112 end
91 //y components must be reals 113 xi = varargin(iy+1)
92 if type(y)<>1 then 114 if xi == [] then
93 error(msprintf(gettext("%s: Wrong type for input argument #%d: Array of reals expected.\n"),"interp1",2)); 115 yi = []
116 return
94 end 117 end
95 //verification of x,y line/column 118 xisize = size(xi)
96 if isvector(x) then 119 if length(size(y))==2 & or(size(y)==1) then // y = 2D vector
97 if find(isnan(x))<>[] then 120 if iscolumn(xisize) then
98 error(msprintf(gettext("%s: Wrong value for input argument #%d: Reals expected but some NaN found.\n"),"interp1",1)); 121 xisize = length(xi)
99 end
100 if isvector(y) then
101 if size(x,"*")<>size(y,"*") then
102 error(msprintf(gettext("%s: Wrong size for input arguments #%d and #%d: Same size expected.\n"),"interp1",1,2));
103 end
104 elseif size(size(y),"*")>=2 then
105 if size(x,"*")<>size(y,1) then
106 error(msprintf(gettext("%s: Wrong size for input arguments #%d and #%d: Same size expected.\n"),"interp1",1,2));
107 end
108 else
109 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector or matrix expected.\n"),"interp1",2));
110 end 122 end
111 else 123 else
112 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector expected.\n"),"interp1",1)); 124 if length(size(xi))==2 & or(size(xi)==1) then // xi = scalar or 2D vector
113 end 125 xisize = length(xi)
114
115 // xi : increase order sorting (for xi)
116 [xtemp,p]=gsort(matrix(x,1,-1),"c","i")
117 x=matrix(xtemp,size(x))
118 x=matrix(x,1,-1)
119 if isvector(y) then
120 y=y(p)
121 elseif size(size(y),"*") then
122 for l=1:size(y,"*")/size(y,1)
123 y(:,l)=y(p,l)
124 end 126 end
125 else
126 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector or matrix expected.\n"),"interp1",2));
127 end
128
129 //-------------------------------------------------
130 // CASE : 3 inputs arguments : yi=interp1(x,y,xi)
131 //-------------------------------------------------
132
133 //default method : linear method is used
134 if size(varargin)==3 then
135 yi=interp1(x,y,xi,"linear",%nan)
136 end 127 end
137 128
138 //-------------------------------------------------- 129 // Checking method
139 // CASE : 4 inputs arguments : yi=interp1(x,y,xi,method) 130 // ---------------
140 //-------------------------------------------------- 131 if rhs > iy+1 & type(varargin(iy+2))<>0 then
141 132 method = varargin(iy+2)
142 if size(varargin)==4 then 133 if size(method,"*") <> 1
143 select part(varargin(4),1) 134 msg = _("%s: Argument #%d: Scalar (1 element) expected.\n")
144 //------------------------------------------- 135 error(msprintf(msg, "interp1", iy+2))
145 // Linear method : yi=linear(x,y,xi,'linear') 136 end
146 //------------------------------------------- 137 if type(method) <> 10 | ~or(convstr(method)==["linear","nearest","spline"])
147 // the values of extrapolation points are nan for linear method 138 msg = _("%s: Argument #%d: Must be in the set {%s}.\n")
148 case "l" 139 error(msprintf(msg, "interp1", iy+2, """linear"",""spline"",""nearest"""))
149
150 yi=interp1(x,y,xi,"linear",%nan)
151
152 //-------------------------------------------
153 // Spline method yi=interp1(x,y,xi,'spline')
154 //-------------------------------------------
155 // the extrapolation used the spline method
156 case "s"
157 if xi==[] then
158 yi=[]
159 return
160 end
161
162 yi=interp1(x,y,xi,"spline","extrap")
163 //----------------------------------------------
164 // Nearest method yi=interp1(x,y,xi,'nearest')
165 //----------------------------------------------
166 // the values of extrapolation points are nan for nearest method
167 case "n"
168 if xi==[] then
169 yi=[]
170 return
171 end
172 yi=interp1(x,y,xi,"nearest",%nan)
173 else
174 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or ''%s'' or ''%s'' expected.\n"),"interp1",4,"linear","nearest"));
175 end 140 end
141 method = convstr(method)
142 else
143 method = "linear"
176 end 144 end
177 145
178 //------------------------------------------------------------------------------------------------- 146 // Checking extrapolation method
179 // CASE : 5 inputs arguments : yi=interp1(x,y,xi,method,'extrap') or yi=interp1(x,y,xi,method,extrapval) 147 // -----------------------------
180 //------------------------------------------------------------------------------------------------- 148 extraVal = []
181 149 if rhs > iy+2 then
182 if size(varargin)==5 then 150 extra = varargin(iy+3)
183 select part(varargin(4),1) 151 err = %f
184 //----------------------------------------------------------------------------------- 152 if type(extra)==10 & size(extra,"*")==1
185 // Linear method : yi=linear(x,y,xi,'linear','extrap') or yi=interp1(x,y,xi,method,extrapval) 153 extra = convstr(extra)
186 //------------------------------------------------------------------------------------ 154 if extra=="extrap"
187 case "l" 155 extra = "natural"
188 xitemp=matrix(xi,-1,1) 156 elseif or(method==["linear" "spline"])
189 // y is a vector 157 if extra=="linear" then
190 if isvector(y) then 158 if method=="linear"
191 yi = resize_matrix(0, size(xitemp)) 159 extra = "natural"
192 [x,ind]=gsort(matrix(x,1,-1),"c","i") 160 end
193 if varargin(5)==%nan then 161 elseif extra=="edgevalue"
194 yi=linear_interpn(xitemp,x,y(ind),"by_nan"); 162 extra = "C0"
195 end 163 elseif extra=="periodic"
196 if type(varargin(5))==10 then 164 msg = []
197 if varargin(5)<>"extrap" then 165 if isvector(y) & y($)<>y(1)
198 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap")); 166 msg = _("%s: Argument #%d: periodicity y($)==y(1) expected.\n")
199 else 167 elseif ~isvector(y) & or(y($,:)<>y(1,:))
200 yi=linear_interpn(xitemp,x,y(ind),"natural"); 168 msg = _("%s: Argument #%d: periodicity y($,:)==y(1,:) expected.\n")
201 end 169 end
202 elseif type(varargin(5))==1 then 170 if msg <> []
203 yi=linear_interpn(xitemp,x,y(ind),"by_nan"); 171 error(msprintf(msg, "interp1",iy))
204 if ~isnan(varargin(5)) then
205 k=find(xitemp>max(x)|xitemp<min(x))
206 yi(k)=varargin(5)
207 end 172 end
173 elseif extra<>"natural"
174 err = %t
208 end 175 end
209 if size(xisize,"*")>=2 176 else // nearest
210 yi=matrix(yi,xisize) 177 if or(extra==["edgevalue" "nearest" "natural"])
178 extra = "natural"
211 else 179 else
212 yi=matrix(yi,xisizetemp) 180 err = %t
213 end
214
215 // y is matrix or hypermatrix
216 elseif size(size(y),"*")>=2 then
217 ysize=size(y)
218 ky=ysize(2:$)
219 yi = resize_matrix(0, [size(xitemp) ky])
220 [x,ind]=gsort(matrix(x,1,-1),"c","i")
221 //extrapolation
222 if type(varargin(5))==10 then
223 if varargin(5)<>"extrap" then
224 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap"));
225 else
226 if xitemp==[] then
227 yi=[]
228 return
229 end
230 for l=1:size(y,"*")/size(y,1)
231 ytemp=y(:,l)
232 yi(:,l)=matrix(linear_interpn(xitemp,x,ytemp(ind),"natural"),size(xitemp))
233 end
234 end
235 elseif type(varargin(5))==1 then
236 if xitemp==[] then
237 yi=[]
238 return
239 end
240 for l=1:size(y,"*")/size(y,1)
241 ytemp=y(:,l)
242 yi(:,l)=matrix(linear_interpn(xitemp,x,ytemp(ind),"by_nan"),size(xitemp))
243 end
244 if ~isnan(varargin(5)) then
245 k=find(xitemp>max(x)|xitemp<min(x))
246 yi(k,:)=varargin(5)
247 end
248 end 181 end
249 yi=matrix(yi,[xisize,ky]) 182 end
183 else
184 if type(extra)<>1 | length(extra)<>1
185 err = %t
250 else 186 else
251 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector or matrix expected.\n"),"interp1",2)); 187 extraVal = extra
188 extra = "natural"
252 end 189 end
190 end
191 if err
192 if method=="nearest"
193 tmp = """extrap"",""edgevalue"",""nearest"""
194 else
195 tmp = """extrap"",""linear"",""periodic"",""edgevalue"""
196 end
197 msg = _("%s: Argument #%d: scalar number or one of {%s} expected.\n")
198 error(msprintf(msg, "interp1", iy+3, tmp))
199 end
200 else
201 extra = "natural"
202 if method <> "spline"
203 extraVal = %nan
204 end
205 end
253 206
254 //------------------------------------------------------------------------------------- 207 // ==========
255 // Spline method yi=interp1(x,y,xi,'spline','extrap') or yi=interp1(x,y,xi,'spline',extrapval) 208 // PROCESSING
256 //------------------------------------------------------------------------------------- 209 // ==========
257 case "s" 210 if ~isreal(y) then
258 if xi==[] then 211 if extraVal==[]
259 if varargin(5)=="extrap"|type(varargin(5))==1 then 212 [extraR, extraI] = (extra, extra);
260 yi=[] 213 else
261 return 214 [extraR, extraI] = (real(extraVal), imag(extraVal));
262 else 215 end
263 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap")); 216 yir = interp1(x, real(y), xi, method, extraR)
217 yii = interp1(x, imag(y), xi, method, extraI)
218 yi = complex(yir, yii)
219 return
220 end
221 // REAL CASE
222 // =========
223 select part(method, 1)
224 case "l"
225 // LINEAR
226 // ------
227 xitemp = matrix(xi,-1,1)
228 // y is a vector
229 if isvector(y) then
230 yi = resize_matrix(0, size(xitemp))
231 if extraVal <> [] then
232 yi = linear_interpn(xitemp, x, y, "by_nan");
233 if ~isnan(extraVal) then
234 k = find(xitemp<min(x) | xitemp>max(x))
235 yi(k) = extraVal
264 end 236 end
237 else
238 yi = linear_interpn(xitemp, x, y, extra);
265 end 239 end
266 xitemp=matrix(xi,-1,1) 240 yi = matrix(yi, xisize)
267 //y is a vector 241
268 if isvector(y) then 242 else // y is matrix or hypermatrix
269 yi = resize_matrix(0, size(xitemp)) 243 ky = size(y)(2:$)
270 yi=interp(xitemp,matrix(x,1,-1),matrix(y,1,-1),splin(matrix(x,1,-1),matrix(y,1,-1)),"natural"); 244 ncolumns = prod(ky)
271 if type(varargin(5))==10 then 245 yi = resize_matrix(0, [size(xitemp) ky])
272 if varargin(5)<>"extrap" then 246 if extraVal <> [] then
273 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap")); 247 for ic = 1:ncolumns
274 end 248 tmp = linear_interpn(xitemp, x, y(:,ic), "by_nan")
275 elseif type(varargin(5))==1 then 249 yi(:,ic) = matrix(tmp, size(xitemp))
276 k=find(xitemp>max(x)|xitemp<min(x))
277 yi(k)=varargin(5)
278 end 250 end
279 if size(xisize,"*")>=2 251 if ~isnan(extraVal) then
280 yi=matrix(yi,xisize) 252 k = find(xitemp < min(x) | xitemp > max(x))
281 else 253 yi(k,:) = extraVal
282 yi=matrix(yi,xisizetemp)
283 end 254 end
284 //y is a matrix or a hypermatrix 255 else // extra=="linear"=="natural" | "periodic" | "edgevalue"=="C0"
285 elseif size(size(y),"*")>=2 then 256 for ic = 1:ncolumns
286 ky=size(y) 257 tmp = linear_interpn(xitemp, x, y(:,ic), extra)
287 ky=ky(2:$) 258 yi(:,ic) = matrix(tmp, size(xitemp))
288 yi = resize_matrix(0, [size(xitemp) ky])
289 for l=1:size(y,"*")/size(y,1)
290 yi(:,l)=matrix(interp(matrix(xi,-1,1),matrix(x,-1,1),y(:,l),splin(matrix(x,-1,1),y(:,l)),"natural"),size(xitemp))//les composante de yi
291 end 259 end
292 //extrapolation 260 end
293 if type(varargin(5))==10 then 261 yi = matrix(yi, [xisize, ky])
294 if varargin(5)<>"extrap" then 262 end
295 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap")); 263
296 end 264 case "s"
297 elseif type(varargin(5))==1 then 265 // SPLINE
298 k=find(xitemp>max(x)|xitemp<min(x)) 266 // ------
299 yi(k,:)=varargin(5) 267 xitemp = matrix(xi,-1,1)
268 // y is a vector
269 if isvector(y) then
270 yi = resize_matrix(0, size(xitemp))
271 yrow = matrix(y,1,-1)
272 if extra=="linear"
273 yi = interp(xitemp, x, yrow, splin(x, yrow, "natural"), extra)
274 elseif extra=="periodic"
275 yi = interp(xitemp, x, yrow, splin(x, yrow, extra), extra)
276 else // default not_a_knot
277 yi = interp(xitemp, x, yrow, splin(x, yrow), extra)
278 end
279 if extraVal <> [] then
280 k = find(xitemp<min(x) | xitemp>max(x))
281 yi(k) = extraVal
282 end
283 yi = matrix(yi, xisize)
284
285 else // y is a matrix or a hypermatrix
286 ky = size(y)(2:$)
287 ncolumns = prod(ky)
288 yi = resize_matrix(0, [size(xitemp) ky])
289 x = matrix(x,-1,1)
290 xi = matrix(xi,-1,1)
291 for ic = 1:ncolumns
292 tmp = y(:,ic)
293 tmp = interp(xi, x, tmp, splin(x, tmp), extra)
294 yi(:,ic) = matrix(tmp, size(xitemp))
295 end
296 // extrapolation
297 if extraVal <> [] then
298 k = find(xitemp<min(x) | xitemp>max(x))
299 yi(k,:) = extraVal
300 end
301 yi = matrix(yi,[xisize,ky])
302 end
303
304 case "n"
305 // NEAREST
306 // -------
307 // if all xi values are nan, retuns nan values for yi
308 if and(isnan(xi)) then
309 yi = xi
310 return
311 end
312
313 xitemp = matrix(xi,1,-1)
314 knan = isnan(xitemp)
315 [xitemp,p] = gsort(xitemp(~knan),"g","i")
316 // For each xi, we search for the k index of the x nearest to xi:
317 k = zeros(xitemp);
318 nx = length(x);
319 j = length(xitemp);
320 i = nx;
321 while j>=1 & i>=1
322 if xitemp(j) >= x(i) then
323 if i <> nx then
324 k(j) = i;
300 end 325 end
301 yi=matrix(yi,[xisize,ky]) 326 j = j - 1;
302 else 327 else
303 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector or matrix expected.\n"),"interp1",2)); 328 i = i - 1;
304 end 329 end
330 end
331 k(xitemp < x(1)) = 1;
332 k(xitemp >= x(nx)) = nx - 1;
333 i = find(xitemp >= matrix((x(k)+x(k+1))/2,size(k)));
334 if i~=[]
335 k(i) = k(i)+1;
336 end
305 337
306 //--------------------------------------------------------------------------------------- 338 // y is vector
307 // Nearest method yi=interp1(x,y,xi,'nearest','extrap') or yi=interp1(x,y,xi,'nearest',extrapval) 339 if isvector(y) then
308 //--------------------------------------------------------------------------------------- 340 yi = y(k)
309 case "n" 341 yi = matrix(yi,1,-1)
310 //if all xi values are nan, retuns nan values for yi 342 // extrapolation
311 if size(find(isnan(xi)),"*")==size(xi,"*") then 343 if extraVal <> [] then
312 if varargin(5)=="extrap"|type(varargin(5))==1 then 344 n = find(xitemp<min(x) | xitemp>max(x) | isnan(xitemp))
313 yi=xi 345 yi(n) = extraVal
314 return
315 else
316 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap"));
317 end
318 end 346 end
319 if xi==[] then 347 yitemp = yi
320 if varargin(5)=="extrap"|type(varargin(5))==1 then 348 yi(p) = yitemp
321 yi=[] 349 ytemp = yi
322 return 350 yi = matrix(xi,1,-1)
323 else 351 yi(knan) = %nan
324 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap")); 352 yi(~knan) = ytemp
325 end 353 yi = matrix(yi, xisize)
354
355 else // y is a matrix or a hypermatrix
356 ky = size(y)(2:$)
357 yi = resize_matrix(0, [size(xitemp,"*") ky])
358 ncolumns = prod(ky)
359 for ic = 1:ncolumns
360 yi(:,ic) = y(k,ic)
326 end 361 end
327 //y is vector 362 // extrapolation
328 if isvector(y) then 363 if extraVal <> [] then
329 xitemp=matrix(xi,1,-1) 364 n = find(xitemp < min(x) | xitemp > max(x))
330 knan=find(isnan(xitemp)) 365 for ic = 1:ncolumns
331 knotnan=find(~isnan(xitemp)) 366 yi(n,ic) = extraVal
332 [xitemp,p]=gsort(matrix(xitemp(knotnan),1,-1),"c","i")
333 yi=matrix(xi,1,-1)
334 k=zeros(xitemp)
335 x_size=size(x,"*")
336 j=size(xitemp,"*")
337 i=x_size
338 while j>=1 & i>=1
339 if xitemp(j)>=x(i) then
340 if i<>x_size then
341 k(j)=i
342 end
343 j=j-1
344 else
345 i=i-1
346 end
347 end
348 k(xitemp<x(1)) = 1;
349 k(xitemp>=x(x_size)) = x_size-1;
350 i = find(xitemp >= matrix((x(k)+x(k+1))/2,size(k)));
351 if i~=[]
352 k(i) = k(i)+1;
353 end
354 yi=y(k)
355 yi=matrix(yi,1,-1)
356 //extrapolation
357 if type(varargin(5))==10 then
358 if varargin(5)<>"extrap" then
359 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap"));
360 end
361 elseif type(varargin(5))==1 then
362 n=find(xitemp>max(x)|xitemp<min(x)|isnan(xitemp)|isnan(xitemp))
363 yi(n)=varargin(5)
364 else
365 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap"));
366 end
367 yitemp=yi
368 yi(p)=yitemp
369 ytemp=yi
370 yi=matrix(xi,1,-1)
371 yi(knan)=%nan
372 yi(knotnan)=ytemp
373 if size(xisize,"*")>=2
374 yi=matrix(yi,xisize)
375 else
376 yi=matrix(yi,xisizetemp)
377 end
378 //y is a matrix or a hypermatrix
379 elseif size(size(y),"*")>=2 then
380 xitemp=matrix(xi,1,-1)
381 knan=find(isnan(xitemp))
382 knotnan=find(~isnan(xitemp))
383 [xitemp,p]=gsort(xitemp(knotnan),"c","i")
384 ind=size(y)
385 ind=ind(2:$)
386 yi = resize_matrix(0, [size(xitemp,"*") ind])
387 k=zeros(xitemp)
388 x_size=size(x,"*")
389 j=size(xitemp,"*")
390 i=x_size
391 while j>=1 & i>=1
392 if xitemp(j)>=x(i) then
393 if i<>x_size then
394 k(j)=i
395 end
396 j=j-1
397 else
398 i=i-1
399 end
400 end 367 end
401 k(xitemp<x(1)) = 1;
402 k(xitemp>=x(x_size)) = x_size-1;
403 i = find(xitemp >= matrix((x(k)+x(k+1))/2,size(k)));
404 if i~=[]
405 k(i) = k(i)+1;
406 end
407 for l=1:size(y,"*")/size(y,1)
408 ytemp=matrix(y(:,l),1,-1)
409 yi(:,l) =ytemp(k)
410 end
411 //extrapolation
412 if type(varargin(5))==10 then
413 if varargin(5)<>"extrap" then
414 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap"));
415 end
416 elseif type(varargin(5))==1 then
417 n=find(xitemp>max(x)|xitemp<min(x))
418 for l=1:size(y,"*")/size(y,1)
419 yi(n,l)=varargin(5)
420 end
421 else
422 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap"));
423 end
424 yitemp=yi
425 for l=1:size(y,"*")/size(y,1)
426 yi(p,l)=yitemp(:,l)
427 end
428 yitemp=yi
429 yi = resize_matrix(0, [size(xi,"*") ind])
430 for l=1:size(y,"*")/size(y,1)
431 yi(knan,l)=%nan
432 yi(knotnan,l)=yitemp(:,l)
433 end
434 yi=matrix(yi,[xisize,ind])
435 else
436 error(msprintf(gettext("%s: Wrong size for input argument #%d: Vector or matrix expected.\n"),"interp1",2));
437 end 368 end
438 else 369 yitemp = yi
439 error(msprintf(gettext("%s: Wrong value for input argument #%d: ''%s'' or real expected.\n"),"interp1",5,"extrap")); 370 for ic = 1:ncolumns
371 yi(p,ic) = yitemp(:,ic)
372 end
373 yitemp = yi
374 yi = resize_matrix(0, [size(xi,"*") ky])
375 for ic = 1:ncolumns
376 yi(knan,ic) = %nan
377 yi(~knan,ic) = yitemp(:,ic)
378 end
379 yi = matrix(yi,[xisize,ky])
440 end 380 end
441 end 381 end
442endfunction 382endfunction
diff --git a/scilab/modules/interpolation/tests/unit_tests/interp1.tst b/scilab/modules/interpolation/tests/unit_tests/interp1.tst
new file mode 100644
index 0000000..3a62bc7
--- /dev/null
+++ b/scilab/modules/interpolation/tests/unit_tests/interp1.tst
@@ -0,0 +1,654 @@
1// =============================================================================
2// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3// Copyright (C) 2021 - Samuel GOUGEON - Le Mans Université
4//
5// This file is distributed under the same license as the Scilab package.
6// =============================================================================
7// <-- CLI SHELL MODE -->
8// <-- NO CHECK REF -->
9
10// Checking the size of the result
11// ===============================
12for m = ["linear", "spline", "nearest"]
13 for x = list(1:0.5:4, (1:0.5:4)')
14 // Y vector => size(result) == size(xi)
15 Xi = list(%pi, 1.2:0.2:2, (1.2:0.2:2)', 1+rand(2,5)*4, ..
16 1+rand(1,2,5)*4, 1+rand(2,1,5)*4, 1+rand(2,5,6)*4);
17 Y = list(5:11, (5:11)');
18 for y = Y
19 for xi = Xi
20 r = interp1(x, y, xi, m);
21 assert_checkequal(size(r), size(xi));
22 end
23 end
24
25 // y = matrix or hypermatrix, xi = scalar or vector
26 Y = list(rand(7,3), rand(7,1,3), rand(7,3,4));
27 Xi = list(%pi, 1.2:0.2:2, (1.2:0.2:2)');
28 for y = Y
29 for xi = Xi
30 r = interp1(x, y, xi, m);
31 assert_checkequal(size(r), [length(xi) size(y)(2:$)]);
32 end
33 end
34
35 // Y and xi are matrices or hypermatrices
36 Y = list(rand(7,3), rand(7,1,3), rand(7,3,4));
37 Xi = list(1+rand(2,5)*4, 1+rand(1,2,5)*4, 1+rand(2,1,5)*4, 1+rand(2,5,6)*4);
38 for y = Y
39 for xi = Xi
40 r = interp1(x, y, xi, m);
41 assert_checkequal(size(r), [size(xi) size(y)(2:$)]);
42 end
43 end
44 end
45end
46
47
48// =============
49// y is a vector
50// =============
51x = 1:4;
52// method = "linear" (default)
53// ---------------------------
54// extrapolation "by_nan" (default)
55for y = list(2:5, (2:5)')
56 assert_checkequal(interp1(x,y, 0), %nan);
57 assert_checkequal(interp1(x,y,1.5), 2.5);
58 assert_checkequal(interp1(x,y,[0.5:5]), [%nan 2.5 3.5 4.5 %nan]);
59 assert_checkequal(interp1(x,y,[0.5:5]'), [%nan 2.5 3.5 4.5 %nan]');
60 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5]), [%nan 2.5 3.5;4.5 %nan %nan]);
61 Ref = cat(3, [%nan 2.5 3.5], [4.5 %nan %nan]);
62 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5])), Ref);
63end
64// extrapolation "extrap" = "linear"
65for y = list(2:5, (2:5)')
66 assert_checkequal(interp1(x,y, 0, ,"extrap"), 1);
67 assert_checkequal(interp1(x,y,1.5, ,"extrap"), 2.5);
68 assert_checkequal(interp1(x,y,[0.5:5], ,"extrap"), [1.5 2.5 3.5 4.5 5.5]);
69 assert_checkequal(interp1(x,y,[0.5:5]', ,"extrap"), [1.5 2.5 3.5 4.5 5.5]');
70 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,"extrap"), ..
71 [1.5 2.5 3.5;4.5 5.5 6.5]);
72 Ref = cat(3, [1.5 2.5 3.5], [4.5 5.5 6.5]);
73 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,"extrap"), Ref);
74end
75// extrapolation by padding
76for y = list(2:5, (2:5)')
77 assert_checkequal(interp1(x,y, 0, ,%pi), %pi);
78 assert_checkequal(interp1(x,y,1.5, ,%pi), 2.5);
79 assert_checkequal(interp1(x,y,[0.5:5], ,%pi), [%pi 2.5 3.5 4.5 %pi]);
80 assert_checkequal(interp1(x,y,[0.5:5]', ,%pi), [%pi 2.5 3.5 4.5 %pi]');
81 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,%pi), ..
82 [%pi 2.5 3.5;4.5 %pi %pi]);
83 Ref = cat(3, [%pi 2.5 3.5], [4.5 %pi %pi]);
84 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,%pi), Ref);
85end
86// periodic extrapolation
87for y = list([2 3 4 2], [2 3 4 2]')
88 assert_checkequal(interp1(x,y, 0, ,"periodic"), 4);
89 assert_checkequal(interp1(x,y,1.5, ,"periodic"), 2.5);
90 assert_checkequal(interp1(x,y,[0.5:5], ,"periodic"), [3 2.5 3.5 3 2.5]);
91 assert_checkequal(interp1(x,y,[0.5:5]', ,"periodic"), [3 2.5 3.5 3 2.5]');
92 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,"periodic"), ..
93 [3 2.5 3.5 ; 3 2.5 3.5]);
94 Ref = cat(3, [3 2.5 3.5], [3 2.5 3.5]);
95 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,"periodic"), Ref);
96end
97// extrapolation by edgevalue
98e = "edgevalue";
99for y = list(2:5, (2:5)')
100 assert_checkequal(interp1(x,y, 0, , e), 2);
101 assert_checkequal(interp1(x,y,1.5, , e), 2.5);
102 assert_checkequal(interp1(x,y,[0.5:5], , e), [2 2.5 3.5 4.5 5]);
103 assert_checkequal(interp1(x,y,[0.5:5]', , e), [2 2.5 3.5 4.5 5]');
104 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], , e), [2 2.5 3.5 ; 4.5 5 5]);
105 Ref = cat(3, [2 2.5 3.5], [4.5 5 5]);
106 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), , e), Ref);
107end
108
109// method = "nearest"
110// ------------------
111m = "nearest";
112// extrapolation "by_nan" (default)
113for y = list(2:5, (2:5)')
114 assert_checkequal(interp1(x,y, [0 5.4], m), [%nan %nan]);
115 assert_checkequal(interp1(x,y,1.5, m), 3);
116 assert_checkequal(interp1(x,y,[0.5:5], m), [%nan 3 4 5 %nan]);
117 assert_checkequal(interp1(x,y,[0.5:5]', m), [%nan 3 4 5 %nan]');
118 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m), [%nan 3 4;5 %nan %nan]);
119 Ref = cat(3, [%nan 3 4], [5 %nan %nan]);
120 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m), Ref);
121end
122// extrapolation "extrap" = "nearest" = "edgevalue"
123for e = ["extrap" "edgevalue" "nearest"]
124 for y = list(2:5, (2:5)')
125 assert_checkequal(interp1(x,y, [0 5.4], m, e), [2 5]);
126 assert_checkequal(interp1(x,y,1.5, m, e), 3);
127 assert_checkequal(interp1(x,y,[0.5:5], m, e), [2 3 4 5 5]);
128 assert_checkequal(interp1(x,y,[0.5:5]', m, e), [2 3 4 5 5]');
129 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, e), ..
130 [2 3 4 ; 5 5 5]);
131 Ref = cat(3, [2 3 4], [5 5 5]);
132 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, e), Ref);
133 end
134end
135// extrapolation by padding
136for y = list(2:5, (2:5)')
137 assert_checkequal(interp1(x,y, [0 5.4], m, %pi), [%pi %pi]);
138 assert_checkequal(interp1(x,y, 1.5, m, %pi), 3);
139 assert_checkequal(interp1(x,y, 0.5:5, m, %pi), [%pi 3 4 5 %pi]);
140 assert_checkequal(interp1(x,y,[0.5:5]', m, %pi), [%pi 3 4 5 %pi]');
141 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, %pi), ..
142 [%pi 3 4 ; 5 %pi %pi]);
143 Ref = cat(3, [%pi 3 4], [5 %pi %pi]);
144 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, %pi), Ref);
145end
146
147// method = "spline"
148// -----------------
149// extrapolation default == "spline" (=="natural")
150m = "spline";
151for y = list(2:5, (2:5)')
152 assert_checkalmostequal(interp1(x,y, [0 5.4], m), [1 6.4], 5*%eps);
153 assert_checkequal(interp1(x,y,1.5, m), 2.5);
154 assert_checkalmostequal(interp1(x,y,[0.5:5], m), 1.5:5.5, %eps);
155 assert_checkalmostequal(interp1(x,y,[0.5:5]', m), (1.5:5.5)', %eps);
156 assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m), [1.5 2.5 3.5;4.5 5.5 6.5], 5*%eps);
157 Ref = cat(3, [1.5 2.5 3.5], [4.5 5.5 6.5]);
158 assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m), Ref, 5*%eps);
159end
160// extrapolation "extrap" == "spline"
161for y = list(2:5, (2:5)')
162 assert_checkalmostequal(interp1(x,y, [0 5.4], m, "extrap"), [1 6.4], 5*%eps);
163 assert_checkequal(interp1(x,y,1.5, m, "extrap"), 2.5);
164 assert_checkalmostequal(interp1(x,y,[0.5:5], m, "extrap"), 1.5:5.5, %eps);
165 assert_checkalmostequal(interp1(x,y,[0.5:5]', m, "extrap"), (1.5:5.5)', %eps);
166 Ref = [1.5 2.5 3.5 ; 4.5 5.5 6.5];
167 assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, "extrap"), Ref, 5*%eps);
168 Ref = cat(3, [1.5 2.5 3.5], [4.5 5.5 6.5]);
169 assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, "extrap"), Ref, 5*%eps);
170end
171// extrapolation by padding
172for y = list(2:5, (2:5)')
173 assert_checkequal(interp1(x,y, [0 5.4], m, %pi), [%pi %pi]);
174 assert_checkequal(interp1(x,y, 1.5, m, %pi), 2.5);
175 assert_checkequal(interp1(x,y, 0.5:5, m, %pi), [%pi 2.5 3.5 4.5 %pi]);
176 assert_checkequal(interp1(x,y,[0.5:5]', m, %pi), [%pi 2.5 3.5 4.5 %pi]');
177 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, %pi), ..
178 [%pi 2.5 3.5 ; 4.5 %pi %pi]);
179 Ref = cat(3, [%pi 2.5 3.5], [4.5 %pi %pi]);
180 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, %pi), Ref);
181end
182// periodic extrapolation
183for y = list([2 3 4 2], [2 3 4 2]')
184 assert_checkalmostequal(interp1(x,y, [0 5.4], m, "periodic"), [4 3.736], %eps);
185 assert_checkequal(interp1(x,y, 1.5, m, "periodic"), 2.125);
186 assert_checkalmostequal(interp1(x,y, 0.5:5, m, "periodic"), [3 2.125 3.875 3 2.125], %eps);
187 assert_checkalmostequal(interp1(x,y,[0.5:5]', m, "periodic"), [3 2.125 3.875 3 2.125]', %eps);
188 assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, "periodic"), ..
189 [3 2.125 3.875 ; 3 2.125 3.875], %eps);
190 Ref = cat(3, [3 2.125 3.875], [3 2.125 3.875]);
191 assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, "periodic"), Ref, %eps);
192end
193// extrapolation by edgevalue
194e = "edgevalue";
195for y = list(2:5, (2:5)')
196 assert_checkequal(interp1(x,y, [0 5.4], m, e), [2 5]);
197 assert_checkequal(interp1(x,y, 1.5, m, e), 2.5);
198 assert_checkequal(interp1(x,y, 0.5:5, m, e), [2 2.5 3.5 4.5 5]);
199 assert_checkequal(interp1(x,y,[0.5:5]', m, e), [2 2.5 3.5 4.5 5]');
200 assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, e), ..
201 [2 2.5 3.5 ; 4.5 5 5]);
202 Ref = cat(3, [2 2.5 3.5], [4.5 5 5]);
203 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, e), Ref);
204end
205
206// =============
207// y is a matrix
208// =============
209x = 1:4;
210y = [0 2 4 6 ; 0.5 1.5 2.5 3.5]';
211// method = "linear" (default)
212// ---------------------------
213// extrapolation "by_nan" (default)
214assert_checkequal(interp1(x,y, 0), [%nan %nan]);
215assert_checkequal(interp1(x,y,1.5), [1 1]);
216assert_checkequal(interp1(x,y, 0.5:5 ), [%nan 1 3 5 %nan; %nan 1 2 3 %nan]');
217assert_checkequal(interp1(x,y,(0.5:5)'), [%nan 1 3 5 %nan; %nan 1 2 3 %nan]');
218Ref = cat(3,[%nan 1 3 ; 5 %nan %nan], [%nan 1 2 ; 3 %nan %nan]);
219assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5]), Ref);
220Ref = matrix([%nan 1 3 5 %nan %nan %nan 1 2 3 %nan %nan],[1 3 2 2]);
221assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5])), Ref);
222
223// extrapolation "extrap" = "linear"
224assert_checkequal(interp1(x,y, 0, ,"extrap"), [-2 -0.5]);
225assert_checkequal(interp1(x,y,1.5, ,"extrap"), [1 1]);
226assert_checkequal(interp1(x,y, 0.5:4, ,"extrap"), [-1 1 3 5 ; 0 1 2 3]');
227assert_checkequal(interp1(x,y,(0.5:4)', ,"extrap"), [-1 1 3 5 ; 0 1 2 3]');
228Ref = cat(3, [-1 1 3 ; 5 7 9], [0 1 2 ; 3 4 5]);
229assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,"extrap"), Ref);
230Ref = matrix([-1 1 3 5 7 9 0 1 2 3 4 5],[1 3 2 2]);
231assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,"extrap"), Ref);
232
233// extrapolation by padding
234assert_checkequal(interp1(x,y, 0, ,%pi), [%pi %pi]);
235assert_checkequal(interp1(x,y,1.5, ,%pi), [1 1]);
236assert_checkequal(interp1(x,y, 0.5:5, ,%pi ), [%pi 1 3 5 %pi; %pi 1 2 3 %pi]');
237assert_checkequal(interp1(x,y,(0.5:5)', ,%pi), [%pi 1 3 5 %pi; %pi 1 2 3 %pi]');
238Ref = cat(3,[%pi 1 3 ; 5 %pi %pi], [%pi 1 2 ; 3 %pi %pi]);
239assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,%pi), Ref);
240Ref = matrix([%pi 1 3 5 %pi %pi %pi 1 2 3 %pi %pi],[1 3 2 2]);
241assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,%pi), Ref);
242
243// periodic extrapolation
244y = [0 2 4 0 ; 0.5 1.5 2.5 0.5]';
245assert_checkalmostequal(interp1(x,y, 5, ,"periodic"), [2 1.5], 2*%eps);
246assert_checkalmostequal(interp1(x,y,1.5, ,"periodic"), [1 1], 2*%eps);
247assert_checkalmostequal(interp1(x,y, 0.5:5, ,"periodic" ), [2 1 3 2 1; 1.5 1 2 1.5 1]', 2*%eps);
248assert_checkalmostequal(interp1(x,y,(0.5:5)', ,"periodic"), [2 1 3 2 1; 1.5 1 2 1.5 1]', 2*%eps);
249Ref = cat(3,[2 1 3 ; 2 1 3], [1.5 1 2 ; 1.5 1 2]);
250assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,"periodic"), Ref, 2*%eps);
251Ref = matrix([2 1 3 2 1 3 1.5 1 2 1.5 1 2],[1 3 2 2]);
252assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,"periodic"), Ref, 2*%eps);
253
254// extrapolation by edgevalue
255e = "edgevalue";
256y = [0 2 4 6 ; 0.5 1.5 2.5 3.5]';
257assert_checkequal(interp1(x,y, 0, , e), [0 0.5]);
258assert_checkequal(interp1(x,y,1.5, , e), [1 1]);
259assert_checkequal(interp1(x,y, 0.5:5, , e ), [0 1 3 5 6 ; 0.5 1 2 3 3.5]');
260assert_checkequal(interp1(x,y,(0.5:5)', , e), [0 1 3 5 6 ; 0.5 1 2 3 3.5]');
261Ref = cat(3,[0 1 3 ; 5 6 6], [0.5 1 2 ; 3 3.5 3.5]);
262assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], , e), Ref);
263Ref = matrix([0 1 3 5 6 6 0.5 1 2 3 3.5 3.5],[1 3 2 2]);
264assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,e), Ref);
265
266// method = "nearest"
267// ------------------
268m = "nearest";
269// extrapolation "by_nan" (default)
270assert_checkequal(interp1(x,y, [0 5.4], m), [%nan %nan ; %nan %nan]);
271assert_checkequal(interp1(x,y,1.5, m), [2 1.5]);
272assert_checkequal(interp1(x,y,[0.5:4], m), [%nan %nan ; 2 1.5 ; 4 2.5 ; 6 3.5]);
273assert_checkequal(interp1(x,y,[0.5:4]', m),[%nan %nan ; 2 1.5 ; 4 2.5 ; 6 3.5]);
274Ref = cat(3, [%nan 2 4 ; 6 %nan %nan], [%nan 1.5 2.5 ; 3.5 %nan %nan]);
275assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m), Ref); // OErr
276Ref = matrix([%nan 2 4 6 %nan %nan %nan 1.5 2.5 3.5 %nan %nan],[1 3 2 2]);
277assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m), Ref);
278
279// extrapolation "extrap" = "nearest" = "edgevalue"
280for e = ["extrap" "nearest" "edgevalue"]
281 assert_checkequal(interp1(x,y, [0 5.4], m, e), [0 0.5 ; 6 3.5]);
282 assert_checkequal(interp1(x,y,1.5, m, e), [2 1.5]);
283 assert_checkequal(interp1(x,y,[0.5:4], m, e), y);
284 assert_checkequal(interp1(x,y,[0.5:4]', m, e), y);
285 Ref = cat(3, [0 2 4 ; 6 6 6], [0.5 1.5 2.5 ; 3.5 3.5 3.5]);
286 assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, e), Ref);
287 Ref = matrix([0 2 4 6 6 6 0.5 1.5 2.5 3.5 3.5 3.5],[1 3 2 2]);
288 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, e), Ref);
289end
290
291// extrapolation by padding
292assert_checkequal(interp1(x,y, [0 5.4], m, %pi), [%pi %pi ; %pi %pi]);
293assert_checkequal(interp1(x,y,1.5, m, %pi), [2 1.5]);
294assert_checkequal(interp1(x,y,[0.5:4], m, %pi), [%pi %pi ; 2 1.5 ; 4 2.5 ; 6 3.5]);
295assert_checkequal(interp1(x,y,[0.5:4]', m, %pi),[%pi %pi ; 2 1.5 ; 4 2.5 ; 6 3.5]);
296Ref = cat(3, [%pi 2 4 ; 6 %pi %pi], [%pi 1.5 2.5 ; 3.5 %pi %pi]);
297assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, %pi), Ref); // OErr
298Ref = matrix([%pi 2 4 6 %pi %pi %pi 1.5 2.5 3.5 %pi %pi],[1 3 2 2]);
299assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, %pi), Ref); // OErr
300
301// method = "spline"
302// -----------------
303m = "spline";
304// extrapolation default == "spline" (=="natural")
305assert_checkalmostequal(interp1(x,y, [0 5.4], m), [-2 -0.5 ; 8.8 4.9], 5*%eps);
306assert_checkequal(interp1(x,y,1.5, m), [1 1]);
307assert_checkequal(interp1(x,y,[0.5:4], m), [-1 1 3 5 ; 0 1 2 3]');
308assert_checkequal(interp1(x,y,[0.5:4]', m), [-1 1 3 5 ; 0 1 2 3]');
309Ref = cat(3,[-1 1 3 ; 5 7 9], [0 1 2 ; 3 4 5]);
310assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m), Ref, 5*%eps);
311Ref = matrix([-1 1 3 5 7 9 0 1 2 3 4 5 ],[1 3 2 2]);
312assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m), Ref, 5*%eps);
313
314// extrapolation "extrap" == "spline"
315assert_checkalmostequal(interp1(x,y, [0 5.4], m, "extrap"), [-2 -0.5 ; 8.8 4.9], 5*%eps);
316assert_checkequal(interp1(x,y,1.5, m, "extrap"), [1 1]);
317assert_checkequal(interp1(x,y,[0.5:4], m, "extrap"), [-1 1 3 5 ; 0 1 2 3]');
318assert_checkequal(interp1(x,y,[0.5:4]', m, "extrap"), [-1 1 3 5 ; 0 1 2 3]');
319Ref = cat(3,[-1 1 3 ; 5 7 9], [0 1 2 ; 3 4 5]);
320assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, "extrap"), Ref, 5*%eps);
321Ref = matrix([-1 1 3 5 7 9 0 1 2 3 4 5 ],[1 3 2 2]);
322assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, "extrap"), Ref, 5*%eps);
323
324// extrapolation by padding
325assert_checkalmostequal(interp1(x,y, [0 5.4], m, %pi), [%pi %pi ; %pi %pi], 5*%eps);
326assert_checkequal(interp1(x,y,1.5, m, %pi), [1 1]);
327assert_checkequal(interp1(x,y,[0.5:4], m, %pi), [%pi 1 3 5 ; %pi 1 2 3]');
328assert_checkequal(interp1(x,y,[0.5:4]', m, %pi), [%pi 1 3 5 ; %pi 1 2 3]');
329Ref = cat(3,[%pi 1 3 ; 5 %pi %pi], [%pi 1 2 ; 3 %pi %pi]);
330assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, %pi), Ref, 5*%eps);
331Ref = matrix([%pi 1 3 5 %pi %pi %pi 1 2 3 %pi %pi],[1 3 2 2]);
332assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, %pi), Ref);
333
334// extrapolation by edgevalue
335e = "edgevalue";
336assert_checkequal(interp1(x,y, [0 5.4], m, e), [0 0.5 ; 6 3.5]);
337assert_checkequal(interp1(x,y,1.5, m, e), [1 1]);
338assert_checkequal(interp1(x,y,[0.5:4], m, e), [0 1 3 5 ; 0.5 1 2 3]');
339assert_checkequal(interp1(x,y,[0.5:4]', m, e), [0 1 3 5 ; 0.5 1 2 3]');
340Ref = cat(3, [0 1 3 ; 5 6 6], [0.5 1 2 ; 3 3.5 3.5]);
341assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, e), Ref);
342Ref = matrix([0 1 3 5 6 6 0.5 1 2 3 3.5 3.5], [1 3 2 2]);
343assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, e), Ref);
344
345// periodic extrapolation
346y = [0 2 4 0 ; 0.5 1.5 2.5 0.5]';
347assert_checkalmostequal(interp1(x,y, [0 5.4], m, "periodic"), [4 2.5 ; 3.136 2.068], 2*%eps);
348assert_checkalmostequal(interp1(x,y,1.5, m, "periodic"), [0.625 0.8125], 5*%eps);
349Ref = [3.125 2.0625 ; 0.625 0.8125 ; 3.375 2.1875 ; 3.125 2.0625];
350assert_checkequal(interp1(x,y,[0.5:4], m, "periodic"), Ref);
351assert_checkequal(interp1(x,y,[0.5:4]', m, "periodic"), Ref);
352Ref = cat(3, [3.125 0.625 3.375 ; 3.125 0.625 3.375], ..
353 [2.0625 0.8125 2.1875 ; 2.0625 0.8125 2.1875]);
354assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], m, "periodic"), Ref, 5*%eps);
355Ref = matrix([3.125 0.625 3.375 3.125 0.625 3.375 2.0625 0.8125 2.1875 2.0625 0.8125 2.1875],[1 3 2 2]);
356assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, "periodic"), Ref, 5*%eps);
357
358
359// ===================
360// y is an hypermatrix
361// ===================
362y = cat(4,[0 2 4 6 ; 0.5 1.5 2.5 3.5]', [1:4 ; 2:5]');
363// method = "linear" (default)
364// ---------------------------
365// extrapolation "by_nan" (default)
366assert_checkequal(interp1(x,y, 0), cat(4,[%nan %nan],[%nan %nan])); // OWrong size and values
367assert_checkequal(interp1(x,y,1.5), cat(4,[1 1],[1.5 2.5])); // OWrong size
368Ref = cat(4,[%nan 1 3 5 ; %nan 1 2 3]', [%nan 1.5:3.5 ; %nan 2.5:4.5]');
369assert_checkequal(interp1(x,y, 0.5:4), Ref);
370assert_checkequal(interp1(x,y,(0.5:4)'), Ref);
371Ref = matrix([%nan 5 1 %nan 3 %nan %nan 3 1 %nan 2 %nan %nan 3.5 1.5 %nan ..
372 2.5 %nan %nan 4.5 2.5 %nan 3.5 %nan],[2 3 2 1 2]);
373assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5]), Ref);
374Ref = matrix([%nan 1 3 5 %nan %nan %nan 1 2 3 %nan %nan %nan 1.5 2.5 3.5 ..
375 %nan %nan %nan 2.5 3.5 4.5 %nan %nan],[1 3 2 2 1 2]);
376assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5])), Ref);
377
378// extrapolation "extrap" = "linear"
379assert_checkequal(interp1(x,y, 0, ,"extrap"), cat(4, [-2 -0.5],[0 1])); // OWrong size
380assert_checkequal(interp1(x,y,1.5, ,"extrap"), cat(4, [1 1],[1.5 2.5])); // OWrong size
381Ref = cat(4, [-1 1 3 5 ; 0 1 2 3]',[0.5:3.5 ; 1.5:4.5]');
382assert_checkequal(interp1(x,y, 0.5:4, ,"extrap"), Ref);
383assert_checkequal(interp1(x,y,(0.5:4)', ,"extrap"), Ref);
384
385Ref = matrix([-1 5 1 7 3 9 0 3 1 4 2 5 0.5 3.5 1.5 4.5 2.5 5.5 1.5 4.5 2.5 5.5 3.5 6.5],..
386 [2 3 2 1 2]);
387assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], ,"extrap"), Ref);
388Ref = matrix([-1 1 3 5 7 9 0 1 2 3 4 5 0.5 1.5 2.5 3.5 4.5 5.5 1.5 2.5 3.5 4.5 5.5 6.5],..
389 [1 3 2 2 1 2]);
390assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), ,"extrap"), Ref);
391
392// extrapolation by padding
393assert_checkequal(interp1(x,y, 0, , %pi), cat(4,[%pi %pi],[%pi %pi])); // OWrong size and values
394assert_checkequal(interp1(x,y,1.5, , %pi), cat(4,[1 1],[1.5 2.5]));
395Ref = cat(4,[%pi 1 3 5 ; %pi 1 2 3]', [%pi 1.5:3.5 ; %pi 2.5:4.5]');
396assert_checkequal(interp1(x,y, 0.5:4, , %pi), Ref);
397assert_checkequal(interp1(x,y,(0.5:4)', , %pi), Ref);
398Ref = matrix([%pi 5 1 %pi 3 %pi %pi 3 1 %pi 2 %pi %pi 3.5 1.5 %pi ..
399 2.5 %pi %pi 4.5 2.5 %pi 3.5 %pi],[2 3 2 1 2]);
400assert_checkequal(interp1(x,y,[0.5 1.5 2.5;3.5 4.5 5.5], , %pi), Ref);
401Ref = matrix([%pi 1 3 5 %pi %pi %pi 1 2 3 %pi %pi %pi 1.5 2.5 3.5 ..
402 %pi %pi %pi 2.5 3.5 4.5 %pi %pi],[1 3 2 2 1 2]);
403assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), , %pi), Ref);
404
405// periodic extrapolation
406y = cat(4,[0 2 4 0 ; 0.5 1.5 2.5 0.5]', [1 2 3 1 ; 2 3 4 2]');
407assert_checkequal(interp1(x,y, 0, , "periodic"), cat(4,[4 2.5],[3 4]));
408assert_checkequal(interp1(x,y,1.5, , "periodic"), cat(4,[1 1],[1.5 2.5]));
409Ref = cat(4,[2 1 3 2 ; 1.5 1 2 1.5]', [2 1.5 2.5 2 ; 3 2.5 3.5 3]');
410assert_checkequal(interp1(x,y, 0.5:4, , "periodic"), Ref);
411assert_checkequal(interp1(x,y,(0.5:4)', , "periodic"), Ref);
412Ref = matrix([2 2 1 1 3 3 1.5 1.5 1 1 2 2 2 2 1.5 1.5 ..
413 2.5 2.5 3 3 2.5 2.5 3.5 3.5],[2 3 2 1 2]);
414assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], , "periodic"), Ref, 2*%eps);
415Ref = matrix([2 1 3 2 1 3 1.5 1 2 1.5 1 2 2 1.5 2.5 2 ..
416 1.5 2.5 3 2.5 3.5 3 2.5 3.5],[1 3 2 2 1 2]);
417assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), , "periodic"), Ref, 2*%eps);
418
419// extrapolation by edgevalue
420e = "edgevalue";
421y = cat(4,[0 2 4 6 ; 0.5 1.5 2.5 3.5]', [1:4 ; 2:5]');
422assert_checkequal(interp1(x,y, 0, , e) , cat(4,[0 0.5],[1 2]));
423assert_checkequal(interp1(x,y,1.5, , e), cat(4,[1 1],[1.5 2.5]));
424Ref = cat(4,[0 1 3 5 ; 0.5 1 2 3]', [1 1.5:3.5 ; 2 2.5:4.5]');
425assert_checkequal(interp1(x,y, 0.5:4, , e), Ref);
426assert_checkequal(interp1(x,y,(0.5:4)', , e), Ref);
427Ref = matrix([0 5 1 6 3 6 0.5 3 1 3.5 2 3.5 1 3.5 1.5 4 ..
428 2.5 4 2 4.5 2.5 5 3.5 5],[2 3 2 1 2]);
429assert_checkequal(interp1(x, y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], , e), Ref);
430Ref = matrix([0 1 3 5 6 6 0.5 1 2 3 3.5 3.5 1 1.5 2.5 3.5 ..
431 4 4 2 2.5 3.5 4.5 5 5],[1 3 2 2 1 2]);
432assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), , e), Ref);
433
434// method = "nearest"
435// ------------------
436m = "nearest";
437// extrapolation "by_nan" (default)
438assert_checkequal(interp1(x,y, [0 5.4], m), %nan*ones(2,2,1,2));
439assert_checkequal(interp1(x,y,1.5, m), cat(4,[2 1.5],[2 3]));
440Ref = cat(4, [%nan %nan ; 2 1.5 ; 4 2.5 ; 6 3.5], [%nan %nan ; 2 3 ; 3 4 ; 4 5]);
441assert_checkequal(interp1(x,y,[0.5:4], m), Ref);
442assert_checkequal(interp1(x,y,[0.5:4]', m), Ref);
443Ref = matrix([%nan 6 2 %nan 4 %nan %nan 3.5 1.5 %nan 2.5 %nan %nan 4 2 %nan 3 %nan %nan 5 3 %nan 4 %nan],[2 3 2 1 2]);
444assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m), Ref); // OErr
445Ref = matrix([%nan 2 4 6 %nan %nan %nan 1.5 2.5 3.5 %nan %nan %nan 2 3 4 %nan %nan %nan 3 4 5 %nan%nan],[1 3 2 2 1 2]);
446assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m), Ref); // OErr
447
448// extrapolation "extrap" = "nearest"
449for e = ["extrap" "nearest" "edgevalue"]
450 assert_checkequal(interp1(x,y, [0 5.4], m, e), cat(4,[0 0.5 ; 6 3.5],[1 2 ; 4 5]));
451 assert_checkequal(interp1(x,y,1.5, m, e), cat(4,[2 1.5],[2 3])); // OWrong values
452 assert_checkequal(interp1(x,y,[0.5:4], m, e), y);
453 assert_checkequal(interp1(x,y,[0.5:4]', m, e), y);
454 Ref = matrix([0 6 2 6 4 6 0.5 3.5 1.5 3.5 2.5 3.5 1 4 2 4 3 4 2 5 3 5 4 5], [2 3 2 1 2]);
455 assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, e), Ref);
456 Ref = matrix([0 2 4 6 6 6 0.5 1.5 2.5 3.5 3.5 3.5 1 2 3 4 4 4 2 3 4 5 5 5],[1 3 2 2 1 2]);
457 assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, e), Ref);
458end
459
460// extrapolation by padding
461assert_checkequal(interp1(x,y, [0 5.4], m, %pi), %pi*ones(2,2,1,2));
462assert_checkequal(interp1(x,y,1.5, m, %pi), cat(4,[2 1.5],[2 3]));
463Ref = cat(4, [%pi %pi ; 2 1.5 ; 4 2.5 ; 6 3.5], [%pi %pi ; 2 3 ; 3 4 ; 4 5]);
464assert_checkequal(interp1(x,y,[0.5:4], m, %pi), Ref);
465assert_checkequal(interp1(x,y,[0.5:4]', m, %pi), Ref);
466Ref = matrix([%pi 6 2 %pi 4 %pi %pi 3.5 1.5 %pi 2.5 %pi %pi 4 2 %pi 3 %pi %pi 5 3 %pi 4 %pi],[2 3 2 1 2]);
467assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, %pi), Ref); // OErr
468Ref = matrix([%pi 2 4 6 %pi %pi %pi 1.5 2.5 3.5 %pi %pi %pi 2 3 4 %pi %pi %pi 3 4 5 %pi%pi],[1 3 2 2 1 2]);
469assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, %pi), Ref); // OErr
470
471// method = "spline"
472// -----------------
473m = "spline";
474// extrapolation default == "spline" (=="natural")
475Ref = cat(4, [-2 -0.5 ; 8.8 4.9], [0 1 ; 5.4 6.4]);
476assert_checkalmostequal(interp1(x,y, [0 5.4], m), Ref, 5*%eps);
477assert_checkequal(interp1(x,y,1.5, m), cat(4,[1 1],[1.5 2.5])); // OWrong values
478Ref = cat(4,[-1 1 3 5 ; 0 1 2 3]',[0.5:3.5 ; 1.5:4.5]');
479assert_checkequal(interp1(x,y,[0.5:4], m), Ref);
480assert_checkequal(interp1(x,y,[0.5:4]', m), Ref);
481Ref = matrix([-1 5 1 7 3 9 0 3 1 4 2 5 0.5 3.5 1.5 4.5 2.5 5.5 1.5 4.5 2.5 5.5 3.5 6.5],..
482 [2 3 2 1 2]);
483assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m), Ref, 5*%eps);
484Ref = matrix([-1 1 3 5 7 9 0 1 2 3 4 5 0.5 1.5 2.5 3.5 4.5 5.5 1.5 2.5 3.5 4.5 5.5 6.5],..
485 [1 3 2 2 1 2]);
486assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m), Ref, 5*%eps);
487
488// extrapolation "extrap" == "spline"
489Ref = cat(4, [-2 -0.5 ; 8.8 4.9], [0 1 ; 5.4 6.4]);
490assert_checkalmostequal(interp1(x,y, [0 5.4], m, "extrap"), Ref, 5*%eps); // OWrong size
491assert_checkequal(interp1(x,y,1.5, m, "extrap"), cat(4, [1 1], [1.5 2.5])); // OWrong size
492Ref = cat(4, [-1 1 3 5; 0:3]', [0.5:3.5 ; 1.5:4.5]');
493assert_checkequal(interp1(x,y,[0.5:4], m, "extrap"), Ref);
494assert_checkequal(interp1(x,y,[0.5:4]', m, "extrap"), Ref);
495Ref = matrix([-1 5 1 7 3 9 0 3 1 4 2 5 0.5 3.5 1.5 4.5 2.5 5.5 1.5 4.5 2.5 5.5 3.5 6.5],..
496 [2 3 2 1 2]);
497assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, "extrap"), Ref, 5*%eps); // OWrong size
498Ref = matrix([-1 1 3 5 7 9 0 1 2 3 4 5 0.5 1.5 2.5 3.5 4.5 5.5 1.5 2.5 3.5 4.5 5.5 6.5],..
499 [1 3 2 2 1 2]);
500assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, "extrap"), Ref, 5*%eps);
501
502// extrapolation by padding
503assert_checkalmostequal(interp1(x,y, [0 5.4], m, %pi), %pi*ones(2,2,1,2), 5*%eps);
504assert_checkequal(interp1(x,y,1.5, m, %pi), cat(4,[1 1],[1.5 2.5])); // OWrong values
505Ref = cat(4,[%pi 1 3 5 ; %pi 1 2 3]',[%pi 1.5 2.5 3.5 ; %pi 2.5 3.5 4.5]');
506assert_checkequal(interp1(x,y,[0.5:4], m, %pi), Ref);
507assert_checkequal(interp1(x,y,[0.5:4]', m, %pi), Ref);
508Ref = matrix([%pi 5 1 %pi 3 %pi %pi 3 1 %pi 2 %pi %pi 3.5 1.5 %pi 2.5 %pi %pi 4.5 2.5 %pi 3.5 %pi], [2 3 2 1 2]);
509assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, %pi), Ref, 5*%eps);
510Ref = matrix([%pi 1 3 5 %pi %pi %pi 1 2 3 %pi %pi %pi 1.5 2.5 3.5 %pi %pi %pi 2.5 3.5 4.5 %pi %pi], [1 3 2 2 1 2]);
511assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, %pi), Ref, 5*%eps);
512
513// extrapolation by edgevalue
514e = "edgevalue";
515assert_checkequal(interp1(x,y, [0 5.4], m, e), cat(4, [0 0.5 ; 6 3.5],[1 2 ; 4 5]));
516assert_checkequal(interp1(x,y,1.5, m, e), cat(4,[1 1],[1.5 2.5]));
517Ref = cat(4, [0 1 3 5 ; 0.5 1 2 3]', [1 1.5 2.5 3.5 ; 2 2.5 3.5 4.5]');
518assert_checkequal(interp1(x,y,[0.5:4], m, e), Ref);
519assert_checkequal(interp1(x,y,[0.5:4]', m, e), Ref);
520Ref = matrix([0 5 1 6 3 6 0.5 3 1 3.5 2 3.5 1 3.5 1.5 4 2.5 4 2 4.5 2.5 5 3.5 5], [2 3 2 1 2]);
521assert_checkequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, e), Ref);
522Ref = matrix([0 1 3 5 6 6 0.5 1 2 3 3.5 3.5 1 1.5 2.5 3.5 4 4 2 2.5 3.5 4.5 5 5], [1 3 2 2 1 2]);
523assert_checkequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, e), Ref);
524
525// periodic extrapolation
526y = cat(4,[0 2 4 0 ; 0.5 1.5 2.5 0.5]', [1 2 3 1 ; 2 3 4 2]');
527Ref = cat(4, [4 2.5 ; 3.136 2.068], [3 4 ; 2.568 3.568]);
528assert_checkalmostequal(interp1(x,y, [0 5.4], m, "periodic"), Ref, 5*%eps);
529assert_checkequal(interp1(x,y,1.5, m, "periodic"), cat(4,[0.625 0.8125],[1.3125 2.3125]));
530Ref = cat(4, [3.125 2.0625 ; 0.625 0.8125 ; 3.375 2.1875 ; 3.125 2.0625], ..
531 [2.5625 3.5625; 1.3125 2.3125 ; 2.6875 3.6875 ; 2.5625 3.5625]);
532assert_checkequal(interp1(x,y,[0.5:4], m, "periodic"), Ref);
533assert_checkequal(interp1(x,y,[0.5:4]', m, "periodic"), Ref);
534Ref = matrix([3.125 3.125 0.625 0.625 3.375 3.375 2.0625 2.0625 0.8125 0.8125 ..
535 2.1875 2.1875 2.5625 2.5625 1.3125 1.3125 2.6875 2.6875 3.5625 3.5625 ..
536 2.3125 2.3125 3.6875 3.6875], [2 3 2 1 2]);
537assert_checkalmostequal(interp1(x,y,[0.5 1.5 2.5 ; 3.5 4.5 5.5], m, "periodic"), Ref, 5*%eps);
538Ref = matrix([3.125 0.625 3.375 3.125 0.625 3.375 2.0625 0.8125 2.1875 2.0625 ..
539 0.8125 2.1875 2.5625 1.3125 2.6875 2.5625 1.3125 2.6875 3.5625 2.3125 ..
540 3.6875 3.5625 2.3125 3.6875], [1 3 2 2 1 2]);
541assert_checkalmostequal(interp1(x, y, cat(3,[0.5 1.5 2.5],[3.5 4.5 5.5]), m, "periodic"), Ref, 5*%eps);
542
543
544// y is complex
545// ============
546// method = "linear" (default)
547// ---------------------------
548assert_checkequal(interp1(0:3, (1:4)+%i, 0.5:2.5), (1.5:3.5)+%i);
549// extrapolation linear
550assert_checkequal(interp1(0:3,(1:4)+%i, 0.5:2:5, , "extrap"), [1.5 3.5 5.5]+%i);
551// Default padding with %nan
552assert_checkequal(interp1((1:4)+%i, 0.5:2.5), [%nan, 1.5+%i, 2.5+%i]);
553// Padding real
554assert_checkequal(interp1(1:3, (2:4)+%i, 0.5:2.5,,7), [7, 2.5+%i, 3.5+%i]);
555// Padding complex
556assert_checkequal(interp1(1:3, (2:4)+%i, [0.5 2.5 4.5],,7-%i), [7-%i, 3.5+%i, 7-%i]);
557
558// method = "nearest"
559// ------------------
560assert_checkequal(interp1(0:3, (1:4)+%i, 0.5:2.5, "nearest"), (2:4)+%i);
561// Default padding with %nan
562assert_checkequal(interp1((1:4)+%i, 0.5:2.5, "nearest"), [%nan, 2+%i, 3+%i]);
563// extrapolation with nearest
564assert_checkequal(interp1((1:4)+%i, 0.5:2:5, "nearest", "extrap"), [1 3 4]+%i);
565// Padding real
566assert_checkequal(interp1(1:3,(2:4)+%i, 0.5:2.5, "nearest", 7), [7, 3+%i, 4+%i]);
567// Padding complex
568assert_checkequal(interp1(1:3,(2:4)+%i, [0.5 2.5 4.5], "nearest",7-%i), [7-%i, 4+%i, 7-%i]);
569
570// method = "spline"
571// -----------------
572assert_checkequal(interp1(0:3, (1:4)+%i, 0.5:2.5, "spline"), (1.5:3.5)+%i);
573// Default padding with spline
574assert_checkequal(interp1((1:4)+%i, 0.5:2.5, "spline"), (0.5:2.5)+%i);
575assert_checkequal(interp1((1:4)+%i, 0.5:2.5, "spline", "extrap"), (0.5:2.5)+%i);
576// Padding real
577assert_checkequal(interp1(1:3,(2:4)+%i, 0.5:2.5, "spline", 7), [7, 2.5+%i, 3.5+%i]);
578// Padding complex
579assert_checkequal(interp1(1:3,(2:4)+%i, [0.5 2.5 4.5], "spline",7-%i), [7-%i, 3.5+%i, 7-%i]);
580
581
582// ==============
583// Error messages
584// ==============
585msg = _("%s: Wrong number of input arguments: Must be between %d and %d.\n")
586msg = msprintf(msg, "interp1", 2, 5);
587assert_checkerror("interp1()", msg);
588assert_checkerror("interp1(1)", msg);
589assert_checkerror("interp1(1,2,3,4,5,6)", msg);
590// x
591msg = msprintf(_("%s: Argument #%d: Real numbers expected.\n"),"interp1",1);
592assert_checkerror("interp1([1 %i],[1 2],1.5)", msg);
593msg = msprintf(_("%s: Argument #%d: Vector expected.\n"), "interp1", 1);
594assert_checkerror("interp1([1 2; 3 4],[1 2],1.5)", msg);
595assert_checkerror("interp1(1,[1 2],1.5)", msg);
596msg = msprintf(_("%s: Argument #%d: Nan value forbidden.\n"), "interp1", 1);
597assert_checkerror("interp1([1 2 %nan],[1 2],1.5)", msg);
598// y
599msg = msprintf(_("%s: Argument #%d: vector or at least 2 rows expected.\n"), "interp1", 2);
600assert_checkerror("interp1(1,1.5)", msg);
601msg = msprintf(_("%s: Argument #%d: Decimal or complex numbers expected.\n"), "interp1",1);
602assert_checkerror("interp1([1 %z],1.5)", msg);
603assert_checkerror("interp1([%t %f],1.5)", msg);
604//
605msg = msprintf(_("%s: Arguments #%d and #%d: Same numbers of elements expected.\n"),"interp1", 1, 2);
606assert_checkerror("interp1([1 2],1:3,3)", msg);
607assert_checkerror("interp1(1:3,1:2,3)", msg);
608msg = msprintf(_("%s: Arguments #%d and #%d: Same numbers of rows expected.\n"), "interp1", 1, 2);
609assert_checkerror("interp1(1:3,rand(2,3),3)", msg);
610 // periodicity
611msg = msprintf(_("%s: Argument #%d: periodicity y($)==y(1) expected.\n"), "interp1", 2);
612assert_checkerror("interp1(1:4,1:4,3,,""periodic"")", msg);
613msg = msprintf(_("%s: Argument #%d: periodicity y($)==y(1) expected.\n"), "interp1", 1);
614assert_checkerror("interp1(1:4,3,""linear"",""periodic"")", msg);
615msg = msprintf(_("%s: Argument #%d: periodicity y($,:)==y(1,:) expected.\n"), "interp1", 2);
616assert_checkerror("interp1(1:3,rand(3,2),3,,""periodic"")", msg);
617msg = msprintf(_("%s: Argument #%d: periodicity y($,:)==y(1,:) expected.\n"), "interp1", 1);
618assert_checkerror("interp1(rand(3,2),3,""linear"",""periodic"")", msg);
619
620// xp
621msg = msprintf(_("%s: Argument #%d: Decimal numbers expected.\n"), "interp1",2);
622assert_checkerror("interp1(1:3,[%i 1 1])", msg);
623msg = msprintf(_("%s: Argument #%d: Decimal numbers expected.\n"), "interp1",3);
624assert_checkerror("interp1(1:3,rand(3,2),%i)", msg);
625assert_checkerror("interp1(1:3,2:4,,""linear"")", msg);
626// method
627msg = msprintf(_("%s: Argument #%d: Scalar (1 element) expected.\n"), "interp1", 3);
628assert_checkerror("interp1(1:3,2,[""linear"" ""spline""])", msg);
629msg = msprintf(_("%s: Argument #%d: Scalar (1 element) expected.\n"), "interp1", 4);
630assert_checkerror("interp1(1:3,2:4,5,[""linear"" ""spline""])", msg);
631Set = """linear"",""spline"",""nearest""";
632msg = msprintf(_("%s: Argument #%d: Must be in the set {%s}.\n"), "interp1", 3, Set);
633assert_checkerror("interp1(1:3,5,""splines"")", msg);
634msg = msprintf(_("%s: Argument #%d: Must be in the set {%s}.\n"), "interp1", 4, Set);
635assert_checkerror("interp1(2:4,1:3,5,""splines"")", msg);
636assert_checkerror("interp1(2:4,1:3,5,%z)", msg);
637// extrapolation
638msg0 = _("%s: Argument #%d: scalar number or one of {%s} expected.\n");
639msg = msprintf(msg0, "interp1", 5, """extrap"",""linear"",""periodic"",""edgevalue""");
640assert_checkerror("interp1(1:3,2:4,5,,""abc"")", msg);
641assert_checkerror("interp1(1:3,2:4,5,,%z)", msg);
642assert_checkerror("interp1(1:3,2:4,5,,1:2)", msg);
643assert_checkerror("interp1(1:3,2:4,5,,int8(3))", msg);
644msg = msprintf(msg0, "interp1", 5, """extrap"",""edgevalue"",""nearest""");
645assert_checkerror("interp1(1:3,2:4,5,""nearest"",""linear"")", msg);
646assert_checkerror("interp1(1:3,2:4,5,""nearest"",""periodic"")", msg);
647msg = msprintf(msg0, "interp1", 4, """extrap"",""linear"",""periodic"",""edgevalue""");
648assert_checkerror("interp1(1:3,2:4,""linear"",""abc"")", msg);
649assert_checkerror("interp1(1:3,2:4,""linear"",%z)", msg);
650assert_checkerror("interp1(1:3,2:4,""linear"",1:2)", msg);
651assert_checkerror("interp1(1:3,2:4,""linear"",int8(3))", msg);
652msg = msprintf(msg0, "interp1", 4, """extrap"",""edgevalue"",""nearest""");
653assert_checkerror("interp1(1:3,2:4,""nearest"",""linear"")", msg);
654assert_checkerror("interp1(1:3,2:4,""nearest"",""periodic"")", msg);