faster box blur for radius < 1

src/libImaging/BoxBlur.c

@@ -19,7 +19,7 @@ ImagingLineBoxBlur32(UINT32 *lineOut, UINT32 *lineIn, int lastx, int radius, int
     __m128i wws = _mm_set1_epi32(ww);
     __m128i fws = _mm_set1_epi32(fw);
     __m128i b23 = _mm_set1_epi32(1 << 23);
-    __m128i shiftmask = _mm_set_epi8(0,0,0,0,0,0,0,0,0,0,0,0,3,7,11,15);    
+    __m128i shiftmask = _mm_set_epi8(0,0,0,0,0,0,0,0,0,0,0,0,3,7,11,15);
     
 
     // #define MOVE_ACC(acc, subtract, add)
@@ -123,6 +123,66 @@ ImagingLineBoxBlur32(UINT32 *lineOut, UINT32 *lineIn, int lastx, int radius, int
 }
 
 
+void static inline
+ImagingLineBoxBlurZero32(UINT32 *lineOut, UINT32 *lineIn, int lastx, int edgeA,
+    int edgeB, UINT32 ww, UINT32 fw)
+{
+    int x;
+    __m128i acc;
+    __m128i bulk;
+    __m128i wws = _mm_set1_epi32(ww);
+    __m128i fws = _mm_set1_epi32(fw);
+    __m128i b23 = _mm_set1_epi32(1 << 23);
+    __m128i shiftmask = _mm_set_epi8(0,0,0,0,0,0,0,0,0,0,0,0,3,7,11,15);
+
+    // #define ADD_FAR(bulk, acc, left, right)
+    //     bulk[0] = (acc[0] * ww) + (lineIn[left][0] + lineIn[right][0]) * fw;
+
+    #define ADD_FAR(bulk, acc, left, right) \
+        { \
+        __m128i tmp3 = _mm_cvtepu8_epi32(_mm_cvtsi32_si128(lineIn[left])); \
+        __m128i tmp4 = _mm_cvtepu8_epi32(_mm_cvtsi32_si128(lineIn[right])); \
+        __m128i tmp5 = _mm_mullo_epi32(_mm_add_epi32(tmp3, tmp4), fws); \
+        __m128i tmp6 = _mm_mullo_epi32(acc, wws); \
+        bulk = _mm_add_epi32(tmp5, tmp6); \
+        }
+
+    // #define SAVE(x, bulk)
+    //    lineOut[x][0] = (UINT8)((bulk[0] + (1 << 23)) >> 24);
+
+    #define SAVE(x, bulk) \
+        lineOut[x] = _mm_cvtsi128_si32( \
+            _mm_shuffle_epi8(_mm_add_epi32(bulk, b23), shiftmask) \
+        )
+
+
+    /* Subtract pixel from left ("0").
+       Add pixels from radius. */
+    for (x = 0; x < edgeA; x++) {
+        acc = _mm_cvtepu8_epi32(_mm_cvtsi32_si128(lineIn[0]));
+        ADD_FAR(bulk, acc, 0, x + 1);
+        SAVE(x, bulk);
+    }
+    /* Subtract previous pixel from "-radius".
+       Add pixels from radius. */
+    for (x = edgeA; x < edgeB; x++) {
+        acc = _mm_cvtepu8_epi32(_mm_cvtsi32_si128(lineIn[x]));
+        ADD_FAR(bulk, acc, x - 1, x + 1);
+        SAVE(x, bulk);
+    }
+    /* Subtract previous pixel from "-radius".
+       Add last pixel. */
+    for (x = edgeB; x <= lastx; x++) {
+        acc = _mm_cvtepu8_epi32(_mm_cvtsi32_si128(lineIn[lastx]));
+        ADD_FAR(bulk, acc, x - 1, lastx);
+        SAVE(x, bulk);
+    }
+
+    #undef ADD_FAR
+    #undef SAVE
+}
+
+
 void static inline
 ImagingLineBoxBlur8(UINT8 *lineOut, UINT8 *lineIn, int lastx, int radius, int edgeA,
     int edgeB, UINT32 ww, UINT32 fw)
@@ -189,6 +249,45 @@ ImagingLineBoxBlur8(UINT8 *lineOut, UINT8 *lineIn, int lastx, int radius, int ed
 }
 
 
+void static inline
+ImagingLineBoxBlurZero8(UINT8 *lineOut, UINT8 *lineIn, int lastx, int edgeA,
+    int edgeB, UINT32 ww, UINT32 fw)
+{
+    int x;
+    UINT32 acc;
+    UINT32 bulk;
+
+    #define MOVE_ACC(acc, subtract, add) \
+        acc += lineIn[add] - lineIn[subtract];
+
+    #define ADD_FAR(bulk, acc, left, right) \
+        bulk = (acc * ww) + (lineIn[left] + lineIn[right]) * fw;
+
+    #define SAVE(x, bulk) \
+        lineOut[x] = (UINT8)((bulk + (1 << 23)) >> 24)
+
+    for (x = 0; x < edgeA; x++) {
+        acc = lineIn[0];
+        ADD_FAR(bulk, acc, 0, x + 1);
+        SAVE(x, bulk);
+    }
+    for (x = edgeA; x < edgeB; x++) {
+        acc = lineIn[x];
+        ADD_FAR(bulk, acc, x - 1, x + 1);
+        SAVE(x, bulk);
+    }
+    for (x = edgeB; x <= lastx; x++) {
+        acc = lineIn[lastx];
+        ADD_FAR(bulk, acc, x - 1, lastx);
+        SAVE(x, bulk);
+    }
+
+    #undef MOVE_ACC
+    #undef ADD_FAR
+    #undef SAVE
+}
+
+
 
 Imaging
 ImagingHorizontalBoxBlur(Imaging imOut, Imaging imIn, float floatRadius)
@@ -208,39 +307,72 @@ ImagingHorizontalBoxBlur(Imaging imOut, Imaging imIn, float floatRadius)
     if (lineOut == NULL)
         return ImagingError_MemoryError();
 
-    // printf(">>> %d %d %d\n", radius, ww, fw);
+    // printf(">>> radius: %d edges: %d %d, weights: %d %d\n",
+    //     radius, edgeA, edgeB, ww, fw);
 
     ImagingSectionEnter(&cookie);
 
     if (imIn->image8)
     {
-        for (y = 0; y < imIn->ysize; y++) {
+        if (radius) {
-            ImagingLineBoxBlur8(
+            for (y = 0; y < imIn->ysize; y++) {
-                (imIn == imOut ? (UINT8 *) lineOut : imOut->image8[y]),
+                ImagingLineBoxBlur8(
-                imIn->image8[y],
+                    (imIn == imOut ? (UINT8 *) lineOut : imOut->image8[y]),
-                imIn->xsize - 1,
+                    imIn->image8[y],
-                radius, edgeA, edgeB,
+                    imIn->xsize - 1,
-                ww, fw
+                    radius, edgeA, edgeB,
-            );
+                    ww, fw
-            if (imIn == imOut) {
+                );
-                // Commit.
+                if (imIn == imOut) {
-                memcpy(imOut->image8[y], lineOut, imIn->xsize);
+                    // Commit.
+                    memcpy(imOut->image8[y], lineOut, imIn->xsize);
+                }
+            }
+        } else {
+            for (y = 0; y < imIn->ysize; y++) {
+                ImagingLineBoxBlurZero8(
+                    (imIn == imOut ? (UINT8 *) lineOut : imOut->image8[y]),
+                    imIn->image8[y],
+                    imIn->xsize - 1,
+                    edgeA, edgeB,
+                    ww, fw
+                );
+                if (imIn == imOut) {
+                    // Commit.
+                    memcpy(imOut->image8[y], lineOut, imIn->xsize);
+                }
             }
         }
     }
     else
     {
-        for (y = 0; y < imIn->ysize; y++) {
+        if (radius) {
-            ImagingLineBoxBlur32(
+            for (y = 0; y < imIn->ysize; y++) {
-                imIn == imOut ? (UINT32 *) lineOut : (UINT32 *) imOut->image32[y],
+                ImagingLineBoxBlur32(
-                (UINT32 *) imIn->image32[y],
+                    imIn == imOut ? (UINT32 *) lineOut : (UINT32 *) imOut->image32[y],
-                imIn->xsize - 1,
+                    (UINT32 *) imIn->image32[y],
-                radius, edgeA, edgeB,
+                    imIn->xsize - 1,
-                ww, fw
+                    radius, edgeA, edgeB,
-            );
+                    ww, fw
-            if (imIn == imOut) {
+                );
-                // Commit.
+                if (imIn == imOut) {
-                memcpy(imOut->image32[y], lineOut, imIn->xsize * 4);
+                    // Commit.
+                    memcpy(imOut->image32[y], lineOut, imIn->xsize * 4);
+                }
+            }
+        } else{
+            for (y = 0; y < imIn->ysize; y++) {
+                ImagingLineBoxBlurZero32(
+                    imIn == imOut ? (UINT32 *) lineOut : (UINT32 *) imOut->image32[y],
+                    (UINT32 *) imIn->image32[y],
+                    imIn->xsize - 1,
+                    edgeA, edgeB,
+                    ww, fw
+                );
+                if (imIn == imOut) {
+                    // Commit.
+                    memcpy(imOut->image32[y], lineOut, imIn->xsize * 4);
+                }
             }
         }
     }