Noise 2

3 files changed, 356 insertions, 8 deletions

diff --git a/Makefile b/Makefile
index da5af2f..223a5b8 100644
--- a/Makefile
+++ b/Makefile
@@ -2,7 +2,7 @@ CFLAGS=-Wall
 LDFLAGS=-Wall -lm 
 GTKFLAGS=`pkg-config --cflags --libs gtk+-2.0`
 
-all: filters fft-test bluenoise trans
+all: filters fft-test bluenoise trans noise2
 
 filters: fft.c fft.h filters.c
 	$(CC) -o filters fft.c filters.c $(LDFLAGS)
@@ -15,3 +15,6 @@ trans: pngtrans.c fft.c fft.h
 
 bluenoise: fft.c fft.h bluenoise.c gtk-utils.c
 	$(CC) $(GTKFLAGS) -o bluenoise fft.c bluenoise.c gtk-utils.c $(LDFLAGS)
+
+noise2: noise2.c fft.c fft.h
+	$(CC) $(GTKFLAGS) -o noise2 noise2.c fft.c $(LDFLAGS)
diff --git a/noise2.c b/noise2.c
new file mode 100644
index 0000000..5974e40
--- /dev/null
+++ b/noise2.c
@@ -0,0 +1,350 @@
+#include <stdlib.h>
+#include <stdint.h>
+#include <cairo.h>
+#include <gtk/gtk.h>
+#include <gdk/gdkpixbuf.h>
+#include "fft.h"
+
+typedef struct complex_image_t complex_image_t;
+struct complex_image_t
+{
+    int width;
+    int height;
+    complex_t *red;
+    complex_t *green;
+    complex_t *blue;
+};
+
+static complex_image_t *
+complex_image_new (int width, int height)
+{
+    complex_image_t *image = g_new0 (complex_image_t, 1);
+
+    image->width = width;
+    image->height = height;
+    image->red = g_new0 (complex_t, width * height);
+    image->green = g_new0 (complex_t, width * height);
+    image->blue = g_new0 (complex_t, width * height);
+    
+    return image;
+}
+
+static complex_image_t *
+complex_image_copy (complex_image_t *image)
+{
+    complex_image_t *copy = complex_image_new (image->width, image->height);
+    int width = image->width;
+    int height = image->height;
+
+    copy->red = g_memdup (image->red, width * height * sizeof (complex_t));
+    copy->green = g_memdup (image->green, width * height * sizeof (complex_t));
+    copy->blue = g_memdup (image->blue, width * height * sizeof (complex_t));
+
+    return copy;
+}
+
+static void
+complex_image_subtract (complex_image_t *image, complex_image_t *other)
+{
+    int i, j;
+    int h = image->height;
+    int w = image->width;
+    
+    for (i = 0; i < h; ++i)
+    {
+	for (j = 0; j < w; ++j)
+	{
+	    int idx = i * w + j;
+
+	    image->red[idx] = complex_sub (image->red[idx], other->red[idx]);
+	    image->green[idx] = complex_sub (image->green[idx], other->green[idx]);
+	    image->blue[idx] = complex_sub (image->blue[idx], other->blue[idx]);
+	}
+    }
+}
+
+static complex_image_t *
+complex_image_from_pixbuf (GdkPixbuf *pixbuf)
+{
+    complex_image_t *result;
+    uint8_t *pdata;
+    int w, h, s;
+    int i, j;
+    gboolean has_alpha;
+    int n_channels;
+
+    w = gdk_pixbuf_get_width (pixbuf);
+    h = gdk_pixbuf_get_height (pixbuf);
+    s = gdk_pixbuf_get_rowstride (pixbuf);
+    pdata = gdk_pixbuf_get_pixels (pixbuf);
+    has_alpha = gdk_pixbuf_get_has_alpha (pixbuf);
+    n_channels = 3 + has_alpha;
+
+    printf ("has alpha: %d\n", has_alpha);
+    
+    result = complex_image_new (w, h);
+
+    for (i = 0; i < h; ++i)
+    {
+	for (j = 0; j < w; ++j)
+	{
+	    uint8_t *p = &pdata[i * s + j * n_channels];
+	    int idx = i * w + j;
+
+	    result->red[idx].re = p[0] / 255.0;
+	    result->green[idx].re = p[1] / 255.0;
+	    result->blue[idx].re = p[2] / 255.0;
+	}
+    }
+
+    return result;
+}
+
+typedef uint8_t (* convert_t) (complex_t c);
+
+static GdkPixbuf *
+pixbuf_from_complex_image (complex_image_t *image, convert_t convert)
+{
+    int w = image->width;
+    int h = image->height;
+    GdkPixbuf *pixbuf = gdk_pixbuf_new (GDK_COLORSPACE_RGB, TRUE, 8, w, h);
+    uint8_t *p_bits = (uint8_t *)gdk_pixbuf_get_pixels (pixbuf);
+    int s = gdk_pixbuf_get_rowstride (pixbuf);
+    int i, j;
+
+    for (i = 0; i < h; ++i)
+    {
+	for (j = 0; j < w; ++j)
+	{
+	    uint8_t *p = &(p_bits[i * s + j * 4]);
+	    int idx = i * w + j;
+
+	    p[0] = convert (image->red[idx]);
+	    p[1] = convert (image->green[idx]);
+	    p[2] = convert (image->blue[idx]);
+	    p[3] = 0xff; // convert (image->alpha[idx]);
+	}
+    }
+
+    return pixbuf;
+}
+
+#define SIZE 1024
+#define SIZE 1024
+
+static gboolean
+on_expose (GtkWidget *widget, GdkEventExpose *expose, gpointer data)
+{
+    GdkPixbuf *pixbuf = data;
+    int i, j;
+    int pwidth = gdk_pixbuf_get_width (pixbuf);
+    int pheight = gdk_pixbuf_get_height (pixbuf);
+    int dwidth, dheight;
+
+    gdk_drawable_get_size (widget->window, &dwidth, &dheight);
+
+    for (i = 0; i < dheight; i += pheight)
+    {
+	for (j = 0; j < dwidth; j += pwidth)
+	{
+	    gdk_draw_pixbuf (widget->window, NULL,
+			     pixbuf, 0, 0, j, i, pwidth, pheight,
+			     GDK_RGB_DITHER_NONE,
+			     0, 0);
+	}
+    }
+    return TRUE;
+}
+
+typedef enum
+{
+    DISPLAY_MAG,
+    DISPLAY_RE
+} display_type_t;
+
+static uint8_t
+convert_mag (complex_t d)
+{
+    double m = complex_mag (d);
+
+    m = log (m) / 10.0;
+
+    if (m > 1.0)
+    {
+	printf ("%f\n", m);
+	m = 1.0;
+    }
+    
+    if (m < 0)
+	m = 0;
+    
+    return (uint8_t) (m * 255.0 + 0.5);
+}
+
+static uint8_t
+convert_re (complex_t c)
+{
+    if (c.re > 1.0)
+	c.re = 1.0;
+    if (c.re < 0)
+	c.re = 0;
+    
+    return c.re * 255.0 + 0.5;
+}
+
+static void
+display (const char *name, complex_image_t *image, display_type_t type)
+{
+    GtkWidget *window, *da;
+    GdkPixbuf *pixbuf;
+    int argc;
+    char **argv;
+    char *arg0 = g_strdup (name);
+    convert_t convert;
+
+    argc = 1;
+    argv = (char **)&arg0;
+
+    gtk_init (&argc, &argv);
+
+    window = gtk_window_new (GTK_WINDOW_TOPLEVEL);
+
+    da = gtk_drawing_area_new ();
+
+    gtk_container_add (GTK_CONTAINER (window), da);
+
+    gtk_window_set_default_size (GTK_WINDOW (window), SIZE, SIZE);
+
+    switch (type)
+    {
+    case DISPLAY_MAG: convert = convert_mag; break;
+    case DISPLAY_RE: convert = convert_re; break;
+    default:
+	g_assert_not_reached();
+    }
+	    
+    pixbuf = pixbuf_from_complex_image (image, convert);
+
+    g_signal_connect (da, "expose_event", G_CALLBACK (on_expose), pixbuf);
+    g_signal_connect (window, "delete_event", G_CALLBACK (gtk_main_quit), NULL);
+
+    gtk_widget_show_all (window);
+
+    gtk_main ();
+}
+
+static void
+image_fft (complex_image_t *image)
+{
+    fft_shift_2d (image->red, image->width);
+    fft_shift_2d (image->green, image->width);
+    fft_shift_2d (image->blue, image->width);
+}
+
+static void
+image_ifft (complex_image_t *image)
+{
+    ifft_shift_2d (image->red, image->width);
+    ifft_shift_2d (image->green, image->width);
+    ifft_shift_2d (image->blue, image->width);
+}
+
+static complex_t
+filter (complex_t c, double dist)
+{
+    double m = complex_mag (c);
+    double arg = complex_arg (c);
+
+    return complex_from_mag_arg (m * (1/(pow (sqrt (dist) / 16, 4) + 1)) , arg);
+}
+
+static void
+low_pass (complex_image_t *image, double d)
+{
+    int w = image->width;
+    int h = image->height;
+    int i, j;
+    double d2 = d * d;
+    int c = h / 2;
+
+    for (i = 0; i < h; ++i)
+    {
+	for (j = 0; j < w; ++j)
+	{
+	    int idx = i * w + j;
+	    double dist = (c - i) * (c - i) + (c - j) * (c -j);
+	    double t;
+
+	    if (dist > 2)
+	    {
+		image->red[idx] = filter (image->red[idx], dist);
+		image->green[idx] = filter (image->green[idx], dist);
+		image->blue[idx] = filter (image->blue[idx], dist);
+	    }
+	    else
+	    {
+		image->red[idx].re = image->red[idx].im = 0;
+		image->green[idx].re = image->green[idx].im = 0;
+		image->blue[idx].re = image->blue[idx].im = 0;
+	    }
+	}
+    }
+}
+
+static void
+make_noise (complex_image_t *image)
+{
+    int i, j;
+    int h, w;
+
+    h = image->height;
+    w = image->width;
+
+    for (i = 0; i < h; ++i)
+    {
+	for (j = 0; j < w; ++j)
+	{
+	    int idx = i * w + j;
+	    double n = 0.5 * (drand48() + drand48());
+
+	    image->red[idx].re = n;
+	    image->green[idx].re = n;
+	    image->blue[idx].re = n;
+	}
+    }
+}
+
+int
+main (int argc, char **argv)
+{
+    char *input, *output;
+    GdkPixbuf *pb;
+    complex_image_t *image;
+    complex_image_t *low_passed;
+
+    g_type_init ();
+
+    image = complex_image_new (32, 32);
+
+    make_noise (image);
+    
+    display ("test", image, DISPLAY_RE);
+    
+    image_fft (image);
+
+    low_passed = complex_image_copy (image);
+
+    low_pass (low_passed, 2 * image->width);
+
+    display ("low passed", low_passed, DISPLAY_MAG);
+    
+    complex_image_subtract (image, low_passed);
+
+    display ("ssubtracted", image, DISPLAY_MAG);
+    
+    image_ifft (image);
+
+    display ("test", image, DISPLAY_RE);
+    
+    return 0;
+}
diff --git a/pngtrans.c b/pngtrans.c
index 4fb6b4c..dea47e4 100644
--- a/pngtrans.c
+++ b/pngtrans.c
@@ -278,16 +278,11 @@ main (int argc, char **argv)
 
     image_fft (image);
 
-#if 0
+    display ("test", image, DISPLAY_MAG);
+
     low_pass (image, 2 * image->width);
-#endif
 
     display ("test", image, DISPLAY_MAG);
-
-#if 0
-#endif
-#if 0
-#endif
     
     image_ifft (image);