2424 * FIXME: Implement image processing gradient filters
2525 */
2626
27+ #include <stdint.h>
2728#include "Imaging.h"
2829
2930#define ROUND_UP (f ) ((int)((f) >= 0.0 ? (f) + 0.5F : (f) - 0.5F))
31+ #define INT32_MAX_F 2147483647.0F
3032
3133static inline UINT8
3234clip8 (float in ) {
33- if (in <= 0.0 ) {
34- return 0 ;
35- }
36- if (in >= 255.0 ) {
37- return 255 ;
38- }
35+ // Branchless clamp to [0, 255].
36+ in = in < 0.0f ? 0.0f : in ;
37+ in = in > 255.0f ? 255.0f : in ;
3938 return (UINT8 )in ;
4039}
4140
4241static inline INT32
4342clip32 (float in ) {
44- if (in <= 0.0 ) {
45- return 0 ;
46- }
47- if (in >= pow (2 , 31 ) - 1 ) {
48- return pow (2 , 31 ) - 1 ;
49- }
50- return (INT32 )in ;
43+ // Clamp the low bound branchlessly (maxss).
44+ in = in < 0.0f ? 0.0f : in ;
45+ // Must explicitly return INT32_MAX for the high bound
46+ // to avoid incorrect overflow to INT_MIN.
47+ return in >= INT32_MAX_F ? INT32_MAX : (INT32 )in ;
5148}
5249
5350Imaging
@@ -116,36 +113,39 @@ kernel_i16(int size, UINT8 *in0, int x, const float *kernel, int bigendian) {
116113 int half_size = (size - 1 ) / 2 ;
117114 for (i = 0 ; i < size ; i ++ ) {
118115 int x1 = x + i - half_size ;
119- result += _i2f (
120- in0 [x1 * 2 + (bigendian ? 1 : 0 )] +
121- (in0 [x1 * 2 + (bigendian ? 0 : 1 )] >> 8 )
122- ) *
116+ result += (float )(in0 [x1 * 2 + (bigendian ? 1 : 0 )] +
117+ (in0 [x1 * 2 + (bigendian ? 0 : 1 )] >> 8 )) *
123118 kernel [i ];
124119 }
125120 return result ;
126121}
127122
128- void
123+ static void
129124ImagingFilter3x3 (Imaging imOut , Imaging im , const float * kernel , float offset ) {
130- #define KERNEL1x3 (in0 , x , kernel , d ) \
131- (_i2f( in0[x - d]) * (kernel)[0] + _i2f (in0[x]) * (kernel)[1] + \
132- _i2f (in0[x + d]) * (kernel)[2])
125+ #define KERNEL1x3 (in0 , x , kernel , d ) \
126+ ((float)( in0[x - d]) * (kernel)[0] + (float) (in0[x]) * (kernel)[1] + \
127+ (float) (in0[x + d]) * (kernel)[2])
133128
134129 int x = 0 , y = 0 ;
130+ // restrict safety: assert im and imOut to be separate,
131+ // so we can safely `restrict` the copy pointers below
132+ assert (im -> image != imOut -> image );
133+ // Hoist for loop optimization.
134+ int xsize = im -> xsize , ysize = im -> ysize ;
135135
136136 memcpy (imOut -> image [0 ], im -> image [0 ], im -> linesize );
137137 if (im -> bands == 1 ) {
138138 // Add one time for rounding
139139 offset += 0.5 ;
140140 if (im -> type == IMAGING_TYPE_INT32 ) {
141- for (y = 1 ; y < im -> ysize - 1 ; y ++ ) {
142- INT32 * in_1 = (INT32 * )im -> image [y - 1 ];
143- INT32 * in0 = (INT32 * )im -> image [y ];
144- INT32 * in1 = (INT32 * )im -> image [y + 1 ];
145- INT32 * out = (INT32 * )imOut -> image [y ];
141+ for (y = 1 ; y < ysize - 1 ; y ++ ) {
142+ INT32 * restrict in_1 = (INT32 * )im -> image [y - 1 ];
143+ INT32 * restrict in0 = (INT32 * )im -> image [y ];
144+ INT32 * restrict in1 = (INT32 * )im -> image [y + 1 ];
145+ INT32 * restrict out = (INT32 * )imOut -> image [y ];
146146
147147 out [0 ] = in0 [0 ];
148- for (x = 1 ; x < im -> xsize - 1 ; x ++ ) {
148+ for (x = 1 ; x < xsize - 1 ; x ++ ) {
149149 float ss = offset ;
150150 ss += KERNEL1x3 (in1 , x , & kernel [0 ], 1 );
151151 ss += KERNEL1x3 (in0 , x , & kernel [3 ], 1 );
@@ -166,17 +166,17 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
166166 bigendian = 1 ;
167167 }
168168 }
169- for (y = 1 ; y < im -> ysize - 1 ; y ++ ) {
170- UINT8 * in_1 = (UINT8 * )im -> image [y - 1 ];
171- UINT8 * in0 = (UINT8 * )im -> image [y ];
172- UINT8 * in1 = (UINT8 * )im -> image [y + 1 ];
173- UINT8 * out = (UINT8 * )imOut -> image [y ];
169+ for (y = 1 ; y < ysize - 1 ; y ++ ) {
170+ UINT8 * restrict in_1 = (UINT8 * )im -> image [y - 1 ];
171+ UINT8 * restrict in0 = (UINT8 * )im -> image [y ];
172+ UINT8 * restrict in1 = (UINT8 * )im -> image [y + 1 ];
173+ UINT8 * restrict out = (UINT8 * )imOut -> image [y ];
174174
175175 out [0 ] = in0 [0 ];
176176 if (im -> type == IMAGING_TYPE_SPECIAL ) {
177177 out [1 ] = in0 [1 ];
178178 }
179- for (x = 1 ; x < im -> xsize - 1 ; x ++ ) {
179+ for (x = 1 ; x < xsize - 1 ; x ++ ) {
180180 float ss = offset ;
181181 if (im -> type == IMAGING_TYPE_SPECIAL ) {
182182 ss += kernel_i16 (3 , in1 , x , & kernel [0 ], bigendian );
@@ -203,15 +203,15 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
203203 } else {
204204 // Add one time for rounding
205205 offset += 0.5 ;
206- for (y = 1 ; y < im -> ysize - 1 ; y ++ ) {
207- UINT8 * in_1 = (UINT8 * )im -> image [y - 1 ];
208- UINT8 * in0 = (UINT8 * )im -> image [y ];
209- UINT8 * in1 = (UINT8 * )im -> image [y + 1 ];
210- UINT8 * out = (UINT8 * )imOut -> image [y ];
206+ for (y = 1 ; y < ysize - 1 ; y ++ ) {
207+ UINT8 * restrict in_1 = (UINT8 * )im -> image [y - 1 ];
208+ UINT8 * restrict in0 = (UINT8 * )im -> image [y ];
209+ UINT8 * restrict in1 = (UINT8 * )im -> image [y + 1 ];
210+ UINT8 * restrict out = (UINT8 * )imOut -> image [y ];
211211
212212 memcpy (out , in0 , sizeof (UINT32 ));
213213 if (im -> bands == 2 ) {
214- for (x = 1 ; x < im -> xsize - 1 ; x ++ ) {
214+ for (x = 1 ; x < xsize - 1 ; x ++ ) {
215215 float ss0 = offset ;
216216 float ss3 = offset ;
217217 UINT32 v ;
@@ -225,7 +225,7 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
225225 memcpy (out + x * sizeof (v ), & v , sizeof (v ));
226226 }
227227 } else if (im -> bands == 3 ) {
228- for (x = 1 ; x < im -> xsize - 1 ; x ++ ) {
228+ for (x = 1 ; x < xsize - 1 ; x ++ ) {
229229 float ss0 = offset ;
230230 float ss1 = offset ;
231231 float ss2 = offset ;
@@ -243,7 +243,7 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
243243 memcpy (out + x * sizeof (v ), & v , sizeof (v ));
244244 }
245245 } else if (im -> bands == 4 ) {
246- for (x = 1 ; x < im -> xsize - 1 ; x ++ ) {
246+ for (x = 1 ; x < xsize - 1 ; x ++ ) {
247247 float ss0 = offset ;
248248 float ss1 = offset ;
249249 float ss2 = offset ;
@@ -271,32 +271,38 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
271271 memcpy (imOut -> image [y ], im -> image [y ], im -> linesize );
272272}
273273
274- void
274+ static void
275275ImagingFilter5x5 (Imaging imOut , Imaging im , const float * kernel , float offset ) {
276- #define KERNEL1x5 (in0 , x , kernel , d ) \
277- (_i2f( in0[x - d - d]) * (kernel)[0] + _i2f (in0[x - d]) * (kernel)[1] + \
278- _i2f( in0[x]) * (kernel)[2] + _i2f (in0[x + d]) * (kernel)[3] + \
279- _i2f (in0[x + d + d]) * (kernel)[4])
276+ #define KERNEL1x5 (in0 , x , kernel , d ) \
277+ ((float)( in0[x - d - d]) * (kernel)[0] + (float) (in0[x - d]) * (kernel)[1] + \
278+ (float)( in0[x]) * (kernel)[2] + (float) (in0[x + d]) * (kernel)[3] + \
279+ (float) (in0[x + d + d]) * (kernel)[4])
280280
281281 int x = 0 , y = 0 ;
282282
283+ // restrict safety: assert im and imOut to be separate,
284+ // so we can safely `restrict` the copy pointers below
285+ assert (im -> image != imOut -> image );
286+ // Hoist for loop optimization.
287+ int xsize = im -> xsize , ysize = im -> ysize ;
288+
283289 memcpy (imOut -> image [0 ], im -> image [0 ], im -> linesize );
284290 memcpy (imOut -> image [1 ], im -> image [1 ], im -> linesize );
285291 if (im -> bands == 1 ) {
286292 // Add one time for rounding
287293 offset += 0.5 ;
288294 if (im -> type == IMAGING_TYPE_INT32 ) {
289- for (y = 2 ; y < im -> ysize - 2 ; y ++ ) {
290- INT32 * in_2 = (INT32 * )im -> image [y - 2 ];
291- INT32 * in_1 = (INT32 * )im -> image [y - 1 ];
292- INT32 * in0 = (INT32 * )im -> image [y ];
293- INT32 * in1 = (INT32 * )im -> image [y + 1 ];
294- INT32 * in2 = (INT32 * )im -> image [y + 2 ];
295- INT32 * out = (INT32 * )imOut -> image [y ];
295+ for (y = 2 ; y < ysize - 2 ; y ++ ) {
296+ INT32 * restrict in_2 = (INT32 * )im -> image [y - 2 ];
297+ INT32 * restrict in_1 = (INT32 * )im -> image [y - 1 ];
298+ INT32 * restrict in0 = (INT32 * )im -> image [y ];
299+ INT32 * restrict in1 = (INT32 * )im -> image [y + 1 ];
300+ INT32 * restrict in2 = (INT32 * )im -> image [y + 2 ];
301+ INT32 * restrict out = (INT32 * )imOut -> image [y ];
296302
297303 out [0 ] = in0 [0 ];
298304 out [1 ] = in0 [1 ];
299- for (x = 2 ; x < im -> xsize - 2 ; x ++ ) {
305+ for (x = 2 ; x < xsize - 2 ; x ++ ) {
300306 float ss = offset ;
301307 ss += KERNEL1x5 (in2 , x , & kernel [0 ], 1 );
302308 ss += KERNEL1x5 (in1 , x , & kernel [5 ], 1 );
@@ -320,21 +326,21 @@ ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
320326 bigendian = 1 ;
321327 }
322328 }
323- for (y = 2 ; y < im -> ysize - 2 ; y ++ ) {
324- UINT8 * in_2 = (UINT8 * )im -> image [y - 2 ];
325- UINT8 * in_1 = (UINT8 * )im -> image [y - 1 ];
326- UINT8 * in0 = (UINT8 * )im -> image [y ];
327- UINT8 * in1 = (UINT8 * )im -> image [y + 1 ];
328- UINT8 * in2 = (UINT8 * )im -> image [y + 2 ];
329- UINT8 * out = (UINT8 * )imOut -> image [y ];
329+ for (y = 2 ; y < ysize - 2 ; y ++ ) {
330+ UINT8 * restrict in_2 = (UINT8 * )im -> image [y - 2 ];
331+ UINT8 * restrict in_1 = (UINT8 * )im -> image [y - 1 ];
332+ UINT8 * restrict in0 = (UINT8 * )im -> image [y ];
333+ UINT8 * restrict in1 = (UINT8 * )im -> image [y + 1 ];
334+ UINT8 * restrict in2 = (UINT8 * )im -> image [y + 2 ];
335+ UINT8 * restrict out = (UINT8 * )imOut -> image [y ];
330336
331337 out [0 ] = in0 [0 ];
332338 out [1 ] = in0 [1 ];
333339 if (im -> type == IMAGING_TYPE_SPECIAL ) {
334340 out [2 ] = in0 [2 ];
335341 out [3 ] = in0 [3 ];
336342 }
337- for (x = 2 ; x < im -> xsize - 2 ; x ++ ) {
343+ for (x = 2 ; x < xsize - 2 ; x ++ ) {
338344 float ss = offset ;
339345 if (im -> type == IMAGING_TYPE_SPECIAL ) {
340346 ss += kernel_i16 (5 , in2 , x , & kernel [0 ], bigendian );
@@ -368,17 +374,17 @@ ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
368374 } else {
369375 // Add one time for rounding
370376 offset += 0.5 ;
371- for (y = 2 ; y < im -> ysize - 2 ; y ++ ) {
372- UINT8 * in_2 = (UINT8 * )im -> image [y - 2 ];
373- UINT8 * in_1 = (UINT8 * )im -> image [y - 1 ];
374- UINT8 * in0 = (UINT8 * )im -> image [y ];
375- UINT8 * in1 = (UINT8 * )im -> image [y + 1 ];
376- UINT8 * in2 = (UINT8 * )im -> image [y + 2 ];
377- UINT8 * out = (UINT8 * )imOut -> image [y ];
377+ for (y = 2 ; y < ysize - 2 ; y ++ ) {
378+ UINT8 * restrict in_2 = (UINT8 * )im -> image [y - 2 ];
379+ UINT8 * restrict in_1 = (UINT8 * )im -> image [y - 1 ];
380+ UINT8 * restrict in0 = (UINT8 * )im -> image [y ];
381+ UINT8 * restrict in1 = (UINT8 * )im -> image [y + 1 ];
382+ UINT8 * restrict in2 = (UINT8 * )im -> image [y + 2 ];
383+ UINT8 * restrict out = (UINT8 * )imOut -> image [y ];
378384
379385 memcpy (out , in0 , sizeof (UINT32 ) * 2 );
380386 if (im -> bands == 2 ) {
381- for (x = 2 ; x < im -> xsize - 2 ; x ++ ) {
387+ for (x = 2 ; x < xsize - 2 ; x ++ ) {
382388 float ss0 = offset ;
383389 float ss3 = offset ;
384390 UINT32 v ;
@@ -396,7 +402,7 @@ ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
396402 memcpy (out + x * sizeof (v ), & v , sizeof (v ));
397403 }
398404 } else if (im -> bands == 3 ) {
399- for (x = 2 ; x < im -> xsize - 2 ; x ++ ) {
405+ for (x = 2 ; x < xsize - 2 ; x ++ ) {
400406 float ss0 = offset ;
401407 float ss1 = offset ;
402408 float ss2 = offset ;
@@ -420,7 +426,7 @@ ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
420426 memcpy (out + x * sizeof (v ), & v , sizeof (v ));
421427 }
422428 } else if (im -> bands == 4 ) {
423- for (x = 2 ; x < im -> xsize - 2 ; x ++ ) {
429+ for (x = 2 ; x < xsize - 2 ; x ++ ) {
424430 float ss0 = offset ;
425431 float ss1 = offset ;
426432 float ss2 = offset ;
0 commit comments