?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Bugfix: the defocus node "no zbuffer" settings was automatically set
[plumiferos.git] / source / blender / nodes / intern / CMP_nodes / CMP_defocus.c
blob069c0be3c7ce9866d1bdac482108600684d46ef1
1 /**
2 * $Id: CMP_defocus.c 10741 2007-05-20 15:52:29Z eeshlo $
3 *
4 * ***** BEGIN GPL LICENSE BLOCK *****
5 *
6 * This program is free software; you can redistribute it and/or
7 * modify it under the terms of the GNU General Public License
8 * as published by the Free Software Foundation; either version 2
9 * of the License, or (at your option) any later version.
10 *
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software Foundation,
18 * Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
19 *
20 * The Original Code is Copyright (C) 2006 Blender Foundation.
21 * All rights reserved.
22 *
23 * The Original Code is: all of this file.
24 *
25 * Contributor(s): none yet.
26 *
27 * ***** END GPL LICENSE BLOCK *****
28 */
29
30 #include "../CMP_util.h"
31
32 /* ************ qdn: Defocus node ****************** */
33 static bNodeSocketType cmp_node_defocus_in[]= {
34 { SOCK_RGBA, 1, "Image", 0.8f, 0.8f, 0.8f, 1.0f, 0.0f, 1.0f},
35 { SOCK_VALUE, 1, "Z", 0.8f, 0.8f, 0.8f, 1.0f, 0.0f, 1.0f},
36 { -1, 0, "" }
37 };
38 static bNodeSocketType cmp_node_defocus_out[]= {
39 { SOCK_RGBA, 0, "Image", 0.0f, 0.0f, 0.0f, 1.0f, 0.0f, 1.0f},
40 { -1, 0, "" }
41 };
42
43
44 // line coefs for point sampling & scancon. data.
45 typedef struct BokehCoeffs {
46 float x0, y0, dx, dy;
47 float ls_x, ls_y;
48 float min_x, min_y, max_x, max_y;
49 } BokehCoeffs;
50
51 // returns array of BokehCoeffs
52 // returns length of array in 'len_bkh',
53 // radius squared of inscribed disk in 'inradsq', needed in getWeight() test,
54 // BKH[8] is the data returned for the bokeh shape & bkh_b[4] is it's 2d bound
55 static void makeBokeh(char bktype, char ro, int* len_bkh, float* inradsq, BokehCoeffs BKH[8], float bkh_b[4])
56 {
57 float x0, x1, y0, y1, dx, dy, iDxy;
58 float w = MAX2(1e-5f, ro)*M_PI/180.f; // never reported stangely enough, but a zero offset causes missing center line...
59 float wi = (360.f/bktype)*M_PI/180.f;
60 int i, ov, nv;
61
62 // bktype must be at least 3 & <= 8
63 bktype = (bktype<3) ? 3 : ((bktype>8) ? 8 : bktype);
64 *len_bkh = bktype;
65 *inradsq = -1.f;
66
67 for (i=0; i<(*len_bkh); i++) {
68 x0 = cos(w);
69 y0 = sin(w);
70 w += wi;
71 x1 = cos(w);
72 y1 = sin(w);
73 if ((*inradsq)<0.f) {
74 // radius squared of inscribed disk
75 float idx=(x0+x1)*0.5f, idy=(y0+y1)*0.5f;
76 *inradsq = idx*idx + idy*idy;
77 }
78 BKH[i].x0 = x0;
79 BKH[i].y0 = y0;
80 dx = x1-x0, dy = y1-y0;
81 iDxy = 1.f / sqrt(dx*dx + dy*dy);
82 dx *= iDxy;
83 dy *= iDxy;
84 BKH[i].dx = dx;
85 BKH[i].dy = dy;
86 }
87
88 // precalc scanconversion data
89 // bokeh bound, not transformed, for scanconvert
90 bkh_b[0] = bkh_b[2] = 1e10f; // xmin/ymin
91 bkh_b[1] = bkh_b[3] = -1e10f; // xmax/ymax
92 ov = (*len_bkh) - 1;
93 for (nv=0; nv<(*len_bkh); nv++) {
94 bkh_b[0] = MIN2(bkh_b[0], BKH[nv].x0); // xmin
95 bkh_b[1] = MAX2(bkh_b[1], BKH[nv].x0); // xmax
96 bkh_b[2] = MIN2(bkh_b[2], BKH[nv].y0); // ymin
97 bkh_b[3] = MAX2(bkh_b[3], BKH[nv].y0); // ymax
98 BKH[nv].min_x = MIN2(BKH[ov].x0, BKH[nv].x0);
99 BKH[nv].max_x = MAX2(BKH[ov].x0, BKH[nv].x0);
100 BKH[nv].min_y = MIN2(BKH[ov].y0, BKH[nv].y0);
101 BKH[nv].max_y = MAX2(BKH[ov].y0, BKH[nv].y0);
102 dy = BKH[nv].y0 - BKH[ov].y0;
103 BKH[nv].ls_x = (BKH[nv].x0 - BKH[ov].x0) / ((dy==0.f) ? 1.f : dy);
104 BKH[nv].ls_y = (BKH[nv].ls_x==0.f) ? 1.f : (1.f/BKH[nv].ls_x);
105 ov = nv;
106 }
107 }
108
109 // test if u/v inside shape & returns weight value
110 static float getWeight(BokehCoeffs* BKH, int len_bkh, float u, float v, float rad, float inradsq)
111 {
112 BokehCoeffs* bc = BKH;
113 float cdist, irad = (rad==0.f) ? 1.f : (1.f/rad);
114 u *= irad;
115 v *= irad;
116
117 // early out test1: if point outside outer unit disk, it cannot be inside shape
118 cdist = u*u + v*v;
119 if (cdist>1.f) return 0.f;
120
121 // early out test2: if point inside or on inner disk, point must be inside shape
122 if (cdist<=inradsq) return 1.f;
123
124 while (len_bkh--) {
125 if ((bc->dy*(u - bc->x0) - bc->dx*(v - bc->y0)) > 0.f) return 0.f;
126 bc++;
127 }
128 return 1.f;
129 }
130
131 // QMC.seq. for sampling, A.Keller, EMS
132 static float RI_vdC(unsigned int bits, unsigned int r)
133 {
134 bits = ( bits << 16) | ( bits >> 16);
135 bits = ((bits & 0x00ff00ff) << 8) | ((bits & 0xff00ff00) >> 8);
136 bits = ((bits & 0x0f0f0f0f) << 4) | ((bits & 0xf0f0f0f0) >> 4);
137 bits = ((bits & 0x33333333) << 2) | ((bits & 0xcccccccc) >> 2);
138 bits = ((bits & 0x55555555) << 1) | ((bits & 0xaaaaaaaa) >> 1);
139 bits ^= r;
140 return (float)((double)bits / 4294967296.0);
141 }
142
143 // single channel IIR gaussian filtering
144 // much faster than anything else, constant time independent of width
145 // should extend to multichannel and make this a node, could be useful
146 static void IIR_gauss_single(CompBuf* buf, float sigma)
147 {
148 double q, q2, sc, cf[4], tsM[9], tsu[3], tsv[3];
149 float *X, *Y, *W;
150 int i, x, y, sz;
151
152 // single channel only for now
153 if (buf->type != CB_VAL) return;
154
155 // <0.5 not valid, though can have a possibly useful sort of sharpening effect
156 if (sigma < 0.5) return;
157
158 // see "Recursive Gabor Filtering" by Young/VanVliet
159 // all factors here in double.prec. Required, because for single.prec it seems to blow up if sigma > ~200
160 if (sigma >= 3.556)
161 q = 0.9804*(sigma - 3.556) + 2.5091;
162 else // sigma >= 0.5
163 q = (0.0561*sigma + 0.5784)*sigma - 0.2568;
164 q2 = q*q;
165 sc = (1.1668 + q)*(3.203729649 + (2.21566 + q)*q);
166 // no gabor filtering here, so no complex multiplies, just the regular coefs.
167 // all negated here, so as not to have to recalc Triggs/Sdika matrix
168 cf[1] = q*(5.788961737 + (6.76492 + 3.0*q)*q)/ sc;
169 cf[2] = -q2*(3.38246 + 3.0*q)/sc;
170 // 0 & 3 unchanged
171 cf[3] = q2*q/sc;
172 cf[0] = 1.0 - cf[1] - cf[2] - cf[3];
173
174 // Triggs/Sdika border corrections,
175 // it seems to work, not entirely sure if it is actually totally correct,
176 // Besides J.M.Geusebroek's anigauss.c (see http://www.science.uva.nl/~mark),
177 // found one other implementation by Cristoph Lampert,
178 // but neither seem to be quite the same, result seems to be ok sofar anyway.
179 // Extra scale factor here to not have to do it in filter,
180 // though maybe this had something to with the precision errors
181 sc = cf[0]/((1.0 + cf[1] - cf[2] + cf[3])*(1.0 - cf[1] - cf[2] - cf[3])*(1.0 + cf[2] + (cf[1] - cf[3])*cf[3]));
182 tsM[0] = sc*(-cf[3]*cf[1] + 1.0 - cf[3]*cf[3] - cf[2]);
183 tsM[1] = sc*((cf[3] + cf[1])*(cf[2] + cf[3]*cf[1]));
184 tsM[2] = sc*(cf[3]*(cf[1] + cf[3]*cf[2]));
185 tsM[3] = sc*(cf[1] + cf[3]*cf[2]);
186 tsM[4] = sc*(-(cf[2] - 1.0)*(cf[2] + cf[3]*cf[1]));
187 tsM[5] = sc*(-(cf[3]*cf[1] + cf[3]*cf[3] + cf[2] - 1.0)*cf[3]);
188 tsM[6] = sc*(cf[3]*cf[1] + cf[2] + cf[1]*cf[1] - cf[2]*cf[2]);
189 tsM[7] = sc*(cf[1]*cf[2] + cf[3]*cf[2]*cf[2] - cf[1]*cf[3]*cf[3] - cf[3]*cf[3]*cf[3] - cf[3]*cf[2] + cf[3]);
190 tsM[8] = sc*(cf[3]*(cf[1] + cf[3]*cf[2]));
191
192 #define YVV(L)\
193 {\
194 W[0] = cf[0]*X[0] + cf[1]*X[0] + cf[2]*X[0] + cf[3]*X[0];\
195 W[1] = cf[0]*X[1] + cf[1]*W[0] + cf[2]*X[0] + cf[3]*X[0];\
196 W[2] = cf[0]*X[2] + cf[1]*W[1] + cf[2]*W[0] + cf[3]*X[0];\
197 for (i=3; i<L; i++)\
198 W[i] = cf[0]*X[i] + cf[1]*W[i-1] + cf[2]*W[i-2] + cf[3]*W[i-3];\
199 tsu[0] = W[L-1] - X[L-1];\
200 tsu[1] = W[L-2] - X[L-1];\
201 tsu[2] = W[L-3] - X[L-1];\
202 tsv[0] = tsM[0]*tsu[0] + tsM[1]*tsu[1] + tsM[2]*tsu[2] + X[L-1];\
203 tsv[1] = tsM[3]*tsu[0] + tsM[4]*tsu[1] + tsM[5]*tsu[2] + X[L-1];\
204 tsv[2] = tsM[6]*tsu[0] + tsM[7]*tsu[1] + tsM[8]*tsu[2] + X[L-1];\
205 Y[L-1] = cf[0]*W[L-1] + cf[1]*tsv[0] + cf[2]*tsv[1] + cf[3]*tsv[2];\
206 Y[L-2] = cf[0]*W[L-2] + cf[1]*Y[L-1] + cf[2]*tsv[0] + cf[3]*tsv[1];\
207 Y[L-3] = cf[0]*W[L-3] + cf[1]*Y[L-2] + cf[2]*Y[L-1] + cf[3]*tsv[0];\
208 for (i=L-4; i>=0; i--)\
209 Y[i] = cf[0]*W[i] + cf[1]*Y[i+1] + cf[2]*Y[i+2] + cf[3]*Y[i+3];\
210 }
211
212 // intermediate buffers
213 sz = MAX2(buf->x, buf->y);
214 Y = MEM_callocN(sz*sizeof(float), "IIR_gauss Y buf");
215 W = MEM_callocN(sz*sizeof(float), "IIR_gauss W buf");
216 // H
217 for (y=0; y<buf->y; y++) {
218 X = &buf->rect[y*buf->x];
219 YVV(buf->x);
220 memcpy(X, Y, sizeof(float)*buf->x);
221 }
222 // V
223 X = MEM_callocN(buf->y*sizeof(float), "IIR_gauss X buf");
224 for (x=0; x<buf->x; x++) {
225 for (y=0; y<buf->y; y++)
226 X[y] = buf->rect[x + y*buf->x];
227 YVV(buf->y);
228 for (y=0; y<buf->y; y++)
229 buf->rect[x + y*buf->x] = Y[y];
230 }
231 MEM_freeN(X);
232
233 MEM_freeN(W);
234 MEM_freeN(Y);
235 #undef YVV
236 }
237
238 static void defocus_blur(bNode *node, CompBuf *new, CompBuf *img, CompBuf *zbuf, float inpval, int no_zbuf)
239 {
240 NodeDefocus *nqd = node->storage;
241 CompBuf *wts; // weights buffer
242 CompBuf *crad; // CoC radius buffer
243 BokehCoeffs BKH[8]; // bokeh shape data, here never > 8 pts.
244 float bkh_b[4] = {0}; // shape 2D bound
245 unsigned int p, px, p4, zp, cp, cp4;
246 float *ctcol, u, v, iZ, ct_crad, lwt, wt=0, cR2=0;
247 float dof_sp, maxfgc, bk_hn_theta=0, inradsq=0;
248 float cam_fdist=1, cam_invfdist=1, cam_lens=35;
249 int x, y, sx, sy, len_bkh=0;
250 float aspect, aperture;
251 int minsz;
252 //float bcrad, nmaxc, scf;
253
254 // get some required params from the current scene camera
255 // (ton) this is wrong, needs fixed
256 Object* camob = G.scene->camera;
257 if (camob && camob->type==OB_CAMERA) {
258 Camera* cam = (Camera*)camob->data;
259 cam_lens = cam->lens;
260 cam_fdist = (cam->YF_dofdist==0.f) ? 1e10f : cam->YF_dofdist;
261 cam_invfdist = 1.f/cam_fdist;
262 }
263
264 // guess work here.. best match with raytraced result
265 minsz = MIN2(img->x, img->y);
266 dof_sp = (float)minsz / (16.f / cam_lens); // <- == aspect * MIN2(img->x, img->y) / tan(0.5f * fov);
267
268 // aperture
269 aspect = (img->x > img->y) ? (img->y / (float)img->x) : (img->x / (float)img->y);
270 aperture = 0.5f*(cam_lens / (aspect*32.f)) / nqd->fstop;
271
272 // if not disk, make bokeh coefficients and other needed data
273 if (nqd->bktype!=0) {
274 makeBokeh(nqd->bktype, nqd->rotation, &len_bkh, &inradsq, BKH, bkh_b);
275 bk_hn_theta = 0.5 * nqd->bktype * sin(2.0 * M_PI / nqd->bktype); // weight factor
276 }
277
278 // accumulated weights
279 wts = alloc_compbuf(img->x, img->y, CB_VAL, 1);
280 // CoC radius buffer
281 crad = alloc_compbuf(img->x, img->y, CB_VAL, 1);
282
283 // if 'no_zbuf' flag set (which is always set if input is not an image),
284 // values are instead interpreted directly as blur radius values
285 if (no_zbuf) {
286 // to prevent *reaaallly* big radius values and impossible calculation times,
287 // limit the maximum to half the image width or height, whichever is smaller
288 float maxr = 0.5f*(float)MIN2(img->x, img->y);
289 for (p=0; p<(unsigned int)(img->x*img->y); p++) {
290 crad->rect[p] = zbuf ? (zbuf->rect[p]*nqd->scale) : inpval;
291 // bug #5921, limit minimum
292 crad->rect[p] = MAX2(1e-5f, crad->rect[p]);
293 crad->rect[p] = MIN2(crad->rect[p], maxr);
294 // if maxblur!=0, limit maximum
295 if (nqd->maxblur != 0.f) crad->rect[p] = MIN2(crad->rect[p], nqd->maxblur);
296 }
297 }
298 else {
299 // actual zbuffer.
300 // separate foreground from background CoC's
301 // then blur background and blend in again with foreground,
302 // improves the 'blurred foreground overlapping in-focus midground' sharp boundary problem.
303 // wts buffer here used for blendmask
304 maxfgc = 0.f; // maximum foreground CoC radius
305 for (y=0; y<img->y; y++) {
306 p = y * img->x;
307 for (x=0; x<img->x; x++) {
308 px = p + x;
309 iZ = (zbuf->rect[px]==0.f) ? 0.f : (1.f/zbuf->rect[px]);
310 crad->rect[px] = 0.5f*(aperture*(dof_sp*(cam_invfdist - iZ) - 1.f));
311 if (crad->rect[px] <= 0.f) {
312 wts->rect[px] = 1.f;
313 crad->rect[px] = -crad->rect[px];
314 if (crad->rect[px] > maxfgc) maxfgc = crad->rect[px];
315 }
316 else crad->rect[px] = wts->rect[px] = 0;
317 }
318 }
319
320 // fast blur...
321 // bug #6656 part 1, probably when previous node_composite.c was split into separate files, it was not properly updated
322 // to include recent cvs commits (well, at least not defocus node), so this part was missing...
323 wt = aperture*128.f;
324 IIR_gauss_single(crad, wt);
325 IIR_gauss_single(wts, wt);
326
327 // bug #6656 part 2a, although foreground blur is not based anymore on closest object,
328 // the rescaling op below was still based on that anyway, and unlike the comment in below code,
329 // the difference is therefore not always that small at all...
330 // so for now commented out, not sure if this is going to cause other future problems, lets just wait and see...
331 /*
332 // find new maximum to scale it back to original
333 // (could skip this, not strictly necessary, in general, difference is quite small, but just in case...)
334 nmaxc = 0;
335 for (p=0; p<(img->x*img->y); p++)
336 if (crad->rect[p] > nmaxc) nmaxc = crad->rect[p];
337 // rescale factor
338 scf = (nmaxc==0.f) ? 1.f: (maxfgc / nmaxc);
339 */
340
341 // and blend...
342 for (y=0; y<img->y; y++) {
343 p = y*img->x;
344 for (x=0; x<img->x; x++) {
345 px = p + x;
346 if (zbuf->rect[px]!=0.f) {
347 iZ = (zbuf->rect[px]==0.f) ? 0.f : (1.f/zbuf->rect[px]);
348
349 // bug #6656 part 2b, do not rescale
350 /*
351 bcrad = 0.5f*fabs(aperture*(dof_sp*(cam_invfdist - iZ) - 1.f));
352 // scale crad back to original maximum and blend
353 crad->rect[px] = bcrad + wts->rect[px]*(scf*crad->rect[px] - bcrad);
354 */
355 crad->rect[px] = 0.5f*fabs(aperture*(dof_sp*(cam_invfdist - iZ) - 1.f));
356
357 // 'bug' #6615, limit minimum radius to 1 pixel, not really a solution, but somewhat mitigates the problem
358 crad->rect[px] = MAX2(crad->rect[px], 0.5f);
359 // if maxblur!=0, limit maximum
360 if (nqd->maxblur != 0.f) crad->rect[px] = MIN2(crad->rect[px], nqd->maxblur);
361 }
362 else crad->rect[px] = 0.f;
363 // clear weights for next part
364 wts->rect[px] = 0.f;
365 }
366 // esc set by main calling process
367 if(node->exec & NODE_BREAK)
368 break;
369 }
370 }
371
372 //------------------------------------------------------------------
373 // main loop
374 for (y=0; y<img->y; y++) {
375 // some sort of visual feedback would be nice, or at least this text in the renderwin header
376 // but for now just print some info in the console every 8 scanlines.
377 if (((y & 7)==0) || (y==(img->y-1))) {
378 printf("\rdefocus: Processing Line %d of %d ... ", y+1, img->y);
379 fflush(stdout);
380 }
381 // esc set by main calling process
382 if(node->exec & NODE_BREAK)
383 break;
384
385 zp = y * img->x;
386 for (x=0; x<img->x; x++) {
387 cp = zp + x;
388 cp4 = cp * img->type;
389
390 // Circle of Confusion radius for current pixel
391 cR2 = ct_crad = crad->rect[cp];
392 // skip if zero (border render)
393 if (ct_crad==0.f) {
394 // related to bug #5921, forgot output image when skipping 0 radius values
395 new->rect[cp4] = img->rect[cp4];
396 if (new->type != CB_VAL) {
397 new->rect[cp4+1] = img->rect[cp4+1];
398 new->rect[cp4+2] = img->rect[cp4+2];
399 new->rect[cp4+3] = img->rect[cp4+3];
400 }
401 continue;
402 }
403 cR2 *= cR2;
404
405 // pixel color
406 ctcol = &img->rect[cp4];
407
408 if (!nqd->preview) {
409 int xs, xe, ys, ye;
410 float lwt, wtcol[4] = {0}, aacol[4] = {0};
411
412 // shape weight
413 if (nqd->bktype==0) // disk
414 wt = 1.f/((float)M_PI*cR2);
415 else
416 wt = 1.f/(cR2*bk_hn_theta);
417
418 // weighted color
419 wtcol[0] = wt*ctcol[0];
420 if (new->type != CB_VAL) {
421 wtcol[1] = wt*ctcol[1];
422 wtcol[2] = wt*ctcol[2];
423 wtcol[3] = wt*ctcol[3];
424 }
425
426 // macro for background blur overlap test
427 // unfortunately, since this is done per pixel,
428 // it has a very significant negative impact on processing time...
429 // (eg. aa disk blur without test: 112 sec, vs with test: 176 sec...)
430 // iff center blur radius > threshold
431 // and if overlap pixel in focus, do nothing, else add color/weigbt
432 // (threshold constant is dependant on amount of blur)
433 #define TESTBG1(c, w) {\
434 if (ct_crad > nqd->bthresh) {\
435 if (crad->rect[p] > nqd->bthresh) {\
436 new->rect[p] += c[0];\
437 wts->rect[p] += w;\
438 }\
439 }\
440 else {\
441 new->rect[p] += c[0];\
442 wts->rect[p] += w;\
443 }\
444 }
445 #define TESTBG4(c, w) {\
446 if (ct_crad > nqd->bthresh) {\
447 if (crad->rect[p] > nqd->bthresh) {\
448 new->rect[p4] += c[0];\
449 new->rect[p4+1] += c[1];\
450 new->rect[p4+2] += c[2];\
451 new->rect[p4+3] += c[3];\
452 wts->rect[p] += w;\
453 }\
454 }\
455 else {\
456 new->rect[p4] += c[0];\
457 new->rect[p4+1] += c[1];\
458 new->rect[p4+2] += c[2];\
459 new->rect[p4+3] += c[3];\
460 wts->rect[p] += w;\
461 }\
462 }
463 if (nqd->bktype == 0) {
464 // Disk
465 int _x, i, j, di;
466 float Dj, T;
467 // AA pixel
468 #define AAPIX(a, b) {\
469 int _ny = b;\
470 if ((_ny >= 0) && (_ny < new->y)) {\
471 int _nx = a;\
472 if ((_nx >=0) && (_nx < new->x)) {\
473 p = _ny*new->x + _nx;\
474 if (new->type==CB_VAL) {\
475 TESTBG1(aacol, lwt);\
476 }\
477 else {\
478 p4 = p * new->type;\
479 TESTBG4(aacol, lwt);\
480 }\
481 }\
482 }\
483 }
484 // circle scanline
485 #define CSCAN(a, b) {\
486 int _ny = y + b;\
487 if ((_ny >= 0) && (_ny < new->y)) {\
488 xs = x - a + 1;\
489 if (xs < 0) xs = 0;\
490 xe = x + a;\
491 if (xe > new->x) xe = new->x;\
492 p = _ny*new->x + xs;\
493 if (new->type==CB_VAL) {\
494 for (_x=xs; _x<xe; _x++, p++) TESTBG1(wtcol, wt);\
495 }\
496 else {\
497 p4 = p * new->type;\
498 for (_x=xs; _x<xe; _x++, p++, p4+=new->type) TESTBG4(wtcol, wt);\
499 }\
500 }\
501 }
502 i = ceil(ct_crad);
503 j = 0;
504 T = 0;
505 while (i > j) {
506 Dj = sqrt(cR2 - j*j);
507 Dj -= floor(Dj);
508 di = 0;
509 if (Dj > T) { i--; di = 1; }
510 T = Dj;
511 aacol[0] = wtcol[0]*Dj;
512 if (new->type != CB_VAL) {
513 aacol[1] = wtcol[1]*Dj;
514 aacol[2] = wtcol[2]*Dj;
515 aacol[3] = wtcol[3]*Dj;
516 }
517 lwt = wt*Dj;
518 if (i!=j) {
519 // outer pixels
520 AAPIX(x+j, y+i);
521 AAPIX(x+j, y-i);
522 if (j) {
523 AAPIX(x-j, y+i); // BL
524 AAPIX(x-j, y-i); // TL
525 }
526 if (di) { // only when i changed, interior of outer section
527 CSCAN(j, i); // bottom
528 CSCAN(j, -i); // top
529 }
530 }
531 // lower mid section
532 AAPIX(x+i, y+j);
533 if (i) AAPIX(x-i, y+j);
534 CSCAN(i, j);
535 // upper mid section
536 if (j) {
537 AAPIX(x+i, y-j);
538 if (i) AAPIX(x-i, y-j);
539 CSCAN(i, -j);
540 }
541 j++;
542 }
543 #undef CSCAN
544 #undef AAPIX
545 }
546 else {
547 // n-agonal
548 int ov, nv;
549 float mind, maxd, lwt;
550 ys = MAX2((int)floor(bkh_b[2]*ct_crad + y), 0);
551 ye = MIN2((int)ceil(bkh_b[3]*ct_crad + y), new->y - 1);
552 for (sy=ys; sy<=ye; sy++) {
553 float fxs = 1e10f, fxe = -1e10f;
554 float yf = (sy - y)/ct_crad;
555 int found = 0;
556 ov = len_bkh - 1;
557 mind = maxd = 0;
558 for (nv=0; nv<len_bkh; nv++) {
559 if ((BKH[nv].max_y >= yf) && (BKH[nv].min_y <= yf)) {
560 float tx = BKH[ov].x0 + BKH[nv].ls_x*(yf - BKH[ov].y0);
561 if (tx < fxs) { fxs = tx; mind = BKH[nv].ls_x; }
562 if (tx > fxe) { fxe = tx; maxd = BKH[nv].ls_x; }
563 if (++found == 2) break;
564 }
565 ov = nv;
566 }
567 if (found) {
568 fxs = fxs*ct_crad + x;
569 fxe = fxe*ct_crad + x;
570 xs = (int)floor(fxs), xe = (int)ceil(fxe);
571 // AA hack for first and last x pixel, near vertical edges only
572 if (fabs(mind) <= 1.f) {
573 if ((xs >= 0) && (xs < new->x)) {
574 lwt = 1.f-(fxs - xs);
575 aacol[0] = wtcol[0]*lwt;
576 p = xs + sy*new->x;
577 if (new->type==CB_VAL) {
578 lwt *= wt;
579 TESTBG1(aacol, lwt);
580 }
581 else {
582 p4 = p * new->type;
583 aacol[1] = wtcol[1]*lwt;
584 aacol[2] = wtcol[2]*lwt;
585 aacol[3] = wtcol[3]*lwt;
586 lwt *= wt;
587 TESTBG4(aacol, lwt);
588 }
589 }
590 }
591 if (fabs(maxd) <= 1.f) {
592 if ((xe >= 0) && (xe < new->x)) {
593 lwt = 1.f-(xe - fxe);
594 aacol[0] = wtcol[0]*lwt;
595 p = xe + sy*new->x;
596 if (new->type==CB_VAL) {
597 lwt *= wt;
598 TESTBG1(aacol, lwt);
599 }
600 else {
601 p4 = p * new->type;
602 aacol[1] = wtcol[1]*lwt;
603 aacol[2] = wtcol[2]*lwt;
604 aacol[3] = wtcol[3]*lwt;
605 lwt *= wt;
606 TESTBG4(aacol, lwt);
607 }
608 }
609 }
610 xs = MAX2(xs+1, 0);
611 xe = MIN2(xe, new->x);
612 // remaining interior scanline
613 p = sy*new->x + xs;
614 if (new->type==CB_VAL) {
615 for (sx=xs; sx<xe; sx++, p++) TESTBG1(wtcol, wt);
616 }
617 else {
618 p4 = p * new->type;
619 for (sx=xs; sx<xe; sx++, p++, p4+=new->type) TESTBG4(wtcol, wt);
620 }
621 }
622 }
623
624 // now traverse in opposite direction, y scanlines,
625 // but this time only draw the near horizontal edges,
626 // applying same AA hack as above
627 xs = MAX2((int)floor(bkh_b[0]*ct_crad + x), 0);
628 xe = MIN2((int)ceil(bkh_b[1]*ct_crad + x), img->x - 1);
629 for (sx=xs; sx<=xe; sx++) {
630 float xf = (sx - x)/ct_crad;
631 float fys = 1e10f, fye = -1e10f;
632 int found = 0;
633 ov = len_bkh - 1;
634 mind = maxd = 0;
635 for (nv=0; nv<len_bkh; nv++) {
636 if ((BKH[nv].max_x >= xf) && (BKH[nv].min_x <= xf)) {
637 float ty = BKH[ov].y0 + BKH[nv].ls_y*(xf - BKH[ov].x0);
638 if (ty < fys) { fys = ty; mind = BKH[nv].ls_y; }
639 if (ty > fye) { fye = ty; maxd = BKH[nv].ls_y; }
640 if (++found == 2) break;
641 }
642 ov = nv;
643 }
644 if (found) {
645 fys = fys*ct_crad + y;
646 fye = fye*ct_crad + y;
647 // near horizontal edges only, line slope <= 1
648 if (fabs(mind) <= 1.f) {
649 int iys = (int)floor(fys);
650 if ((iys >= 0) && (iys < new->y)) {
651 lwt = 1.f - (fys - iys);
652 aacol[0] = wtcol[0]*lwt;
653 p = sx + iys*new->x;
654 if (new->type==CB_VAL) {
655 lwt *= wt;
656 TESTBG1(aacol, lwt);
657 }
658 else {
659 p4 = p * new->type;
660 aacol[1] = wtcol[1]*lwt;
661 aacol[2] = wtcol[2]*lwt;
662 aacol[3] = wtcol[3]*lwt;
663 lwt *= wt;
664 TESTBG4(aacol, lwt);
665 }
666 }
667 }
668 if (fabs(maxd) <= 1.f) {
669 int iye = ceil(fye);
670 if ((iye >= 0) && (iye < new->y)) {
671 lwt = 1.f - (iye - fye);
672 aacol[0] = wtcol[0]*lwt;
673 p = sx + iye*new->x;
674 if (new->type==CB_VAL) {
675 lwt *= wt;
676 TESTBG1(aacol, lwt);
677 }
678 else {
679 p4 = p * new->type;
680 aacol[1] = wtcol[1]*lwt;
681 aacol[2] = wtcol[2]*lwt;
682 aacol[3] = wtcol[3]*lwt;
683 lwt *= wt;
684 TESTBG4(aacol, lwt);
685 }
686 }
687 }
688 }
689 }
690
691 }
692 #undef TESTBG4
693 #undef TESTBG1
694
695 }
696 else {
697 // sampled, simple rejection sampling here, good enough
698 unsigned int maxsam, s, ui = BLI_rand()*BLI_rand();
699 float wcor, cpr = BLI_frand();
700 if (no_zbuf)
701 maxsam = nqd->samples; // no zbuffer input, use sample value directly
702 else {
703 // depth adaptive sampling hack, the more out of focus, the more samples taken, 16 minimum.
704 maxsam = (int)(0.5f + nqd->samples*(1.f-(float)exp(-fabs(zbuf->rect[cp] - cam_fdist))));
705 if (maxsam < 16) maxsam = 16;
706 }
707 wcor = 1.f/(float)maxsam;
708 for (s=0; s<maxsam; ++s) {
709 u = ct_crad*(2.f*RI_vdC(s, ui) - 1.f);
710 v = ct_crad*(2.f*(s + cpr)/(float)maxsam - 1.f);
711 sx = (int)(x + u + 0.5f), sy = (int)(y + v + 0.5f);
712 if ((sx<0) || (sx >= new->x) || (sy<0) || (sy >= new->y)) continue;
713 p = sx + sy*new->x;
714 p4 = p * new->type;
715 if (nqd->bktype==0) // Disk
716 lwt = ((u*u + v*v)<=cR2) ? wcor : 0.f;
717 else // AA not needed here
718 lwt = wcor * getWeight(BKH, len_bkh, u, v, ct_crad, inradsq);
719 // prevent background bleeding onto in-focus pixels, user-option
720 if (ct_crad > nqd->bthresh) { // if center blur > threshold
721 if (crad->rect[p] > nqd->bthresh) { // if overlap pixel in focus, do nothing, else add color/weigbt
722 new->rect[p4] += ctcol[0] * lwt;
723 if (new->type != CB_VAL) {
724 new->rect[p4+1] += ctcol[1] * lwt;
725 new->rect[p4+2] += ctcol[2] * lwt;
726 new->rect[p4+3] += ctcol[3] * lwt;
727 }
728 wts->rect[p] += lwt;
729 }
730 }
731 else {
732 new->rect[p4] += ctcol[0] * lwt;
733 if (new->type != CB_VAL) {
734 new->rect[p4+1] += ctcol[1] * lwt;
735 new->rect[p4+2] += ctcol[2] * lwt;
736 new->rect[p4+3] += ctcol[3] * lwt;
737 }
738 wts->rect[p] += lwt;
739 }
740 }
741 }
742
743 }
744 }
745
746 // finally, normalize
747 for (y=0; y<new->y; y++) {
748 p = y * new->x;
749 p4 = p * new->type;
750 for (x=0; x<new->x; x++) {
751 float dv = (wts->rect[p]==0.f) ? 1.f : (1.f/wts->rect[p]);
752 new->rect[p4] *= dv;
753 if (new->type!=CB_VAL) {
754 new->rect[p4+1] *= dv;
755 new->rect[p4+2] *= dv;
756 new->rect[p4+3] *= dv;
757 }
758 p++;
759 p4 += new->type;
760 }
761 }
762
763 free_compbuf(crad);
764 free_compbuf(wts);
765
766 printf("Done\n");
767 }
768
769
770 static void node_composit_exec_defocus(void *data, bNode *node, bNodeStack **in, bNodeStack **out)
771 {
772 CompBuf *new, *old, *zbuf_use = NULL, *img = in[0]->data, *zbuf = in[1]->data;
773 NodeDefocus *nqd = node->storage;
774 int no_zbuf = nqd->no_zbuf;
775
776 if ((img==NULL) || (out[0]->hasoutput==0)) return;
777
778 // if image not valid type or fstop==infinite (128), nothing to do, pass in to out
779 if (((img->type!=CB_RGBA) && (img->type!=CB_VAL)) || ((no_zbuf==0) && (nqd->fstop==128.f))) {
780 out[0]->data = pass_on_compbuf(img);
781 return;
782 }
783
784 if (zbuf!=NULL) {
785 // Zbuf input, check to make sure, single channel, same size
786 // doesn't have to be actual zbuffer, but must be value type
787 if ((zbuf->x != img->x) || (zbuf->y != img->y)) {
788 // could do a scale here instead...
789 printf("Z input must be same size as image !\n");
790 return;
791 }
792 zbuf_use = typecheck_compbuf(zbuf, CB_VAL);
793 }
794 else no_zbuf = 1; // no zbuffer input
795
796 // ok, process
797 old = img;
798 if (nqd->gamco) {
799 // gamma correct, blender func is simplified, fixed value & RGBA only,
800 // should make user param. also depremul and premul afterwards, gamma
801 // correction can't work with premul alpha
802 old = dupalloc_compbuf(img);
803 premul_compbuf(old, 1);
804 gamma_correct_compbuf(old, 0);
805 premul_compbuf(old, 0);
806 }
807
808 new = alloc_compbuf(old->x, old->y, old->type, 1);
809 defocus_blur(node, new, old, zbuf_use, in[1]->vec[0]*nqd->scale, no_zbuf);
810
811 if (nqd->gamco) {
812 premul_compbuf(new, 1);
813 gamma_correct_compbuf(new, 1);
814 premul_compbuf(new, 0);
815 free_compbuf(old);
816 }
817 if(node->exec & NODE_BREAK) {
818 free_compbuf(new);
819 new= NULL;
820 }
821 out[0]->data = new;
822 if (zbuf_use && (zbuf_use != zbuf)) free_compbuf(zbuf_use);
823 }
824
825 static void node_composit_init_defocus(bNode* node)
826 {
827 /* qdn: defocus node */
828 NodeDefocus *nbd = MEM_callocN(sizeof(NodeDefocus), "node defocus data");
829 nbd->bktype = 0;
830 nbd->rotation = 0.f;
831 nbd->preview = 1;
832 nbd->gamco = 0;
833 nbd->samples = 16;
834 nbd->fstop = 128.f;
835 nbd->maxblur = 0;
836 nbd->bthresh = 1.f;
837 nbd->scale = 1.f;
838 nbd->no_zbuf = 1;
839 node->storage = nbd;
840 }
841
842 bNodeType cmp_node_defocus = {
843 /* *next,*prev */ NULL, NULL,
844 /* type code */ CMP_NODE_DEFOCUS,
845 /* name */ "Defocus",
846 /* width+range */ 150, 120, 200,
847 /* class+opts */ NODE_CLASS_OP_FILTER, NODE_OPTIONS,
848 /* input sock */ cmp_node_defocus_in,
849 /* output sock */ cmp_node_defocus_out,
850 /* storage */ "NodeDefocus",
851 /* execfunc */ node_composit_exec_defocus,
852 /* butfunc */ NULL,
853 /* initfunc */ node_composit_init_defocus,
854 /* freestoragefunc */ node_free_standard_storage,
855 /* copystoragefunc */ node_copy_standard_storage,
856 /* id */ NULL
857 };
858
859
Diego's blender repository