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
|
#define IMAGER_NO_CONTEXT
#include "imager.h"
#include <math.h>
static double
gauss(int x, double std) {
return 1.0/(sqrt(2.0*PI)*std)*exp(-(double)(x)*(double)(x)/(2*std*std));
}
/* Counters are as follows
l: lines
i: columns
c: filter coeffs
ch: channels
pc: coeff equalization
*/
int
i_gaussian(i_img *im, double stddev) {
int i, c, ch;
i_img_dim x, y;
double pc;
double *coeff;
double res[MAXCHANNELS];
i_img *timg;
int radius, diameter;
dIMCTXim(im);
im_log((aIMCTX, 1,"i_gaussian(im %p, stdev %.2f)\n",im,stddev));
i_clear_error();
if (stddev <= 0) {
i_push_error(0, "stddev must be positive");
return 0;
}
/* totally silly cutoff */
if (stddev > 1000) {
stddev = 1000;
}
timg = i_sametype(im, im->xsize, im->ysize);
if (im->bits <= 8)
radius = ceil(2 * stddev);
else
radius = ceil(3 * stddev);
diameter = 1 + radius * 2;
coeff = mymalloc(sizeof(double) * diameter);
for(i=0;i <= radius;i++)
coeff[radius + i]=coeff[radius - i]=gauss(i, stddev);
pc=0;
for(i=0; i < diameter; i++)
pc+=coeff[i];
for(i=0;i < diameter;i++)
coeff[i] /= pc;
#code im->bits <= 8
IM_COLOR rcolor;
for(y = 0; y < im->ysize; y++) {
for(x = 0; x < im->xsize; x++) {
pc=0.0;
for(ch=0;ch<im->channels;ch++)
res[ch]=0;
for(c = 0;c < diameter; c++)
if (IM_GPIX(im,x+c-radius,y,&rcolor)!=-1) {
for(ch=0;ch<im->channels;ch++)
res[ch]+= rcolor.channel[ch] * coeff[c];
pc+=coeff[c];
}
for(ch=0;ch<im->channels;ch++) {
double value = res[ch] / pc;
rcolor.channel[ch] = value > IM_SAMPLE_MAX ? IM_SAMPLE_MAX : IM_ROUND(value);
}
IM_PPIX(timg, x, y, &rcolor);
}
}
for(x = 0;x < im->xsize; x++) {
for(y = 0; y < im->ysize; y++) {
pc=0.0;
for(ch=0; ch<im->channels; ch++)
res[ch]=0;
for(c=0; c < diameter; c++)
if (IM_GPIX(timg, x, y+c-radius, &rcolor)!=-1) {
for(ch=0;ch<im->channels;ch++)
res[ch]+= rcolor.channel[ch] * coeff[c];
pc+=coeff[c];
}
for(ch=0;ch<im->channels;ch++) {
double value = res[ch]/pc;
rcolor.channel[ch] = value > IM_SAMPLE_MAX ? IM_SAMPLE_MAX : IM_ROUND(value);
}
IM_PPIX(im, x, y, &rcolor);
}
}
#/code
myfree(coeff);
i_img_destroy(timg);
return 1;
}
|
