1
/* This file was taken from modified dcraw published by Paul Lee
2
on January 23, 2009, taking dcraw ver.8.90/rev.1.417
4
http://sites.google.com/site/demosaicalgorithms/modified-dcraw
6
As modified dcraw source code was published, the release under
7
GPL Version 2 or later option could be applied, so this file
8
is taken under this premise.
12
Color Demosaicing Using Variance of Color Differences
13
by K.H Chung and Y.H Chan
15
#define PIX_SORT(a,b) { if ((a)>(b)) {temp=(a);(a)=(b);(b)=temp;} }
19
If ahd == 0 then its just VCD
22
void CLASS vcd_interpolate(int ahd_cutoff)
24
int row, col, indx, c, d, i, j;
26
int w1, w2, w3, w4, w5, w6, v0, v1, v2, v3, v4, LH, LV, var_dir;
27
int AHD_cnt, LH_cnt, LV_cnt, varH_cnt, varV_cnt, varD_cnt, T_cnt;
28
int dC0, dC1, dC2, dC3, dC4, temp;
29
double e, T, d1, d3, d5, d7, d9, varH, varV, varD, var0;
33
if (verbose) fprintf(stderr,_("VCD interpolation...\n"));
36
/* assume VCD's T value is based on gamma=2.22 test images */
38
AHD_cnt = LH_cnt = LV_cnt = varH_cnt = varV_cnt = varD_cnt = 0;
45
border_interpolate(6);
47
/* run AHD interpolation on Green channel with threshold */
48
if (ahd_cutoff > 0) ahd_partial_interpolate(ahd_cutoff);
50
/* Interpolate green pixels on red/blue pixel locations */
51
for (row=6; row < height-6; row++)
52
for (col=6+(FC(row,6) & 1), c=FC(row,col); col < width-6; col+=2) {
56
if (image[indx][1] > 0) {
62
ABS(pix[-2-w2][c]-pix[-w2][c]) + ABS(pix[-w2][c]-pix[2-w2][c]) +
63
ABS(pix[-2 ][c]-pix[ 0][c]) + ABS(pix[ 0][c]-pix[2 ][c]) +
64
ABS(pix[-2+w2][c]-pix[ w2][c]) + ABS(pix[ w2][c]-pix[2+w2][c]) +
65
ABS(pix[-2-w1][1]-pix[-w1][1]) + ABS(pix[-w1][1]-pix[2-w1][1]) +
66
ABS(pix[-2+w1][1]-pix[ w1][1]) + ABS(pix[ w1][1]-pix[2+w1][1]) +
67
ABS(pix[-1-w2][1]-pix[-w2][c]) + ABS(pix[-w2][c]-pix[1-w2][1]) +
68
ABS(pix[-1 ][1]-pix[ 0][c]) + ABS(pix[ 0][c]-pix[1 ][1]) +
69
ABS(pix[-1+w2][1]-pix[ w2][c]) + ABS(pix[ w2][c]-pix[1+w2][1]) +
70
ABS(pix[-1-w1][d]-pix[-w1][1]) + ABS(pix[-w1][1]-pix[1-w1][d]) +
71
ABS(pix[-1+w1][d]-pix[ w1][1]) + ABS(pix[ w1][1]-pix[1+w1][d]);
74
ABS(pix[-2-w2][c]-pix[-2][c]) + ABS(pix[-2][c]-pix[-2+w2][c]) +
75
ABS(pix[ -w2][c]-pix[ 0][c]) + ABS(pix[ 0][c]-pix[ w2][c]) +
76
ABS(pix[ 2-w2][c]-pix[ 2][c]) + ABS(pix[ 2][c]-pix[ 2+w2][c]) +
77
ABS(pix[-1-w2][1]-pix[-1][1]) + ABS(pix[-1][1]-pix[-1+w2][1]) +
78
ABS(pix[ 1-w2][1]-pix[ 1][1]) + ABS(pix[ 1][1]-pix[ 1+w2][1]) +
79
ABS(pix[-2-w1][1]-pix[-2][c]) + ABS(pix[-2][c]-pix[-2+w1][1]) +
80
ABS(pix[ -w1][1]-pix[ 0][c]) + ABS(pix[ 0][c]-pix[ w1][1]) +
81
ABS(pix[ 2-w1][1]-pix[ 2][c]) + ABS(pix[ 2][c]-pix[ 2+w1][1]) +
82
ABS(pix[-1-w1][d]-pix[-1][1]) + ABS(pix[-1][1]-pix[-1+w1][d]) +
83
ABS(pix[ 1-w1][d]-pix[ 1][1]) + ABS(pix[ 1][1]-pix[ 1+w1][d]);
85
e = (double)(LH) / (double)(LV);
96
/* varH, varV, varD: Eq. (11)~(18) */
99
d1 = (double)(pix[-6][c]-2*(pix[-5][1]-pix[-4][c]+pix[-3][1])+pix[-2][c])/65535.0;
100
d3 = (double)(pix[-4][c]-2*(pix[-3][1]-pix[-2][c]+pix[-1][1])+pix[ 0][c])/65535.0;
101
d5 = (double)(pix[-2][c]-2*(pix[-1][1]-pix[ 0][c]+pix[ 1][1])+pix[ 2][c])/65535.0;
102
d7 = (double)(pix[ 0][c]-2*(pix[ 1][1]-pix[ 2][c]+pix[ 3][1])+pix[ 4][c])/65535.0;
103
d9 = (double)(pix[ 2][c]-2*(pix[ 3][1]-pix[ 4][c]+pix[ 5][1])+pix[ 6][c])/65535.0;
104
/* variance assuming d2 =(d1+d3)/2, d4=(d3+d5)/2, and etc */
106
d1*(18*d1 - 3*d3 - 12*d5 - 12*d7 - 9*d9) +
107
d3*(19*d3 - 7*d5 - 16*d7 - 12*d9) +
108
d5*(19*d5 - 7*d7 - 12*d9) +
114
d1 = (double)(pix[-w6][c]-2*(pix[-w5][1]-pix[-w4][c]+pix[-w3][1])+pix[-w2][c])/65535.0;
115
d3 = (double)(pix[-w4][c]-2*(pix[-w3][1]-pix[-w2][c]+pix[-w1][1])+pix[ 0][c])/65535.0;
116
d5 = (double)(pix[-w2][c]-2*(pix[-w1][1]-pix[ 0][c]+pix[ w1][1])+pix[ w2][c])/65535.0;
117
d7 = (double)(pix[ 0][c]-2*(pix[ w1][1]-pix[ w2][c]+pix[ w3][1])+pix[ w4][c])/65535.0;
118
d9 = (double)(pix[ w2][c]-2*(pix[ w3][1]-pix[ w4][c]+pix[ w5][1])+pix[ w6][c])/65535.0;
120
d1*(18*d1 - 3*d3 - 12*d5 - 12*d7 - 9*d9) +
121
d3*(19*d3 - 7*d5 - 16*d7 - 12*d9) +
122
d5*(19*d5 - 7*d7 - 12*d9) +
129
d1 = (double)(4*pix[-4 ][c]+pix[-4-w2][c]+pix[-6][c]+pix[-2 ][c]+pix[-4+w2][c]
130
-2*(pix[-4-w1][1]+pix[-5 ][1]+pix[-3][1]+pix[-4+w1][1]))/65535.0;
131
d3 = (double)(4*pix[-2 ][c]+pix[-2-w2][c]+pix[-4][c]+pix[ 0 ][c]+pix[-2+w2][c]
132
-2*(pix[-2-w1][1]+pix[-3 ][1]+pix[-1][1]+pix[-2+w1][1]))/65535.0;
133
d5 = (double)(4*pix[ 0 ][c]+pix[ -w2][c]+pix[-2][c]+pix[ 2 ][c]+pix[ w2][c]
134
-2*(pix[ -w1][1]+pix[-1 ][1]+pix[ 1][1]+pix[ w1][1]))/65535.0;
135
d7 = (double)(4*pix[ 2 ][c]+pix[ 2-w2][c]+pix[ 0][c]+pix[ 4 ][c]+pix[ 2+w2][c]
136
-2*(pix[ 2-w1][1]+pix[ 1 ][1]+pix[ 3][1]+pix[ 2+w1][1]))/65535.0;
137
d9 = (double)(4*pix[ 4 ][c]+pix[ 4-w2][c]+pix[ 2][c]+pix[ 6 ][c]+pix[ 4+w2][c]
138
-2*(pix[ 4-w1][1]+pix[ 3 ][1]+pix[ 5][1]+pix[ 4+w1][1]))/65535.0;
140
d1*(18*d1 - 3*d3 - 12*d5 - 12*d7 - 9*d9) +
141
d3*(19*d3 - 7*d5 - 16*d7 - 12*d9) +
142
d5*(19*d5 - 7*d7 - 12*d9) +
146
d1 = (double)(4*pix[-w4][c]+pix[-w6 ][c]+pix[-w4-2][c]+pix[-w4+2][c]+pix[-w2][c]
147
-2*(pix[-w5][1]+pix[-w4-1][1]+pix[-w4+1][1]+pix[-w3 ][1]))/65535.0;
148
d3 = (double)(4*pix[-w2][c]+pix[-w4 ][c]+pix[-w2-2][c]+pix[-w2+2][c]+pix[ 0][c]
149
-2*(pix[-w3][1]+pix[-w2-1][1]+pix[-w2+1][1]+pix[-w1 ][1]))/65535.0;
150
d7 = (double)(4*pix[ w2][c]+pix[ 0][c]+pix[ w2-2][c]+pix[ w2+2][c]+pix[ w4][c]
151
-2*(pix[ w1][1]+pix[ w2-1][1]+pix[ w2+1][1]+pix[ w3 ][1]))/65535.0;
152
d9 = (double)(4*pix[ w4][c]+pix[ w2 ][c]+pix[ w4-2][c]+pix[ w4+2][c]+pix[ w6][c]
153
-2*(pix[ w3][1]+pix[ w4-1][1]+pix[ w4+1][1]+pix[ w5 ][1]))/65535.0;
155
d1*(18*d1 - 3*d3 - 12*d5 - 12*d7 - 9*d9) +
156
d3*(19*d3 - 7*d5 - 16*d7 - 12*d9) +
157
d5*(19*d5 - 7*d7 - 12*d9) +
160
/* scale varD to equalize to varH and varV */
162
if (varD < var0) var_dir = 3;
166
else if (var_dir == 2)
168
else if (var_dir == 3)
172
/* limit values within surrounding green values
173
for overshooted pixel values, revert to linear interpolation */
175
/* Eq.(3) - Horizontal */
177
(2*(pix[ -1][1]+pix[0][c]+pix[ 1][1])-pix[ -2][c]-pix[ 2][c]+2) >> 2;
182
v3 = MAX(2*v1 - v2,0);
183
v4 = MIN(2*v2 - v1,65535);
184
if (v0 < v3 || v0 > v4) {
185
v0 = (pix[-3][1] + pix[3][1] +
186
18*(2*pix[0][c] - pix[-2][c] - pix[2][c]) +
187
63*(pix[-1][1] + pix[1][1]) + 64) >> 7;
188
if (v0 < v3 || v0 > v4) {
189
v0 = (4*(v1 + v2) + 2*pix[0][c]-pix[-2][c]-pix[2][c] + 4) >> 3;
190
/* Bi-linear if anti-aliasing overshoots */
191
if (v0 < v3 || v0 > v4) v0 = (v1 + v2 + 1) >> 1; } }
193
else if (var_dir == 2) {
194
/* Eq.(4) - Vertical */
196
(2*(pix[-w1][1]+pix[0][c]+pix[w1][1])-pix[-w2][c]-pix[w2][c]+2) >> 2;
201
v3 = MAX(2*v1 - v2,0);
202
v4 = MIN(2*v2 - v1,65535);
203
if (v0 < v3 || v0 > v4) {
204
v0 = (pix[-w3][1] + pix[w3][1] +
205
18*(2*pix[0][c] - pix[-w2][c] - pix[w2][c]) +
206
63*(pix[-w1][1] + pix[w1][1]) + 64) >> 7;
207
if (v0 < v3 || v0 > v4) {
208
v0 = (4*(v1 + v2) + 2*pix[0][c]-pix[-w2][c]-pix[w2][c] + 4) >> 3;
209
/* Bi-linear if anti-aliasing overshoots */
210
if (v0 < v3 || v0 > v4) v0 = (v1 + v2 + 1) >> 1; }}
212
else if (var_dir == 3) {
213
/* Eq.(5) - Diagonal */
214
v0 = 2*(pix[ -1][1]+pix[0][c]+pix[ 1][1])-pix[ -2][c]-pix[ 2][c]+2;
215
v0 += 2*(pix[-w1][1]+pix[0][c]+pix[w1][1])-pix[-w2][c]-pix[w2][c]+2;
217
v1 = MIN(MIN(pix[-1][1],pix[1][1]),MIN(pix[-w1][1],pix[w1][1]));
218
v2 = MAX(MAX(pix[-1][1],pix[1][1]),MAX(pix[-w1][1],pix[w1][1]));
219
v3 = MAX(2*v1 - v2,0);
220
v4 = MIN(2*v2 - v1,65535);
221
if (v0 < v3 || v0 > v4)
222
v0 = (pix[-w1][1] + pix[ -1][1] + pix[ 1][1] + pix[w1][1] + 2) >> 2;
224
else if (var_dir == 4) {
231
Interpolote red/blue pixels on BLUE/RED pixel locations
232
using pattern regcognition on differential color plane
234
for (row=1; row < height-1; row++)
235
for (col=1+(FC(row,1) & 1), c=2-FC(row,col); col < width-1; col+=2) {
236
indx = row*width+col;
238
if (pix[0][c] == 0) {
239
dC1 = pix[-w1-1][1]-pix[-w1-1][c];
240
dC2 = pix[-w1+1][1]-pix[-w1+1][c];
241
dC3 = pix[ w1-1][1]-pix[ w1-1][c];
242
dC4 = pix[ w1+1][1]-pix[ w1+1][c];
243
dC0 = dC1 + dC2 + dC3 + dC4;
248
j = (dC1 > dC0) + (dC2 > dC0) + (dC3 > dC0) + (dC4 > dC0);
249
if (j == 3 || j == 1) {
250
/* edge-corner pattern: median of color differential values */
257
/* stripe pattern: average along diagonal */
258
v1 = ABS(pix[-w1-1][c]-pix[w1+1][c]);
259
v2 = ABS(pix[-w1+1][c]-pix[w1-1][c]);
266
v0 = (((int)(pix[0][1]) << 3) - dC0 + 4) >> 3;
267
/* apply anti-aliasing if overshoot */
268
if (v0 < 0 || v0 > 65535)
269
v0 = (pix[-w1-1][c]+pix[-w1+1][c]+pix[w1-1][c]+pix[w1+1][c]+2) >> 2;
274
Interpolote red/blue pixels on GREEN pixel locations
275
using pattern regcognition on differential color plane
277
for (row=1; row < height-1; row++)
278
for (col=1+(FC(row,2) & 1), c=FC(row,col+1); col < width-1; col+=2) {
279
indx = row*width+col;
281
for (i=0; i < 2; c=2-c, i++) {
282
if (pix[0][c] == 0) {
283
dC1 = pix[-w1][1]-pix[-w1][c];
284
dC2 = pix[ -1][1]-pix[ -1][c];
285
dC3 = pix[ 1][1]-pix[ 1][c];
286
dC4 = pix[ w1][1]-pix[ w1][c];
287
dC0 = dC1 + dC2 + dC3 + dC4;
292
j = (dC1 > dC0) + (dC2 > dC0) + (dC3 > dC0) + (dC4 > dC0);
293
if (j == 3 || j == 1) {
294
/* edge-corner pattern: median of color differential values */
301
/* stripe pattern: average along diagonal */
302
v1 = ABS(pix[-w1][c]-pix[w1][c]);
303
v2 = ABS(pix[ -1][c]-pix[ 1][c]);
310
v0 = (((int)(pix[0][1]) << 3) - dC0 + 4) >> 3;
311
/* apply anti-aliasing if overshoot */
312
if (v0 < 0 || v0 > 65535) {
314
v0 = (pix[ -1][c]+pix[ 1][c]+1) >> 1;
316
v0 = (pix[-w1][c]+pix[w1][c]+1) >> 1; }
322
/* Compute statistics */
325
if (ahd_cutoff > 0) {
326
T_cnt = AHD_cnt + LH_cnt + LV_cnt + varH_cnt + varV_cnt + varD_cnt;
328
_("\tAHD, LH, LV, varH, varV varD = %4.2f, %4.2f, %4.2f, %4.2f, %4.2f, %4.2f (%%)\n"),
329
100*(double)AHD_cnt/(double)T_cnt,
330
100*(double)LH_cnt/(double)T_cnt,
331
100*(double)LV_cnt/(double)T_cnt,
332
100*(double)varH_cnt/(double)T_cnt,
333
100*(double)varV_cnt/(double)T_cnt,
334
100*(double)varD_cnt/(double)T_cnt); }
336
T_cnt = LH_cnt + LV_cnt + varH_cnt + varV_cnt + varD_cnt;
338
_("\tLH, LV, varH, varV varD = %4.2f, %4.2f, %4.2f, %4.2f, %4.2f (%%)\n"),
339
100*(double)LH_cnt/(double)T_cnt,
340
100*(double)LV_cnt/(double)T_cnt,
341
100*(double)varH_cnt/(double)T_cnt,
342
100*(double)varV_cnt/(double)T_cnt,
343
100*(double)varD_cnt/(double)T_cnt); }
348
dt = ((double)(t2-t1)) / CLOCKS_PER_SEC;
349
if (verbose) fprintf(stderr,_("\telapsed time = %5.3fs\n"),dt);
added
removed
libraw/demosaic-pack-GPL2/vcd_interpolate.c
