* median.c: add a median filter. Highly unoptimised for now.

git-svn-id: file:///srv/caca.zoy.org/var/lib/svn/libpipi/trunk@2744 92316355-f0b4-4df1-b90c-862c8a59935f

pipi/Makefile.am

@@ -54,7 +54,8 @@ filter_sources = \
 	filter/autocontrast.c \
 	filter/blur.c \
 	filter/convolution.c filter/convolution_template.h \
-	filter/color.c
+	filter/color.c \
+	filter/median.c
 
 quantize_sources = \
 	quantize/reduce.c

pipi/context.c

@@ -164,6 +164,29 @@ int pipi_command(pipi_context_t *ctx, char const *cmd, ...)
         pipi_free(src);
         ctx->images[ctx->nimages - 1] = dst;
     }
+    else if(!strcmp(cmd, "median"))
+    {
+        pipi_image_t *src, *dst;
+        char const *arg;
+        va_list ap;
+        double w, h;
+
+        if(ctx->nimages < 1)
+            return -1;
+        va_start(ap, cmd);
+        arg = va_arg(ap, char const *);
+        va_end(ap);
+        w = h = atof(arg);
+        arg = strchr(arg, 'x');
+        if(arg)
+            h = atof(arg + 1);
+        src = ctx->images[ctx->nimages - 1];
+        dst = pipi_median_ext(src, w, h);
+        if(dst == NULL)
+            return -1;
+        pipi_free(src);
+        ctx->images[ctx->nimages - 1] = dst;
+    }
     else if(!strcmp(cmd, "geometry"))
     {
         pipi_image_t *src, *dst;

pipi/filter/median.c

@@ -0,0 +1,147 @@
+/*
+ *  libpipi       Proper image processing implementation library
+ *  Copyright (c) 2004-2008 Sam Hocevar <sam@zoy.org>
+ *                All Rights Reserved
+ *
+ *  $Id$
+ *
+ *  This library is free software. It comes without any warranty, to
+ *  the extent permitted by applicable law. You can redistribute it
+ *  and/or modify it under the terms of the Do What The Fuck You Want
+ *  To Public License, Version 2, as published by Sam Hocevar. See
+ *  http://sam.zoy.org/wtfpl/COPYING for more details.
+ */
+
+/*
+ * median.c: median filter functions
+ */
+
+#include "config.h"
+#include "common.h"
+
+#include <stdlib.h>
+#include <stdio.h>
+#include <string.h>
+#include <math.h>
+
+#include "pipi.h"
+#include "pipi_internals.h"
+
+pipi_image_t *pipi_median(pipi_image_t *src, int radius)
+{
+    return pipi_median_ext(src, radius, radius);
+}
+
+/* FIXME: this is in desperate want of optimisation. Here is what could
+ * be done to improve the performance:
+ *  - prefetch the neighbourhood; most neighbours are the same as the
+ *    previous pixels.
+ *  - use a better sort algorithm; bubble sort is ridiculous
+ *  - even better, use state-of-the art median selection, ie. O(3n), for
+ *    most common combinations (9, 25, 49, 81). */
+pipi_image_t *pipi_median_ext(pipi_image_t *src, int rx, int ry)
+{
+    pipi_image_t *dst;
+    pipi_pixels_t *srcp, *dstp;
+    float *srcdata, *dstdata;
+    double *list;
+    int x, y, w, h, i, j, k, size, gray;
+
+    w = src->w;
+    h = src->h;
+    size = (2 * rx + 1) * (2 * ry + 1);
+
+    gray = (src->last_modified == PIPI_PIXELS_Y_F);
+
+    srcp = gray ? pipi_getpixels(src, PIPI_PIXELS_Y_F)
+                : pipi_getpixels(src, PIPI_PIXELS_RGBA_F);
+    srcdata = (float *)srcp->pixels;
+
+    dst = pipi_new(w, h);
+    dstp = gray ? pipi_getpixels(dst, PIPI_PIXELS_Y_F)
+                : pipi_getpixels(dst, PIPI_PIXELS_RGBA_F);
+    dstdata = (float *)dstp->pixels;
+
+    list = malloc(size * (gray ? 1 : 4) * sizeof(double));
+
+    for(y = 0; y < h; y++)
+    {
+        for(x = 0; x < w; x++)
+        {
+            double tmp;
+            double *parser = list;
+
+            /* Make a list of neighbours */
+            for(j = -ry; j <= ry; j++)
+            {
+                int y2 = y + j;
+                if(y2 < 0) y2 = h - 1 - ((-y2 - 1) % h);
+                else if(y2 > 0) y2 = y2 % h;
+
+                for(i = -rx; i <= rx; i++)
+                {
+                    int x2 = x + i;
+                    if(x2 < 0) x2 = w - 1 - ((-x2 - 1) % w);
+                    else if(x2 > 0) x2 = x2 % w;
+
+                    if(gray)
+                    {
+                        *parser++ = srcdata[y2 * w + x2];
+                    }
+                    else
+                    {
+                        *parser++ = srcdata[4 * (y2 * w + x2)];
+                        *parser++ = srcdata[4 * (y2 * w + x2) + 1];
+                        *parser++ = srcdata[4 * (y2 * w + x2) + 2];
+                        *parser++ = srcdata[4 * (y2 * w + x2) + 3];
+                    }
+                }
+            }
+
+            /* Sort the list */
+            for(i = 0; i < size; i++)
+            {
+                for(j = 0; j < i; j++)
+                {
+                    if(gray)
+                    {
+                        if(list[i] < list[j])
+                        {
+                            tmp = list[i];
+                            list[i] = list[j];
+                            list[j] = tmp;
+                        }
+                    }
+                    else
+                    {
+                        for(k = 0; k < 4; k++)
+                        {
+                            if(list[4 * i + k] < list[4 * j + k])
+                            {
+                                tmp = list[4 * i + k];
+                                list[4 * i + k] = list[4 * j + k];
+                                list[4 * j + k] = tmp;
+                            }
+                        }
+                    }
+                }
+            }
+
+            /* Store the median value */
+            if(gray)
+            {
+                dstdata[y * w + x] = list[0 * size / 2];
+            }
+            else
+            {
+                dstdata[4 * (y * w + x)] = list[size / 8 * 4];
+                dstdata[4 * (y * w + x) + 1] = list[size / 8 * 4 + 1];
+                dstdata[4 * (y * w + x) + 2] = list[size / 8 * 4 + 2];
+                dstdata[4 * (y * w + x) + 3] = list[size / 8 * 4 + 3];
+            }
+        }
+    }
+
+    return dst;
+}
+

pipi/pipi.h

@@ -126,6 +126,8 @@ extern pipi_image_t *pipi_brightness(pipi_image_t *, double);
 extern pipi_image_t *pipi_contrast(pipi_image_t *, double);
 extern pipi_image_t *pipi_autocontrast(pipi_image_t *);
 extern pipi_image_t *pipi_invert(pipi_image_t *);
+extern pipi_image_t *pipi_median(pipi_image_t *, int);
+extern pipi_image_t *pipi_median_ext(pipi_image_t *, int, int);
 
 extern pipi_image_t *pipi_tile(pipi_image_t *, int, int);
 extern int pipi_flood_fill(pipi_image_t *,

src/pipi.c

@@ -73,6 +73,14 @@ int main(int argc, char *argv[])
                 return EXIT_FAILURE;
             argv++;
         }
+        else if(!strcmp(argv[0], "--median"))
+        {
+            if(argv[1] == NULL)
+                return EXIT_FAILURE;
+            if(pipi_command(ctx, "median", argv[1]) != 0)
+                return EXIT_FAILURE;
+            argv++;
+        }
         else if(!strcmp(argv[0], "--gray"))
         {
             if(pipi_command(ctx, "gray") != 0)