Bilinear interpolation: Difference between revisions
m (→{{header|J}}) |
(→{{header|J}}: generalization) |
||
Line 105: | Line 105: | ||
<lang J> |
<lang J> |
||
Note 'FEA' |
Note 'FEA' |
||
Here we develop a general method to generate isoparametric interpolants. |
|||
⚫ | |||
⚫ | |||
⚫ | |||
⚫ | |||
at the coordinates within the element with the known values at the nodes. |
|||
⚫ | |||
⚫ | |||
Line 115: | Line 118: | ||
| | |
| | |
||
| (0,0) | |
| (0,0) | |
||
| * | |
| * | |
||
| | |
| | |
||
| | |
| | |
||
| | |
| | |
||
+---------------+ |
+---------------+ |
||
Line 126: | Line 129: | ||
f1( 1,-1) = 1, f1(all other corners) is 0. |
f1( 1,-1) = 1, f1(all other corners) is 0. |
||
... |
... |
||
The coordinates (xi,eta) are a vector length 2 as argument y |
|||
Supply to interpolate the function values at nodes 0 1 2 and 3 as the x argument. |
|||
⚫ | |||
Choose a shape function. |
|||
f0 =: 1r4 * [: */ 1 1&- NB. f0(xi,eta) = (1-xi)(1-eta)/4 |
|||
Use shape functions C0 + C1*xi + C2*eta + C3*xi*eta . |
|||
Given (xi,eta) as the vector y form a vector of the |
|||
f2 =: _1r4 * [: */ 1 _1&- NB. ... |
|||
coefficients of the constants (1, xi, eta, and their product) |
|||
f3 =: 1r4 * [: */ _1 _1&- |
|||
shape_function =: 1 , {. , {: , */ |
|||
functions =: f0`f1`f2`f3`:0 |
|||
CORNERS NB. are the ordered coordinates of the corners |
|||
⚫ | |||
_1 _1 |
|||
NB. prove the function values at the corners are correct. |
|||
1 _1 |
|||
assert (=i.4) -: functions"1 CORNERS NB. 4 by 4 identity matrix |
|||
_1 1 |
|||
1 1 |
|||
(=i.4) NB. rows of the identity matrix are the values of each shape functions at each corner |
|||
1 0 0 0 |
|||
0 1 0 0 |
|||
0 0 1 0 |
|||
0 0 0 1 |
|||
(=i.4x) %. shape_function"1 x: CORNERS NB. Compute the values of the constants as rational numbers. |
|||
1r4 1r4 1r4 1r4 |
|||
_1r4 1r4 _1r4 1r4 |
|||
_1r4 _1r4 1r4 1r4 |
|||
1r4 _1r4 _1r4 1r4 |
|||
This method extends to higher order interpolants having more nodes or to other dimensions. |
|||
⚫ | |||
⚫ | |||
mp =: +/ .* NB. matrix product |
|||
⚫ | |||
shape_function =: 1 , {. , {: , */ |
|||
COEFFICIENTS =: (=i.4) %. shape_function"1 CORNERS |
|||
shape_functions =: COEFFICIENTS mp shape_function |
|||
⚫ | |||
</lang> |
</lang> |
||
<pre> |
<pre> |
||
Note |
Note 'demonstrate the interpolant with a saddle' |
||
lower left has value 1, |
lower left has value 1, |
||
lower right: 2 |
lower right: 2 |
||
upper left: 2.2 |
upper left: 2.2 |
||
upper right: 0.7 |
upper right: 0.7 |
||
GRID is an array of coordinate pairs (hence rank 3) |
|||
ordered so that _1 _1 is in the lower left of the picture, |
|||
directed with 1 _1 to lower right. |
|||
) |
) |
||
Line 158: | Line 178: | ||
viewmat SADDLE |
viewmat SADDLE |
||
</pre> |
</pre> |
||
[[Image:J_bilinear_interpolant.jpg |
[[Image:J_bilinear_interpolant.jpg]] |
||
=={{header|Racket}}== |
=={{header|Racket}}== |
Revision as of 22:29, 9 December 2014
Bilinear interpolation is linear interpolation in 2 dimensions, and is typically used for image scaling and for 2D finite element analysis.
C
<lang c>#include <stdint.h> typedef struct {
uint32_t *pixels; unsigned int w; unsigned int h;
} image_t;
- define getByte(value, n) (value >> (n*8) & 0xFF)
uint32_t getpixel(image_t *image, unsigned int x, unsigned int y){
return image->pixels[(y*image->w)+x];
} float lerp(float s, float e, float t){return s+(e-s)*t;} float blerp(float c00, float c10, float c01, float c11, float tx, float ty){
return lerp(lerp(c00, c10, tx), lerp(c01, c11, tx), ty);
} void putpixel(image_t *image, unsigned int x, unsigned int y, uint32_t color){
image->pixels[(y*image->w) + x] = color;
} void scale(image_t *src, image_t *dst, float scalex, float scaley){
int newWidth = (int)src->w*scalex; int newHeight= (int)src->h*scaley; int x, y; for(x= 0, y=0; y < newHeight; x++){ if(x > newWidth){ x = 0; y++; } float gx = x / (float)(newWidth) * (src->w-1); float gy = y / (float)(newHeight) * (src->h-1); int gxi = (int)gx; int gyi = (int)gy; uint32_t result=0; uint32_t c00 = getpixel(src, gxi, gyi); uint32_t c10 = getpixel(src, gxi+1, gyi); uint32_t c01 = getpixel(src, gxi, gyi+1); uint32_t c11 = getpixel(src, gxi+1, gyi+1); uint8_t i; for(i = 0; i < 3; i++){ //((uint8_t*)&result)[i] = blerp( ((uint8_t*)&c00)[i], ((uint8_t*)&c10)[i], ((uint8_t*)&c01)[i], ((uint8_t*)&c11)[i], gxi - gx, gyi - gy); // this is shady result |= (uint8_t)blerp(getByte(c00, i), getByte(c10, i), getByte(c01, i), getByte(c11, i), gx - gxi, gy -gyi) << (8*i); } putpixel(dst,x, y, result); }
}</lang>
D
This uses the module from the Grayscale Image task.
<lang d>import grayscale_image;
/// Currently this accepts only a Grayscale image, for simplicity. Image!Gray rescaleGray(in Image!Gray src, in float scaleX, in float scaleY) pure nothrow @safe in {
assert(src !is null, "Input Image is null."); assert(src.nx > 1 && src.ny > 1, "Minimal input image size is 2x2."); assert(cast(uint)(src.nx * scaleX) > 0, "Output image width must be > 0."); assert(cast(uint)(src.ny * scaleY) > 0, "Output image height must be > 0.");
} body {
alias FP = float; static FP lerp(in FP s, in FP e, in FP t) pure nothrow @safe @nogc { return s + (e - s) * t; }
static FP blerp(in FP c00, in FP c10, in FP c01, in FP c11, in FP tx, in FP ty) pure nothrow @safe @nogc { return lerp(lerp(c00, c10, tx), lerp(c01, c11, tx), ty); }
immutable newWidth = cast(uint)(src.nx * scaleX); immutable newHeight = cast(uint)(src.ny * scaleY); auto result = new Image!Gray(newWidth, newHeight, true);
foreach (immutable y; 0 .. newHeight) foreach (immutable x; 0 .. newWidth) { immutable FP gx = x / FP(newWidth) * (src.nx - 1); immutable FP gy = y / FP(newHeight) * (src.ny - 1); immutable gxi = cast(uint)gx; immutable gyi = cast(uint)gy;
immutable c00 = src[gxi, gyi ]; immutable c10 = src[gxi + 1, gyi ]; immutable c01 = src[gxi, gyi + 1]; immutable c11 = src[gxi + 1, gyi + 1];
immutable pixel = blerp(c00, c10, c01, c11, gx - gxi, gy - gyi); result[x, y] = Gray(cast(ubyte)pixel); }
return result;
}
void main() {
const im = loadPGM!Gray(null, "lena.pgm"); im.rescaleGray(0.3, 0.1).savePGM("lena_smaller.pgm"); im.rescaleGray(1.3, 1.8).savePGM("lena_larger.pgm");
}</lang>
J
<lang J> Note 'FEA'
Here we develop a general method to generate isoparametric interpolants.
The interpolant is the dot product of the four shape function values evaluated at the coordinates within the element with the known values at the nodes. The sum of four shape functions of two variables (xi, eta) is 1 at each of four nodes. Let the base element have nodal coordinates (xi, eta) of +/-1.
2 3 (1,1) +---------------+ | | | | | (0,0) | | * | | | | | | | +---------------+ 0 1
determine f0(xi,eta), ..., f3(xi,eta). f0(-1,-1) = 1, f0(all other corners) is 0. f1( 1,-1) = 1, f1(all other corners) is 0. ...
Choose a shape function. Use shape functions C0 + C1*xi + C2*eta + C3*xi*eta . Given (xi,eta) as the vector y form a vector of the coefficients of the constants (1, xi, eta, and their product)
shape_function =: 1 , {. , {: , */
CORNERS NB. are the ordered coordinates of the corners _1 _1 1 _1 _1 1 1 1 (=i.4) NB. rows of the identity matrix are the values of each shape functions at each corner 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 (=i.4x) %. shape_function"1 x: CORNERS NB. Compute the values of the constants as rational numbers. 1r4 1r4 1r4 1r4 _1r4 1r4 _1r4 1r4 _1r4 _1r4 1r4 1r4 1r4 _1r4 _1r4 1r4
This method extends to higher order interpolants having more nodes or to other dimensions.
)
mp =: +/ .* NB. matrix product
CORNERS =: 21 A.-.+:#:i.4 shape_function =: 1 , {. , {: , */ COEFFICIENTS =: (=i.4) %. shape_function"1 CORNERS shape_functions =: COEFFICIENTS mp shape_function interpolate =: mp shape_functions </lang>
Note 'demonstrate the interpolant with a saddle' lower left has value 1, lower right: 2 upper left: 2.2 upper right: 0.7 ) GRID =: |.,~"0/~(%~i:)100 SADDLE =: 1 2 2.2 0.7 interpolate"_ 1 GRID viewmat SADDLE
File:J bilinear interpolant.jpg
Racket
This mimics the Wikipedia example. <lang racket>#lang racket (require images/flomap)
(define fm
(draw-flomap (λ (dc) (define (pixel x y color) (send dc set-pen color 1 'solid) (send dc draw-point (+ x .5) (+ y 0.5))) (send dc set-alpha 1) (pixel 0 0 "blue") (pixel 0 1 "red") (pixel 1 0 "red") (pixel 1 1 "green")) 2 2))
(flomap->bitmap
(build-flomap 4 250 250 (λ (k x y) (flomap-bilinear-ref fm k (+ 1/2 (/ x 250)) (+ 1/2 (/ y 250))))))</lang>