1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
|
/*
* Polynomial interpolation template
* Copyright 2019 Anton Khirnov <anton@khirnov.net>
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#if CONTIGUOUS
# define CONT cont
#else
# define CONT generic
#endif
#define JOIN3(a, b, c) a ## _ ## b ## _ ## c
#define FUNC2(name, cont, stencil) JOIN3(name, cont, stencil)
#define FUNC(name) FUNC2(name, CONT, STENCIL)
static void
FUNC(interp2d_transfer_line)(double *dst, ptrdiff_t dst_len,
const double *src, ptrdiff_t src_stride,
const ptrdiff_t *idx_x, const double *fact_x, const double *fact_y
#if !CONTIGUOUS
, ptrdiff_t dst_stride0, ptrdiff_t src_stride0
# define SSTRIDE1 src_stride0
# define DSTRIDE1 dst_stride0
#else
# define SSTRIDE1 1
# define DSTRIDE1 1
#endif
)
{
for (size_t x = 0; x < dst_len; x++) {
const double *src_data = src + SSTRIDE1 * idx_x[x];
double val = 0.0;
for (ptrdiff_t idx0 = 0; idx0 < STENCIL; idx0++) {
double tmp = 0.0;
for (ptrdiff_t idx1 = 0; idx1 < STENCIL; idx1++)
tmp += src_data[idx0 * src_stride + idx1 * SSTRIDE1] * fact_x[STENCIL * x + idx1];
val += tmp * fact_y[idx0];
}
dst[x * DSTRIDE1] = val;
}
}
static void
FUNC(interp1d_transfer_line)(double *dst, ptrdiff_t dst_len,
const double *src, const ptrdiff_t *idx_x, const double *fact_x
#if !CONTIGUOUS
, ptrdiff_t dst_stride0, ptrdiff_t src_stride0
#endif
)
{
for (size_t x = 0; x < dst_len; x++) {
const double *src_data = src + SSTRIDE1 * idx_x[x];
double val = 0.0;
for (ptrdiff_t idx = 0; idx < STENCIL; idx++)
val += src_data[idx * SSTRIDE1] * fact_x[STENCIL * x + idx];
dst[x * DSTRIDE1] = val;
}
}
#undef SSTRIDE1
#undef DSTRIDE1
#undef CONT
#undef JOIN3
#undef FUNC2
#undef FUNC
|
