source: libpipi/trunk/pipi/resample/bicubic.c @ 4696

Last change on this file since 4696 was 4696, checked in by Sam Hocevar, 9 years ago

Implement bicubic resampling. Lacks some blurring in the pre-pass, maybe.

  • Property svn:keywords set to Id
File size: 5.5 KB
Line 
1/*
2 *  libpipi       Pathetic image processing interface library
3 *  Copyright (c) 2004-2010 Sam Hocevar <sam@hocevar.net>
4 *                All Rights Reserved
5 *
6 *  $Id: bicubic.c 4696 2010-10-17 00:57:27Z sam $
7 *
8 *  This library is free software. It comes without any warranty, to
9 *  the extent permitted by applicable law. You can redistribute it
10 *  and/or modify it under the terms of the Do What The Fuck You Want
11 *  To Public License, Version 2, as published by Sam Hocevar. See
12 *  http://sam.zoy.org/wtfpl/COPYING for more details.
13 */
14
15/*
16 * bicubic.c: Bicubic image resizing functions
17 */
18
19#include "config.h"
20
21#include <stdlib.h>
22#include <string.h>
23
24#include "pipi.h"
25#include "pipi_internals.h"
26
27static inline int min_int(int a, int b) { return a < b ? a : b; }
28static inline int max_int(int a, int b) { return a > b ? a : b; }
29
30pipi_image_t *pipi_resize_bicubic(pipi_image_t *src, int w, int h)
31{
32    float *srcdata, *dstdata, *p0, *p1, *p2, *p3;
33    pipi_image_t *dst;
34    pipi_pixels_t *srcp, *dstp;
35    int x, y, i, sw, dw, sh, dh, i0, i1, i2, i3;
36    float scalex, scaley;
37
38    srcp = pipi_get_pixels(src, PIPI_PIXELS_RGBA_F32);
39    srcdata = (float *)srcp->pixels;
40
41    dst = pipi_new(w, h);
42    dstp = pipi_get_pixels(dst, PIPI_PIXELS_RGBA_F32);
43    dstdata = (float *)dstp->pixels;
44
45    sw = src->w; sh = src->h;
46    dw = dst->w; dh = dst->h;
47
48    scalex = dw > 1 ? (float)(sw - 1) / (dw - 1) : 1.0f;
49    scaley = dh > 1 ? (float)(sh - 1) / (dh - 1) : 1.0f;
50
51    for(y = 0; y < dh; y++)
52    {
53        float yfloat = scaley * y;
54        int yint = (int)yfloat;
55        float y1 = yfloat - yint;
56
57        p0 = srcdata + 4 * sw * min_int(max_int(0, yint - 1), sh - 1);
58        p1 = srcdata + 4 * sw * min_int(max_int(0, yint    ), sh - 1);
59        p2 = srcdata + 4 * sw * min_int(max_int(0, yint + 1), sh - 1);
60        p3 = srcdata + 4 * sw * min_int(max_int(0, yint + 2), sh - 1);
61
62        for (x = 0; x < dw; x++)
63        {
64            float xfloat = scalex * x;
65            int xint = (int)xfloat;
66            float x1 = xfloat - xint;
67
68            i0 = 4 * min_int(max_int(0, xint - 1), sw - 1);
69            i1 = 4 * min_int(max_int(0, xint    ), sw - 1);
70            i2 = 4 * min_int(max_int(0, xint + 1), sw - 1);
71            i3 = 4 * min_int(max_int(0, xint + 2), sw - 1);
72
73            for (i = 0; i < 4; i++, i0++, i1++, i2++, i3++)
74            {
75                float a00 = p1[i1];
76                float a01 = .5f * (p2[i1] - p0[i1]);
77                float a02 = p0[i1] - 2.5f * p1[i1]
78                            + 2.f * p2[i1] - .5f * p3[i1];
79                float a03 = .5f * (p3[i1] - p0[i1]) + 1.5f * (p1[i1] - p2[i1]);
80
81                float a10 = .5f * (p1[i2] - p1[i0]);
82                float a11 = .25f * (p0[i0] - p2[i0] - p0[i2] + p2[i2]);
83                float a12 = .5f * (p0[i2] - p0[i0]) + 1.25f * (p1[i0] - p1[i2])
84                            + .25f * (p3[i0] - p3[i2]) + p2[i2] - p2[i0];
85                float a13 = .25f * (p0[i0] - p3[i0] - p0[i2] + p3[i2])
86                            + .75f * (p2[i0] - p1[i0] + p1[i2] - p2[i2]);
87
88                float a20 = p1[i0] - 2.5f * p1[i1]
89                            + 2.f * p1[i2] - .5f * p1[i3];
90                float a21 = .5f * (p2[i0] - p0[i0]) + 1.25f * (p0[i1] - p2[i1])
91                            + .25f * (p0[i3] - p2[i3]) - p0[i2] + p2[i2];
92                float a22 = p0[i0] - p3[i2] - 2.5f * (p1[i0] + p0[i1])
93                            + 2.f * (p2[i0] + p0[i2]) - .5f * (p3[i0] + p0[i3])
94                            + 6.25f * p1[i1] - 5.f * (p2[i1] + p1[i2])
95                            + 1.25f * (p3[i1] + p1[i3])
96                            + 4.f * p2[i2] - p2[i3] + .25f * p3[i3];
97                float a23 = 1.5f * (p1[i0] - p2[i0]) + .5f * (p3[i0] - p0[i0])
98                            + 1.25f * (p0[i1] - p3[i1])
99                            + 3.75f * (p2[i1] - p1[i1]) + p3[i2] - p0[i2]
100                            + 3.f * (p1[i2] - p2[i2]) + .25f * (p0[i3] - p3[i3])
101                            + .75f * (p2[i3] - p1[i3]);
102
103                float a30 = .5f * (p1[i3] - p1[i0]) + 1.5f * (p1[i1] - p1[i2]);
104                float a31 = .25f * (p0[i0] - p2[i0]) + .25f * (p2[i3] - p0[i3])
105                            + .75f * (p2[i1] - p0[i1] + p0[i2] - p2[i2]);
106                float a32 = -.5f * p0[i0] + 1.25f * p1[i0] - p2[i0]
107                            + .25f * p3[i0] + 1.5f * p0[i1] - 3.75f * p1[i1]
108                            + 3.f * p2[i1] - .75f * p3[i1] - 1.5f * p0[i2]
109                            + 3.75f * p1[i2] - 3.f * p2[i2] + .75f * p3[i2]
110                            + .5f * p0[i3] - 1.25f * p1[i3] + p2[i3]
111                            - .25f * p3[i3];
112                float a33 = .25f * p0[i0] - .75f * p1[i0] + .75f * p2[i0]
113                            - .25f * p3[i0] - .75f * p0[i1] + 2.25f * p1[i1]
114                            - 2.25f * p2[i1] + .75f * p3[i1] + .75f * p0[i2]
115                            - 2.25f * p1[i2] + 2.25f * p2[i2] - .75f * p3[i2]
116                            - .25f * p0[i3] + .75f * p1[i3] - .75f * p2[i3]
117                            + .25f * p3[i3];
118
119                float y2 = y1 * y1; float y3 = y2 * y1;
120                float x2 = x1 * x1; float x3 = x2 * x1;
121
122                float p = a00 + a01 * y1 + a02 * y2 + a03 * y3 +
123                   + a10 * x1 + a11 * x1 * y1 + a12 * x1 * y2 + a13 * x1 * y3
124                   + a20 * x2 + a21 * x2 * y1 + a22 * x2 * y2 + a23 * x2 * y3
125                   + a30 * x3 + a31 * x3 * y1 + a32 * x3 * y2 + a33 * x3 * y3;
126                if (p < 0.0f) p = 0.0f;
127                if (p > 1.0f) p = 1.0f;
128
129                dstdata[(y * dw + x) * 4 + i] = p;
130            }
131        }
132    }
133
134    return dst;
135}
136
Note: See TracBrowser for help on using the repository browser.