summaryrefslogtreecommitdiff
path: root/noise2.c
blob: 606bd2382cd3c3d24490ad1772ae944a6918b8d3 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
#include <stdlib.h>
#include <stdint.h>
#include <cairo.h>
#include <gtk/gtk.h>
#include <gdk/gdkpixbuf.h>
#include "fft.h"
#include "image.h"

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);
    
    complex_image_show ("test", image, CONVERT_RE);
    
    image_fft (image);

    low_passed = complex_image_copy (image);

    low_pass (low_passed, 2 * image->width);

    complex_image_show ("low passed", low_passed, CONVERT_MAG);
    
    complex_image_subtract (image, low_passed);

    complex_image_show ("ssubtracted", image, CONVERT_MAG);
    
    image_ifft (image);

    complex_image_show ("test", image, CONVERT_RE);
    
    return 0;
}