Pillow

1
/*
2
 * The Python Imaging Library
3
 * $Id$
4
 *
5
 * apply convolution kernel to image
6
 *
7
 * history:
8
 * 1995-11-26 fl   Created, supports 3x3 kernels
9
 * 1995-11-27 fl   Added 5x5 kernels, copy border
10
 * 1999-07-26 fl   Eliminated a few compiler warnings
11
 * 2002-06-09 fl   Moved kernel definitions to Python
12
 * 2002-06-11 fl   Support floating point kernels
13
 * 2003-09-15 fl   Added ImagingExpand helper
14
 *
15
 * Copyright (c) Secret Labs AB 1997-2002.  All rights reserved.
16
 * Copyright (c) Fredrik Lundh 1995.
17
 *
18
 * See the README file for information on usage and redistribution.
19
 */
20

21
/*
22
 * FIXME: Support RGB and RGBA/CMYK modes as well
23
 * FIXME: Expand image border (current version leaves border as is)
24
 * FIXME: Implement image processing gradient filters
25
 */
26

27
#include "Imaging.h"
28

29
static inline UINT8
30
clip8(float in) {
31
    if (in <= 0.0) {
32
        return 0;
33
    }
34
    if (in >= 255.0) {
35
        return 255;
36
    }
37
    return (UINT8)in;
38
}
39

40
static inline INT32
41
clip32(float in) {
42
    if (in <= 0.0) {
43
        return 0;
44
    }
45
    if (in >= pow(2, 31) - 1) {
46
        return pow(2, 31) - 1;
47
    }
48
    return (INT32)in;
49
}
50

51
Imaging
52
ImagingExpand(Imaging imIn, int xmargin, int ymargin) {
53
    Imaging imOut;
54
    int x, y;
55
    ImagingSectionCookie cookie;
56

57
    if (xmargin < 0 && ymargin < 0) {
58
        return (Imaging)ImagingError_ValueError("bad kernel size");
59
    }
60

61
    imOut = ImagingNewDirty(
62
        imIn->mode, imIn->xsize + 2 * xmargin, imIn->ysize + 2 * ymargin
63
    );
64
    if (!imOut) {
65
        return NULL;
66
    }
67

68
#define EXPAND_LINE(type, image, yin, yout)                        \
69
    {                                                              \
70
        for (x = 0; x < xmargin; x++) {                            \
71
            imOut->image[yout][x] = imIn->image[yin][0];           \
72
        }                                                          \
73
        for (x = 0; x < imIn->xsize; x++) {                        \
74
            imOut->image[yout][x + xmargin] = imIn->image[yin][x]; \
75
        }                                                          \
76
        for (x = 0; x < xmargin; x++) {                            \
77
            imOut->image[yout][xmargin + imIn->xsize + x] =        \
78
                imIn->image[yin][imIn->xsize - 1];                 \
79
        }                                                          \
80
    }
81

82
#define EXPAND(type, image)                                                       \
83
    {                                                                             \
84
        for (y = 0; y < ymargin; y++) {                                           \
85
            EXPAND_LINE(type, image, 0, y);                                       \
86
        }                                                                         \
87
        for (y = 0; y < imIn->ysize; y++) {                                       \
88
            EXPAND_LINE(type, image, y, y + ymargin);                             \
89
        }                                                                         \
90
        for (y = 0; y < ymargin; y++) {                                           \
91
            EXPAND_LINE(type, image, imIn->ysize - 1, ymargin + imIn->ysize + y); \
92
        }                                                                         \
93
    }
94

95
    ImagingSectionEnter(&cookie);
96
    if (imIn->image8) {
97
        EXPAND(UINT8, image8);
98
    } else {
99
        EXPAND(INT32, image32);
100
    }
101
    ImagingSectionLeave(&cookie);
102

103
    ImagingCopyPalette(imOut, imIn);
104

105
    return imOut;
106
}
107

108
void
109
ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
110
#define KERNEL1x3(in0, x, kernel, d)                               \
111
    (_i2f(in0[x - d]) * (kernel)[0] + _i2f(in0[x]) * (kernel)[1] + \
112
     _i2f(in0[x + d]) * (kernel)[2])
113

114
    int x = 0, y = 0;
115

116
    memcpy(imOut->image[0], im->image[0], im->linesize);
117
    if (im->bands == 1) {
118
        // Add one time for rounding
119
        offset += 0.5;
120
        if (im->type == IMAGING_TYPE_INT32) {
121
            for (y = 1; y < im->ysize - 1; y++) {
122
                INT32 *in_1 = (INT32 *)im->image[y - 1];
123
                INT32 *in0 = (INT32 *)im->image[y];
124
                INT32 *in1 = (INT32 *)im->image[y + 1];
125
                INT32 *out = (INT32 *)imOut->image[y];
126

127
                out[0] = in0[0];
128
                for (x = 1; x < im->xsize - 1; x++) {
129
                    float ss = offset;
130
                    ss += KERNEL1x3(in1, x, &kernel[0], 1);
131
                    ss += KERNEL1x3(in0, x, &kernel[3], 1);
132
                    ss += KERNEL1x3(in_1, x, &kernel[6], 1);
133
                    out[x] = clip32(ss);
134
                }
135
                out[x] = in0[x];
136
            }
137
        } else {
138
            for (y = 1; y < im->ysize - 1; y++) {
139
                UINT8 *in_1 = (UINT8 *)im->image[y - 1];
140
                UINT8 *in0 = (UINT8 *)im->image[y];
141
                UINT8 *in1 = (UINT8 *)im->image[y + 1];
142
                UINT8 *out = (UINT8 *)imOut->image[y];
143

144
                out[0] = in0[0];
145
                for (x = 1; x < im->xsize - 1; x++) {
146
                    float ss = offset;
147
                    ss += KERNEL1x3(in1, x, &kernel[0], 1);
148
                    ss += KERNEL1x3(in0, x, &kernel[3], 1);
149
                    ss += KERNEL1x3(in_1, x, &kernel[6], 1);
150
                    out[x] = clip8(ss);
151
                }
152
                out[x] = in0[x];
153
            }
154
        }
155
    } else {
156
        // Add one time for rounding
157
        offset += 0.5;
158
        for (y = 1; y < im->ysize - 1; y++) {
159
            UINT8 *in_1 = (UINT8 *)im->image[y - 1];
160
            UINT8 *in0 = (UINT8 *)im->image[y];
161
            UINT8 *in1 = (UINT8 *)im->image[y + 1];
162
            UINT8 *out = (UINT8 *)imOut->image[y];
163

164
            memcpy(out, in0, sizeof(UINT32));
165
            if (im->bands == 2) {
166
                for (x = 1; x < im->xsize - 1; x++) {
167
                    float ss0 = offset;
168
                    float ss3 = offset;
169
                    UINT32 v;
170
                    ss0 += KERNEL1x3(in1, x * 4 + 0, &kernel[0], 4);
171
                    ss3 += KERNEL1x3(in1, x * 4 + 3, &kernel[0], 4);
172
                    ss0 += KERNEL1x3(in0, x * 4 + 0, &kernel[3], 4);
173
                    ss3 += KERNEL1x3(in0, x * 4 + 3, &kernel[3], 4);
174
                    ss0 += KERNEL1x3(in_1, x * 4 + 0, &kernel[6], 4);
175
                    ss3 += KERNEL1x3(in_1, x * 4 + 3, &kernel[6], 4);
176
                    v = MAKE_UINT32(clip8(ss0), 0, 0, clip8(ss3));
177
                    memcpy(out + x * sizeof(v), &v, sizeof(v));
178
                }
179
            } else if (im->bands == 3) {
180
                for (x = 1; x < im->xsize - 1; x++) {
181
                    float ss0 = offset;
182
                    float ss1 = offset;
183
                    float ss2 = offset;
184
                    UINT32 v;
185
                    ss0 += KERNEL1x3(in1, x * 4 + 0, &kernel[0], 4);
186
                    ss1 += KERNEL1x3(in1, x * 4 + 1, &kernel[0], 4);
187
                    ss2 += KERNEL1x3(in1, x * 4 + 2, &kernel[0], 4);
188
                    ss0 += KERNEL1x3(in0, x * 4 + 0, &kernel[3], 4);
189
                    ss1 += KERNEL1x3(in0, x * 4 + 1, &kernel[3], 4);
190
                    ss2 += KERNEL1x3(in0, x * 4 + 2, &kernel[3], 4);
191
                    ss0 += KERNEL1x3(in_1, x * 4 + 0, &kernel[6], 4);
192
                    ss1 += KERNEL1x3(in_1, x * 4 + 1, &kernel[6], 4);
193
                    ss2 += KERNEL1x3(in_1, x * 4 + 2, &kernel[6], 4);
194
                    v = MAKE_UINT32(clip8(ss0), clip8(ss1), clip8(ss2), 0);
195
                    memcpy(out + x * sizeof(v), &v, sizeof(v));
196
                }
197
            } else if (im->bands == 4) {
198
                for (x = 1; x < im->xsize - 1; x++) {
199
                    float ss0 = offset;
200
                    float ss1 = offset;
201
                    float ss2 = offset;
202
                    float ss3 = offset;
203
                    UINT32 v;
204
                    ss0 += KERNEL1x3(in1, x * 4 + 0, &kernel[0], 4);
205
                    ss1 += KERNEL1x3(in1, x * 4 + 1, &kernel[0], 4);
206
                    ss2 += KERNEL1x3(in1, x * 4 + 2, &kernel[0], 4);
207
                    ss3 += KERNEL1x3(in1, x * 4 + 3, &kernel[0], 4);
208
                    ss0 += KERNEL1x3(in0, x * 4 + 0, &kernel[3], 4);
209
                    ss1 += KERNEL1x3(in0, x * 4 + 1, &kernel[3], 4);
210
                    ss2 += KERNEL1x3(in0, x * 4 + 2, &kernel[3], 4);
211
                    ss3 += KERNEL1x3(in0, x * 4 + 3, &kernel[3], 4);
212
                    ss0 += KERNEL1x3(in_1, x * 4 + 0, &kernel[6], 4);
213
                    ss1 += KERNEL1x3(in_1, x * 4 + 1, &kernel[6], 4);
214
                    ss2 += KERNEL1x3(in_1, x * 4 + 2, &kernel[6], 4);
215
                    ss3 += KERNEL1x3(in_1, x * 4 + 3, &kernel[6], 4);
216
                    v = MAKE_UINT32(clip8(ss0), clip8(ss1), clip8(ss2), clip8(ss3));
217
                    memcpy(out + x * sizeof(v), &v, sizeof(v));
218
                }
219
            }
220
            memcpy(out + x * sizeof(UINT32), in0 + x * sizeof(UINT32), sizeof(UINT32));
221
        }
222
    }
223
    memcpy(imOut->image[y], im->image[y], im->linesize);
224
}
225

226
void
227
ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
228
#define KERNEL1x5(in0, x, kernel, d)                                       \
229
    (_i2f(in0[x - d - d]) * (kernel)[0] + _i2f(in0[x - d]) * (kernel)[1] + \
230
     _i2f(in0[x]) * (kernel)[2] + _i2f(in0[x + d]) * (kernel)[3] +         \
231
     _i2f(in0[x + d + d]) * (kernel)[4])
232

233
    int x = 0, y = 0;
234

235
    memcpy(imOut->image[0], im->image[0], im->linesize);
236
    memcpy(imOut->image[1], im->image[1], im->linesize);
237
    if (im->bands == 1) {
238
        // Add one time for rounding
239
        offset += 0.5;
240
        if (im->type == IMAGING_TYPE_INT32) {
241
            for (y = 2; y < im->ysize - 2; y++) {
242
                INT32 *in_2 = (INT32 *)im->image[y - 2];
243
                INT32 *in_1 = (INT32 *)im->image[y - 1];
244
                INT32 *in0 = (INT32 *)im->image[y];
245
                INT32 *in1 = (INT32 *)im->image[y + 1];
246
                INT32 *in2 = (INT32 *)im->image[y + 2];
247
                INT32 *out = (INT32 *)imOut->image[y];
248

249
                out[0] = in0[0];
250
                out[1] = in0[1];
251
                for (x = 2; x < im->xsize - 2; x++) {
252
                    float ss = offset;
253
                    ss += KERNEL1x5(in2, x, &kernel[0], 1);
254
                    ss += KERNEL1x5(in1, x, &kernel[5], 1);
255
                    ss += KERNEL1x5(in0, x, &kernel[10], 1);
256
                    ss += KERNEL1x5(in_1, x, &kernel[15], 1);
257
                    ss += KERNEL1x5(in_2, x, &kernel[20], 1);
258
                    out[x] = clip32(ss);
259
                }
260
                out[x + 0] = in0[x + 0];
261
                out[x + 1] = in0[x + 1];
262
            }
263
        } else {
264
            for (y = 2; y < im->ysize - 2; y++) {
265
                UINT8 *in_2 = (UINT8 *)im->image[y - 2];
266
                UINT8 *in_1 = (UINT8 *)im->image[y - 1];
267
                UINT8 *in0 = (UINT8 *)im->image[y];
268
                UINT8 *in1 = (UINT8 *)im->image[y + 1];
269
                UINT8 *in2 = (UINT8 *)im->image[y + 2];
270
                UINT8 *out = (UINT8 *)imOut->image[y];
271

272
                out[0] = in0[0];
273
                out[1] = in0[1];
274
                for (x = 2; x < im->xsize - 2; x++) {
275
                    float ss = offset;
276
                    ss += KERNEL1x5(in2, x, &kernel[0], 1);
277
                    ss += KERNEL1x5(in1, x, &kernel[5], 1);
278
                    ss += KERNEL1x5(in0, x, &kernel[10], 1);
279
                    ss += KERNEL1x5(in_1, x, &kernel[15], 1);
280
                    ss += KERNEL1x5(in_2, x, &kernel[20], 1);
281
                    out[x] = clip8(ss);
282
                }
283
                out[x + 0] = in0[x + 0];
284
                out[x + 1] = in0[x + 1];
285
            }
286
        }
287
    } else {
288
        // Add one time for rounding
289
        offset += 0.5;
290
        for (y = 2; y < im->ysize - 2; y++) {
291
            UINT8 *in_2 = (UINT8 *)im->image[y - 2];
292
            UINT8 *in_1 = (UINT8 *)im->image[y - 1];
293
            UINT8 *in0 = (UINT8 *)im->image[y];
294
            UINT8 *in1 = (UINT8 *)im->image[y + 1];
295
            UINT8 *in2 = (UINT8 *)im->image[y + 2];
296
            UINT8 *out = (UINT8 *)imOut->image[y];
297

298
            memcpy(out, in0, sizeof(UINT32) * 2);
299
            if (im->bands == 2) {
300
                for (x = 2; x < im->xsize - 2; x++) {
301
                    float ss0 = offset;
302
                    float ss3 = offset;
303
                    UINT32 v;
304
                    ss0 += KERNEL1x5(in2, x * 4 + 0, &kernel[0], 4);
305
                    ss3 += KERNEL1x5(in2, x * 4 + 3, &kernel[0], 4);
306
                    ss0 += KERNEL1x5(in1, x * 4 + 0, &kernel[5], 4);
307
                    ss3 += KERNEL1x5(in1, x * 4 + 3, &kernel[5], 4);
308
                    ss0 += KERNEL1x5(in0, x * 4 + 0, &kernel[10], 4);
309
                    ss3 += KERNEL1x5(in0, x * 4 + 3, &kernel[10], 4);
310
                    ss0 += KERNEL1x5(in_1, x * 4 + 0, &kernel[15], 4);
311
                    ss3 += KERNEL1x5(in_1, x * 4 + 3, &kernel[15], 4);
312
                    ss0 += KERNEL1x5(in_2, x * 4 + 0, &kernel[20], 4);
313
                    ss3 += KERNEL1x5(in_2, x * 4 + 3, &kernel[20], 4);
314
                    v = MAKE_UINT32(clip8(ss0), 0, 0, clip8(ss3));
315
                    memcpy(out + x * sizeof(v), &v, sizeof(v));
316
                }
317
            } else if (im->bands == 3) {
318
                for (x = 2; x < im->xsize - 2; x++) {
319
                    float ss0 = offset;
320
                    float ss1 = offset;
321
                    float ss2 = offset;
322
                    UINT32 v;
323
                    ss0 += KERNEL1x5(in2, x * 4 + 0, &kernel[0], 4);
324
                    ss1 += KERNEL1x5(in2, x * 4 + 1, &kernel[0], 4);
325
                    ss2 += KERNEL1x5(in2, x * 4 + 2, &kernel[0], 4);
326
                    ss0 += KERNEL1x5(in1, x * 4 + 0, &kernel[5], 4);
327
                    ss1 += KERNEL1x5(in1, x * 4 + 1, &kernel[5], 4);
328
                    ss2 += KERNEL1x5(in1, x * 4 + 2, &kernel[5], 4);
329
                    ss0 += KERNEL1x5(in0, x * 4 + 0, &kernel[10], 4);
330
                    ss1 += KERNEL1x5(in0, x * 4 + 1, &kernel[10], 4);
331
                    ss2 += KERNEL1x5(in0, x * 4 + 2, &kernel[10], 4);
332
                    ss0 += KERNEL1x5(in_1, x * 4 + 0, &kernel[15], 4);
333
                    ss1 += KERNEL1x5(in_1, x * 4 + 1, &kernel[15], 4);
334
                    ss2 += KERNEL1x5(in_1, x * 4 + 2, &kernel[15], 4);
335
                    ss0 += KERNEL1x5(in_2, x * 4 + 0, &kernel[20], 4);
336
                    ss1 += KERNEL1x5(in_2, x * 4 + 1, &kernel[20], 4);
337
                    ss2 += KERNEL1x5(in_2, x * 4 + 2, &kernel[20], 4);
338
                    v = MAKE_UINT32(clip8(ss0), clip8(ss1), clip8(ss2), 0);
339
                    memcpy(out + x * sizeof(v), &v, sizeof(v));
340
                }
341
            } else if (im->bands == 4) {
342
                for (x = 2; x < im->xsize - 2; x++) {
343
                    float ss0 = offset;
344
                    float ss1 = offset;
345
                    float ss2 = offset;
346
                    float ss3 = offset;
347
                    UINT32 v;
348
                    ss0 += KERNEL1x5(in2, x * 4 + 0, &kernel[0], 4);
349
                    ss1 += KERNEL1x5(in2, x * 4 + 1, &kernel[0], 4);
350
                    ss2 += KERNEL1x5(in2, x * 4 + 2, &kernel[0], 4);
351
                    ss3 += KERNEL1x5(in2, x * 4 + 3, &kernel[0], 4);
352
                    ss0 += KERNEL1x5(in1, x * 4 + 0, &kernel[5], 4);
353
                    ss1 += KERNEL1x5(in1, x * 4 + 1, &kernel[5], 4);
354
                    ss2 += KERNEL1x5(in1, x * 4 + 2, &kernel[5], 4);
355
                    ss3 += KERNEL1x5(in1, x * 4 + 3, &kernel[5], 4);
356
                    ss0 += KERNEL1x5(in0, x * 4 + 0, &kernel[10], 4);
357
                    ss1 += KERNEL1x5(in0, x * 4 + 1, &kernel[10], 4);
358
                    ss2 += KERNEL1x5(in0, x * 4 + 2, &kernel[10], 4);
359
                    ss3 += KERNEL1x5(in0, x * 4 + 3, &kernel[10], 4);
360
                    ss0 += KERNEL1x5(in_1, x * 4 + 0, &kernel[15], 4);
361
                    ss1 += KERNEL1x5(in_1, x * 4 + 1, &kernel[15], 4);
362
                    ss2 += KERNEL1x5(in_1, x * 4 + 2, &kernel[15], 4);
363
                    ss3 += KERNEL1x5(in_1, x * 4 + 3, &kernel[15], 4);
364
                    ss0 += KERNEL1x5(in_2, x * 4 + 0, &kernel[20], 4);
365
                    ss1 += KERNEL1x5(in_2, x * 4 + 1, &kernel[20], 4);
366
                    ss2 += KERNEL1x5(in_2, x * 4 + 2, &kernel[20], 4);
367
                    ss3 += KERNEL1x5(in_2, x * 4 + 3, &kernel[20], 4);
368
                    v = MAKE_UINT32(clip8(ss0), clip8(ss1), clip8(ss2), clip8(ss3));
369
                    memcpy(out + x * sizeof(v), &v, sizeof(v));
370
                }
371
            }
372
            memcpy(
373
                out + x * sizeof(UINT32), in0 + x * sizeof(UINT32), sizeof(UINT32) * 2
374
            );
375
        }
376
    }
377
    memcpy(imOut->image[y], im->image[y], im->linesize);
378
    memcpy(imOut->image[y + 1], im->image[y + 1], im->linesize);
379
}
380

381
Imaging
382
ImagingFilter(Imaging im, int xsize, int ysize, const FLOAT32 *kernel, FLOAT32 offset) {
383
    Imaging imOut;
384
    ImagingSectionCookie cookie;
385

386
    if (im->type != IMAGING_TYPE_UINT8 && im->type != IMAGING_TYPE_INT32) {
387
        return (Imaging)ImagingError_ModeError();
388
    }
389

390
    if (im->xsize < xsize || im->ysize < ysize) {
391
        return ImagingCopy(im);
392
    }
393

394
    if ((xsize != 3 && xsize != 5) || xsize != ysize) {
395
        return (Imaging)ImagingError_ValueError("bad kernel size");
396
    }
397

398
    imOut = ImagingNewDirty(im->mode, im->xsize, im->ysize);
399
    if (!imOut) {
400
        return NULL;
401
    }
402

403
    ImagingSectionEnter(&cookie);
404
    if (xsize == 3) {
405
        /* 3x3 kernel. */
406
        ImagingFilter3x3(imOut, im, kernel, offset);
407
    } else {
408
        /* 5x5 kernel. */
409
        ImagingFilter5x5(imOut, im, kernel, offset);
410
    }
411
    ImagingSectionLeave(&cookie);
412
    return imOut;
413
}
414