Color Management, Stage 2: Switch color pipeline to use OpenColorIO
[blender.git] / source / blender / nodes / composite / node_composite_util.c
1 /*
2  * ***** BEGIN GPL LICENSE BLOCK *****
3  *
4  * This program is free software; you can redistribute it and/or
5  * modify it under the terms of the GNU General Public License
6  * as published by the Free Software Foundation; either version 2
7  * of the License, or (at your option) any later version. 
8  *
9  * This program is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with this program; if not, write to the Free Software Foundation,
16  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
17  *
18  * The Original Code is Copyright (C) 2006 Blender Foundation.
19  * All rights reserved.
20  *
21  * The Original Code is: all of this file.
22  *
23  * Contributor(s): none yet.
24  *
25  * ***** END GPL LICENSE BLOCK *****
26  */
27
28 /** \file blender/nodes/composite/node_composite_util.c
29  *  \ingroup nodes
30  */
31
32 #include "node_composite_util.h"
33
34 #ifdef WITH_COMPOSITOR_LEGACY
35
36 #include <limits.h>
37
38 CompBuf *alloc_compbuf(int sizex, int sizey, int type, int alloc)
39 {
40         CompBuf *cbuf= MEM_callocN(sizeof(CompBuf), "compbuf");
41         
42         cbuf->x= sizex;
43         cbuf->y= sizey;
44         cbuf->xrad= sizex/2;
45         cbuf->yrad= sizey/2;
46         
47         cbuf->type= type;
48         if (alloc) {
49                 if (cbuf->type==CB_RGBA)
50                         cbuf->rect= MEM_mapallocN(4*sizeof(float)*sizex*sizey, "compbuf RGBA rect");
51                 else if (cbuf->type==CB_VEC3)
52                         cbuf->rect= MEM_mapallocN(3*sizeof(float)*sizex*sizey, "compbuf Vector3 rect");
53                 else if (cbuf->type==CB_VEC2)
54                         cbuf->rect= MEM_mapallocN(2*sizeof(float)*sizex*sizey, "compbuf Vector2 rect");
55                 else
56                         cbuf->rect= MEM_mapallocN(sizeof(float)*sizex*sizey, "compbuf Fac rect");
57                 cbuf->malloc= 1;
58         }
59         cbuf->disprect.xmin = 0;
60         cbuf->disprect.ymin = 0;
61         cbuf->disprect.xmax = sizex;
62         cbuf->disprect.ymax = sizey;
63         
64         return cbuf;
65 }
66
67 CompBuf *dupalloc_compbuf(CompBuf *cbuf)
68 {
69         CompBuf *dupbuf= alloc_compbuf(cbuf->x, cbuf->y, cbuf->type, 1);
70         if (dupbuf) {
71                 memcpy(dupbuf->rect, cbuf->rect, cbuf->type*sizeof(float)*cbuf->x*cbuf->y);
72         
73                 dupbuf->xof= cbuf->xof;
74                 dupbuf->yof= cbuf->yof;
75         }       
76         return dupbuf;
77 }
78
79 /* instead of reference counting, we create a list */
80 CompBuf *pass_on_compbuf(CompBuf *cbuf)
81 {
82         CompBuf *dupbuf= (cbuf)? alloc_compbuf(cbuf->x, cbuf->y, cbuf->type, 0): NULL;
83         CompBuf *lastbuf;
84         
85         if (dupbuf) {
86                 dupbuf->rect= cbuf->rect;
87                 dupbuf->xof= cbuf->xof;
88                 dupbuf->yof= cbuf->yof;
89                 dupbuf->malloc= 0;
90                 
91                 /* get last buffer in list, and append dupbuf */
92                 for (lastbuf= cbuf; lastbuf; lastbuf= lastbuf->next)
93                         if (lastbuf->next==NULL)
94                                 break;
95                 lastbuf->next= dupbuf;
96                 dupbuf->prev= lastbuf;
97         }       
98         return dupbuf;
99 }
100
101
102 void free_compbuf(CompBuf *cbuf)
103 {
104         /* check referencing, then remove from list and set malloc tag */
105         if (cbuf->prev || cbuf->next) {
106                 if (cbuf->prev)
107                         cbuf->prev->next= cbuf->next;
108                 if (cbuf->next)
109                         cbuf->next->prev= cbuf->prev;
110                 if (cbuf->malloc) {
111                         if (cbuf->prev)
112                                 cbuf->prev->malloc= 1;
113                         else
114                                 cbuf->next->malloc= 1;
115                         cbuf->malloc= 0;
116                 }
117         }
118         
119         if (cbuf->malloc && cbuf->rect)
120                 MEM_freeN(cbuf->rect);
121
122         MEM_freeN(cbuf);
123 }
124
125 void print_compbuf(char *str, CompBuf *cbuf)
126 {
127         printf("Compbuf %s %d %d %p\n", str, cbuf->x, cbuf->y, (void *)cbuf->rect);
128         
129 }
130
131 void compbuf_set_node(CompBuf *cbuf, bNode *node)
132 {
133         if (cbuf) cbuf->node = node;
134 }
135
136
137 CompBuf *get_cropped_compbuf(rcti *drect, float *rectf, int rectx, int recty, int type)
138 {
139         CompBuf *cbuf;
140         rcti disprect= *drect;
141         float *outfp;
142         int dx, y;
143         
144         if (disprect.xmax>rectx) disprect.xmax = rectx;
145         if (disprect.ymax>recty) disprect.ymax = recty;
146         if (disprect.xmin>= disprect.xmax) return NULL;
147         if (disprect.ymin>= disprect.ymax) return NULL;
148         
149         cbuf= alloc_compbuf(BLI_RCT_SIZE_X(&disprect), BLI_RCT_SIZE_Y(&disprect), type, 1);
150         outfp= cbuf->rect;
151         rectf += type*(disprect.ymin*rectx + disprect.xmin);
152         dx= type*cbuf->x;
153         for (y=cbuf->y; y>0; y--, outfp+=dx, rectf+=type*rectx)
154                 memcpy(outfp, rectf, sizeof(float)*dx);
155         
156         return cbuf;
157 }
158
159 CompBuf *scalefast_compbuf(CompBuf *inbuf, int newx, int newy)
160 {
161         CompBuf *outbuf; 
162         float *rectf, *newrectf, *rf;
163         int x, y, c, pixsize= inbuf->type;
164         int ofsx, ofsy, stepx, stepy;
165         
166         if (inbuf->x==newx && inbuf->y==newy)
167                 return dupalloc_compbuf(inbuf);
168         
169         outbuf= alloc_compbuf(newx, newy, inbuf->type, 1);
170         newrectf= outbuf->rect;
171         
172         stepx = (65536.0 * (inbuf->x - 1.0) / (newx - 1.0)) + 0.5;
173         stepy = (65536.0 * (inbuf->y - 1.0) / (newy - 1.0)) + 0.5;
174         ofsy = 32768;
175         
176         for (y = newy; y > 0 ; y--) {
177                 rectf = inbuf->rect;
178                 rectf += pixsize * (ofsy >> 16) * inbuf->x;
179
180                 ofsy += stepy;
181                 ofsx = 32768;
182                 
183                 for (x = newx ; x>0 ; x--) {
184                         
185                         rf= rectf + pixsize*(ofsx >> 16);
186                         for (c=0; c<pixsize; c++)
187                                 newrectf[c] = rf[c];
188                         
189                         newrectf+= pixsize;
190                         
191                         ofsx += stepx;
192                 }
193         }
194         
195         return outbuf;
196 }
197
198 void typecheck_compbuf_color(float *out, float *in, int outtype, int intype)
199 {
200         if (intype == outtype) {
201                 memcpy(out, in, sizeof(float)*outtype);
202         }
203         else if (outtype==CB_VAL) {
204                 if (intype==CB_VEC2) {
205                         *out= 0.5f*(in[0]+in[1]);
206                 }
207                 else if (intype==CB_VEC3) {
208                         *out= 0.333333f*(in[0]+in[1]+in[2]);
209                 }
210                 else if (intype==CB_RGBA) {
211                         *out = rgb_to_bw(in);
212                 }
213         }
214         else if (outtype==CB_VEC2) {
215                 if (intype==CB_VAL) {
216                         out[0]= in[0];
217                         out[1]= in[0];
218                 }
219                 else if (intype==CB_VEC3) {
220                         out[0]= in[0];
221                         out[1]= in[1];
222                 }
223                 else if (intype==CB_RGBA) {
224                         out[0]= in[0];
225                         out[1]= in[1];
226                 }
227         }
228         else if (outtype==CB_VEC3) {
229                 if (intype==CB_VAL) {
230                         out[0]= in[0];
231                         out[1]= in[0];
232                         out[2]= in[0];
233                 }
234                 else if (intype==CB_VEC2) {
235                         out[0]= in[0];
236                         out[1]= in[1];
237                         out[2]= 0.0f;
238                 }
239                 else if (intype==CB_RGBA) {
240                         out[0]= in[0];
241                         out[1]= in[1];
242                         out[2]= in[2];
243                 }
244         }
245         else if (outtype==CB_RGBA) {
246                 if (intype==CB_VAL) {
247                         out[0]= in[0];
248                         out[1]= in[0];
249                         out[2]= in[0];
250                         out[3]= 1.0f;
251                 }
252                 else if (intype==CB_VEC2) {
253                         out[0]= in[0];
254                         out[1]= in[1];
255                         out[2]= 0.0f;
256                         out[3]= 1.0f;
257                 }
258                 else if (intype==CB_VEC3) {
259                         out[0]= in[0];
260                         out[1]= in[1];
261                         out[2]= in[2];
262                         out[3]= 1.0f;
263                 }
264         }
265 }
266
267 CompBuf *typecheck_compbuf(CompBuf *inbuf, int type)
268 {
269         if (inbuf && inbuf->type!=type) {
270                 CompBuf *outbuf;
271                 float *inrf, *outrf;
272                 int x;
273
274                 outbuf= alloc_compbuf(inbuf->x, inbuf->y, type, 1); 
275
276                 /* warning note: xof and yof are applied in pixelprocessor, but should be copied otherwise? */
277                 outbuf->xof= inbuf->xof;
278                 outbuf->yof= inbuf->yof;
279
280                 if (inbuf->rect_procedural) {
281                         outbuf->rect_procedural= inbuf->rect_procedural;
282                         copy_v3_v3(outbuf->procedural_size, inbuf->procedural_size);
283                         copy_v3_v3(outbuf->procedural_offset, inbuf->procedural_offset);
284                         outbuf->procedural_type= inbuf->procedural_type;
285                         outbuf->node= inbuf->node;
286                         return outbuf;
287                 }
288
289                 inrf= inbuf->rect;
290                 outrf= outbuf->rect;
291                 x= inbuf->x*inbuf->y;
292                 
293                 if (type==CB_VAL) {
294                         if (inbuf->type==CB_VEC2) {
295                                 for (; x>0; x--, outrf+= 1, inrf+= 2)
296                                         *outrf= 0.5f*(inrf[0]+inrf[1]);
297                         }
298                         else if (inbuf->type==CB_VEC3) {
299                                 for (; x>0; x--, outrf+= 1, inrf+= 3)
300                                         *outrf= 0.333333f*(inrf[0]+inrf[1]+inrf[2]);
301                         }
302                         else if (inbuf->type==CB_RGBA) {
303                                 for (; x>0; x--, outrf+= 1, inrf+= 4)
304                                         *outrf = rgb_to_bw(inrf);
305                         }
306                 }
307                 else if (type==CB_VEC2) {
308                         if (inbuf->type==CB_VAL) {
309                                 for (; x>0; x--, outrf+= 2, inrf+= 1) {
310                                         outrf[0]= inrf[0];
311                                         outrf[1]= inrf[0];
312                                 }
313                         }
314                         else if (inbuf->type==CB_VEC3) {
315                                 for (; x>0; x--, outrf+= 2, inrf+= 3) {
316                                         outrf[0]= inrf[0];
317                                         outrf[1]= inrf[1];
318                                 }
319                         }
320                         else if (inbuf->type==CB_RGBA) {
321                                 for (; x>0; x--, outrf+= 2, inrf+= 4) {
322                                         outrf[0]= inrf[0];
323                                         outrf[1]= inrf[1];
324                                 }
325                         }
326                 }
327                 else if (type==CB_VEC3) {
328                         if (inbuf->type==CB_VAL) {
329                                 for (; x>0; x--, outrf+= 3, inrf+= 1) {
330                                         outrf[0]= inrf[0];
331                                         outrf[1]= inrf[0];
332                                         outrf[2]= inrf[0];
333                                 }
334                         }
335                         else if (inbuf->type==CB_VEC2) {
336                                 for (; x>0; x--, outrf+= 3, inrf+= 2) {
337                                         outrf[0]= inrf[0];
338                                         outrf[1]= inrf[1];
339                                         outrf[2]= 0.0f;
340                                 }
341                         }
342                         else if (inbuf->type==CB_RGBA) {
343                                 for (; x>0; x--, outrf+= 3, inrf+= 4) {
344                                         outrf[0]= inrf[0];
345                                         outrf[1]= inrf[1];
346                                         outrf[2]= inrf[2];
347                                 }
348                         }
349                 }
350                 else if (type==CB_RGBA) {
351                         if (inbuf->type==CB_VAL) {
352                                 for (; x>0; x--, outrf+= 4, inrf+= 1) {
353                                         outrf[0]= inrf[0];
354                                         outrf[1]= inrf[0];
355                                         outrf[2]= inrf[0];
356                                         outrf[3]= 1.0f;
357                                 }
358                         }
359                         else if (inbuf->type==CB_VEC2) {
360                                 for (; x>0; x--, outrf+= 4, inrf+= 2) {
361                                         outrf[0]= inrf[0];
362                                         outrf[1]= inrf[1];
363                                         outrf[2]= 0.0f;
364                                         outrf[3]= 1.0f;
365                                 }
366                         }
367                         else if (inbuf->type==CB_VEC3) {
368                                 for (; x>0; x--, outrf+= 4, inrf+= 3) {
369                                         outrf[0]= inrf[0];
370                                         outrf[1]= inrf[1];
371                                         outrf[2]= inrf[2];
372                                         outrf[3]= 1.0f;
373                                 }
374                         }
375                 }
376                 
377                 return outbuf;
378         }
379         return inbuf;
380 }
381
382 float *compbuf_get_pixel(CompBuf *cbuf, float *defcol, float *use, int x, int y, int xrad, int yrad)
383 {
384         if (cbuf) {
385                 if (cbuf->rect_procedural) {
386                         cbuf->rect_procedural(cbuf, use, (float)x/(float)xrad, (float)y/(float)yrad);
387                         return use;
388                 }
389                 else {
390                         static float col[4]= {0.0f, 0.0f, 0.0f, 0.0f};
391                         
392                         /* map coords */
393                         x-= cbuf->xof;
394                         y-= cbuf->yof;
395                         
396                         if (y<-cbuf->yrad || y>= -cbuf->yrad+cbuf->y) return col;
397                         if (x<-cbuf->xrad || x>= -cbuf->xrad+cbuf->x) return col;
398                         
399                         return cbuf->rect + cbuf->type*( (cbuf->yrad+y)*cbuf->x + (cbuf->xrad+x) );
400                 }
401         }
402         else return defcol;
403 }
404
405 /* **************************************************** */
406
407 static CompBuf *composit_check_compbuf(CompBuf *cbuf, int type, CompBuf *outbuf)
408 {
409         /* check type */
410         CompBuf *dbuf= typecheck_compbuf(cbuf, type);
411
412         /* if same as output and translated, duplicate so pixels don't interfere */
413         if (dbuf == outbuf && !dbuf->rect_procedural && (dbuf->xof || dbuf->yof))
414                 dbuf= dupalloc_compbuf(dbuf);
415         
416         return dbuf;
417 }
418
419 /* Pixel-to-Pixel operation, 1 Image in, 1 out */
420 void composit1_pixel_processor(bNode *node, CompBuf *out, CompBuf *src_buf, float *src_col,
421                                                                           void (*func)(bNode *, float *, float *), 
422                                                                           int src_type)
423 {
424         CompBuf *src_use;
425         float *outfp=out->rect, *srcfp;
426         float color[4]; /* local color if compbuf is procedural */
427         int xrad, yrad, x, y;
428         
429         src_use= composit_check_compbuf(src_buf, src_type, out);
430         
431         xrad= out->xrad;
432         yrad= out->yrad;
433         
434         for (y= -yrad; y<-yrad+out->y; y++) {
435                 for (x= -xrad; x<-xrad+out->x; x++, outfp+=out->type) {
436                         srcfp= compbuf_get_pixel(src_use, src_col, color, x, y, xrad, yrad);
437                         func(node, outfp, srcfp);
438                 }
439         }
440         
441         if (src_use!=src_buf)
442                 free_compbuf(src_use);
443 }
444
445 /* Pixel-to-Pixel operation, 2 Images in, 1 out */
446 void composit2_pixel_processor(bNode *node, CompBuf *out, CompBuf *src_buf, float *src_col,
447                                                                           CompBuf *fac_buf, float *fac, void (*func)(bNode *, float *, float *, float *), 
448                                                                           int src_type, int fac_type)
449 {
450         CompBuf *src_use, *fac_use;
451         float *outfp=out->rect, *srcfp, *facfp;
452         float color[4]; /* local color if compbuf is procedural */
453         int xrad, yrad, x, y;
454         
455         src_use= composit_check_compbuf(src_buf, src_type, out);
456         fac_use= composit_check_compbuf(fac_buf, fac_type, out);
457
458         xrad= out->xrad;
459         yrad= out->yrad;
460         
461         for (y= -yrad; y<-yrad+out->y; y++) {
462                 for (x= -xrad; x<-xrad+out->x; x++, outfp+=out->type) {
463                         srcfp= compbuf_get_pixel(src_use, src_col, color, x, y, xrad, yrad);
464                         facfp= compbuf_get_pixel(fac_use, fac, color, x, y, xrad, yrad);
465                         
466                         func(node, outfp, srcfp, facfp);
467                 }
468         }
469         if (src_use!=src_buf)
470                 free_compbuf(src_use);
471         if (fac_use!=fac_buf)
472                 free_compbuf(fac_use);
473 }
474
475 /* Pixel-to-Pixel operation, 3 Images in, 1 out */
476 void composit3_pixel_processor(bNode *node, CompBuf *out, CompBuf *src1_buf, float *src1_col, CompBuf *src2_buf, float *src2_col, 
477                                                                           CompBuf *fac_buf, float *fac, void (*func)(bNode *, float *, float *, float *, float *), 
478                                                                           int src1_type, int src2_type, int fac_type)
479 {
480         CompBuf *src1_use, *src2_use, *fac_use;
481         float *outfp=out->rect, *src1fp, *src2fp, *facfp;
482         float color[4]; /* local color if compbuf is procedural */
483         int xrad, yrad, x, y;
484         
485         src1_use= composit_check_compbuf(src1_buf, src1_type, out);
486         src2_use= composit_check_compbuf(src2_buf, src2_type, out);
487         fac_use= composit_check_compbuf(fac_buf, fac_type, out);
488         
489         xrad= out->xrad;
490         yrad= out->yrad;
491         
492         for (y= -yrad; y<-yrad+out->y; y++) {
493                 for (x= -xrad; x<-xrad+out->x; x++, outfp+=out->type) {
494                         src1fp= compbuf_get_pixel(src1_use, src1_col, color, x, y, xrad, yrad);
495                         src2fp= compbuf_get_pixel(src2_use, src2_col, color, x, y, xrad, yrad);
496                         facfp= compbuf_get_pixel(fac_use, fac, color, x, y, xrad, yrad);
497                         
498                         func(node, outfp, src1fp, src2fp, facfp);
499                 }
500         }
501         
502         if (src1_use!=src1_buf)
503                 free_compbuf(src1_use);
504         if (src2_use!=src2_buf)
505                 free_compbuf(src2_use);
506         if (fac_use!=fac_buf)
507                 free_compbuf(fac_use);
508 }
509
510 /* Pixel-to-Pixel operation, 4 Images in, 1 out */
511 void composit4_pixel_processor(bNode *node, CompBuf *out, CompBuf *src1_buf, float *src1_col, CompBuf *fac1_buf, float *fac1, 
512                                                                           CompBuf *src2_buf, float *src2_col, CompBuf *fac2_buf, float *fac2, 
513                                                                           void (*func)(bNode *, float *, float *, float *, float *, float *), 
514                                                                           int src1_type, int fac1_type, int src2_type, int fac2_type)
515 {
516         CompBuf *src1_use, *src2_use, *fac1_use, *fac2_use;
517         float *outfp=out->rect, *src1fp, *src2fp, *fac1fp, *fac2fp;
518         float color[4]; /* local color if compbuf is procedural */
519         int xrad, yrad, x, y;
520         
521         src1_use= composit_check_compbuf(src1_buf, src1_type, out);
522         src2_use= composit_check_compbuf(src2_buf, src2_type, out);
523         fac1_use= composit_check_compbuf(fac1_buf, fac1_type, out);
524         fac2_use= composit_check_compbuf(fac2_buf, fac2_type, out);
525         
526         xrad= out->xrad;
527         yrad= out->yrad;
528         
529         for (y= -yrad; y<-yrad+out->y; y++) {
530                 for (x= -xrad; x<-xrad+out->x; x++, outfp+=out->type) {
531                         src1fp= compbuf_get_pixel(src1_use, src1_col, color, x, y, xrad, yrad);
532                         src2fp= compbuf_get_pixel(src2_use, src2_col, color, x, y, xrad, yrad);
533                         fac1fp= compbuf_get_pixel(fac1_use, fac1, color, x, y, xrad, yrad);
534                         fac2fp= compbuf_get_pixel(fac2_use, fac2, color, x, y, xrad, yrad);
535                         
536                         func(node, outfp, src1fp, fac1fp, src2fp, fac2fp);
537                 }
538         }
539         
540         if (src1_use!=src1_buf)
541                 free_compbuf(src1_use);
542         if (src2_use!=src2_buf)
543                 free_compbuf(src2_use);
544         if (fac1_use!=fac1_buf)
545                 free_compbuf(fac1_use);
546         if (fac2_use!=fac2_buf)
547                 free_compbuf(fac2_use);
548 }
549
550
551 CompBuf *valbuf_from_rgbabuf(CompBuf *cbuf, int channel)
552 {
553         CompBuf *valbuf= alloc_compbuf(cbuf->x, cbuf->y, CB_VAL, 1);
554         float *valf, *rectf;
555         int tot;
556         
557         /* warning note: xof and yof are applied in pixelprocessor, but should be copied otherwise? */
558         valbuf->xof= cbuf->xof;
559         valbuf->yof= cbuf->yof;
560         
561         valf= valbuf->rect;
562
563         /* defaults to returning alpha channel */
564         if ((channel < CHAN_R) || (channel > CHAN_A)) channel = CHAN_A;
565
566         rectf= cbuf->rect + channel;
567         
568         for (tot= cbuf->x*cbuf->y; tot>0; tot--, valf++, rectf+=4)
569                 *valf= *rectf;
570         
571         return valbuf;
572 }
573
574 void valbuf_to_rgbabuf(CompBuf *valbuf, CompBuf *cbuf, int channel)
575 {
576         float *valf, *rectf;
577         int tot;
578
579         valf= valbuf->rect;
580
581         /* defaults to returning alpha channel */
582         if ((channel < CHAN_R) || (channel > CHAN_A)) channel = CHAN_A;
583
584         rectf = cbuf->rect + channel;
585
586         for (tot= cbuf->x*cbuf->y; tot>0; tot--, valf++, rectf+=4)
587                 *rectf = *valf;
588 }
589
590 static CompBuf *generate_procedural_preview(CompBuf *cbuf, int newx, int newy)
591 {
592         CompBuf *outbuf;
593         float *outfp;
594         int xrad, yrad, x, y;
595         
596         outbuf= alloc_compbuf(newx, newy, CB_RGBA, 1);
597
598         outfp= outbuf->rect;
599         xrad= outbuf->xrad;
600         yrad= outbuf->yrad;
601         
602         for (y= -yrad; y<-yrad+outbuf->y; y++)
603                 for (x= -xrad; x<-xrad+outbuf->x; x++, outfp+=outbuf->type)
604                         cbuf->rect_procedural(cbuf, outfp, (float)x/(float)xrad, (float)y/(float)yrad);
605
606         return outbuf;
607 }
608
609 /* OCIO_TODO: this function is only used by legacy compositor system only, which would likely be removed soon,
610  *            keep check for old color management flag for now
611  */
612 void generate_preview(void *data, bNode *node, CompBuf *stackbuf)
613 {
614         RenderData *rd= data;
615         bNodePreview *preview= node->preview;
616         int xsize, ysize;
617         int profile_from= (rd->color_mgt_flag & R_COLOR_MANAGEMENT)? IB_PROFILE_LINEAR_RGB: IB_PROFILE_SRGB;
618         int predivide= (rd->color_mgt_flag & R_COLOR_MANAGEMENT_PREDIVIDE);
619         int dither= 0;
620         unsigned char *rect;
621         
622         if (preview && stackbuf) {
623                 CompBuf *cbuf, *stackbuf_use;
624                 
625                 if (stackbuf->rect==NULL && stackbuf->rect_procedural==NULL) return;
626                 
627                 stackbuf_use= typecheck_compbuf(stackbuf, CB_RGBA);
628
629                 if (stackbuf->x > stackbuf->y) {
630                         xsize= 140;
631                         ysize= (140*stackbuf->y)/stackbuf->x;
632                 }
633                 else {
634                         ysize= 140;
635                         xsize= (140*stackbuf->x)/stackbuf->y;
636                 }
637                 
638                 if (stackbuf_use->rect_procedural)
639                         cbuf= generate_procedural_preview(stackbuf_use, xsize, ysize);
640                 else
641                         cbuf= scalefast_compbuf(stackbuf_use, xsize, ysize);
642
643                 /* convert to byte for preview */
644                 rect= MEM_callocN(sizeof(unsigned char)*4*xsize*ysize, "bNodePreview.rect");
645
646                 IMB_buffer_byte_from_float(rect, cbuf->rect,
647                         4, dither, IB_PROFILE_SRGB, profile_from, predivide, 
648                         xsize, ysize, xsize, xsize);
649                 
650                 free_compbuf(cbuf);
651                 if (stackbuf_use!=stackbuf)
652                         free_compbuf(stackbuf_use);
653
654                 // BLI_lock_thread(LOCK_PREVIEW);
655
656                 if (preview->rect)
657                         MEM_freeN(preview->rect);
658                 preview->xsize= xsize;
659                 preview->ysize= ysize;
660                 preview->rect= rect;
661
662                 // BLI_unlock_thread(LOCK_PREVIEW);
663         }
664 }
665
666 void do_rgba_to_yuva(bNode *UNUSED(node), float *out, float *in)
667 {
668         rgb_to_yuv(in[0], in[1], in[2], &out[0], &out[1], &out[2]);
669         out[3]=in[3];
670 }
671
672 void do_rgba_to_hsva(bNode *UNUSED(node), float *out, float *in)
673 {
674         rgb_to_hsv(in[0], in[1], in[2], &out[0], &out[1], &out[2]);
675         out[3]=in[3];
676 }
677
678 void do_rgba_to_ycca(bNode *UNUSED(node), float *out, float *in)
679 {
680         rgb_to_ycc(in[0], in[1], in[2], &out[0], &out[1], &out[2], BLI_YCC_ITU_BT601);
681         out[3]=in[3];
682 }
683
684 void do_yuva_to_rgba(bNode *UNUSED(node), float *out, float *in)
685 {
686         yuv_to_rgb(in[0], in[1], in[2], &out[0], &out[1], &out[2]);
687         out[3]=in[3];
688 }
689
690 void do_hsva_to_rgba(bNode *UNUSED(node), float *out, float *in)
691 {
692         hsv_to_rgb(in[0], in[1], in[2], &out[0], &out[1], &out[2]);
693         out[3]=in[3];
694 }
695
696 void do_ycca_to_rgba(bNode *UNUSED(node), float *out, float *in)
697 {
698         ycc_to_rgb(in[0], in[1], in[2], &out[0], &out[1], &out[2], BLI_YCC_ITU_BT601);
699         out[3]=in[3];
700 }
701
702 void do_copy_rgba(bNode *UNUSED(node), float *out, float *in)
703 {
704         copy_v4_v4(out, in);
705 }
706
707 void do_copy_rgb(bNode *UNUSED(node), float *out, float *in)
708 {
709         copy_v3_v3(out, in);
710         out[3]= 1.0f;
711 }
712
713 void do_copy_value(bNode *UNUSED(node), float *out, float *in)
714 {
715         out[0]= in[0];
716 }
717
718 void do_copy_a_rgba(bNode *UNUSED(node), float *out, float *in, float *fac)
719 {
720         copy_v3_v3(out, in);
721         out[3]= *fac;
722 }
723
724 /* only accepts RGBA buffers */
725 void gamma_correct_compbuf(CompBuf *img, int inversed)
726 {
727         float *drect;
728         int x;
729
730         if (img->type!=CB_RGBA) return;
731
732         drect= img->rect;
733         if (inversed) {
734                 for (x=img->x*img->y; x>0; x--, drect+=4) {
735                         if (drect[0]>0.0f) drect[0]= sqrt(drect[0]); else drect[0]= 0.0f;
736                         if (drect[1]>0.0f) drect[1]= sqrt(drect[1]); else drect[1]= 0.0f;
737                         if (drect[2]>0.0f) drect[2]= sqrt(drect[2]); else drect[2]= 0.0f;
738                 }
739         }
740         else {
741                 for (x=img->x*img->y; x>0; x--, drect+=4) {
742                         if (drect[0]>0.0f) drect[0]*= drect[0]; else drect[0]= 0.0f;
743                         if (drect[1]>0.0f) drect[1]*= drect[1]; else drect[1]= 0.0f;
744                         if (drect[2]>0.0f) drect[2]*= drect[2]; else drect[2]= 0.0f;
745                 }
746         }
747 }
748
749 void premul_compbuf(CompBuf *img, int inversed)
750 {
751         float *drect;
752         int x;
753
754         if (img->type!=CB_RGBA) return;
755
756         drect= img->rect;
757         if (inversed) {
758                 for (x=img->x*img->y; x>0; x--, drect+=4) {
759                         if (fabsf(drect[3]) < 1e-5f) {
760                                 drect[0]= 0.0f;
761                                 drect[1]= 0.0f;
762                                 drect[2]= 0.0f;
763                         }
764                         else {
765                                 drect[0] /= drect[3];
766                                 drect[1] /= drect[3];
767                                 drect[2] /= drect[3];
768                         }
769                 }
770         }
771         else {
772                 for (x=img->x*img->y; x>0; x--, drect+=4) {
773                         drect[0] *= drect[3];
774                         drect[1] *= drect[3];
775                         drect[2] *= drect[3];
776                 }
777         }
778 }
779
780
781
782 /*
783  *  2D Fast Hartley Transform, used for convolution
784  */
785
786 typedef float fREAL;
787
788 // returns next highest power of 2 of x, as well it's log2 in L2
789 static unsigned int nextPow2(unsigned int x, unsigned int* L2)
790 {
791         unsigned int pw, x_notpow2 = x & (x-1);
792         *L2 = 0;
793         while (x>>=1) ++(*L2);
794         pw = 1 << (*L2);
795         if (x_notpow2) { (*L2)++;  pw<<=1; }
796         return pw;
797 }
798
799 //------------------------------------------------------------------------------
800
801 // from FXT library by Joerg Arndt, faster in order bitreversal
802 // use: r = revbin_upd(r, h) where h = N>>1
803 static unsigned int revbin_upd(unsigned int r, unsigned int h)
804 {
805         while (!((r^=h)&h)) h >>= 1;
806         return r;
807 }
808 //------------------------------------------------------------------------------
809 static void FHT(fREAL* data, unsigned int M, unsigned int inverse)
810 {
811         double tt, fc, dc, fs, ds, a = M_PI;
812         fREAL t1, t2;
813         int n2, bd, bl, istep, k, len = 1 << M, n = 1;
814
815         int i, j = 0;
816         unsigned int Nh = len >> 1;
817         for (i=1;i<(len-1);++i) {
818                 j = revbin_upd(j, Nh);
819                 if (j>i) {
820                         t1 = data[i];
821                         data[i] = data[j];
822                         data[j] = t1;
823                 }
824         }
825
826         do {
827                 fREAL* data_n = &data[n];
828
829                 istep = n << 1;
830                 for (k=0; k<len; k+=istep) {
831                         t1 = data_n[k];
832                         data_n[k] = data[k] - t1;
833                         data[k] += t1;
834                 }
835
836                 n2 = n >> 1;
837                 if (n>2) {
838                         fc = dc = cos(a);
839                         fs = ds = sqrt(1.0 - fc*fc); //sin(a);
840                         bd = n-2;
841                         for (bl=1; bl<n2; bl++) {
842                                 fREAL* data_nbd = &data_n[bd];
843                                 fREAL* data_bd = &data[bd];
844                                 for (k=bl; k<len; k+=istep) {
845                                         t1 = fc*data_n[k] + fs*data_nbd[k];
846                                         t2 = fs*data_n[k] - fc*data_nbd[k];
847                                         data_n[k] = data[k] - t1;
848                                         data_nbd[k] = data_bd[k] - t2;
849                                         data[k] += t1;
850                                         data_bd[k] += t2;
851                                 }
852                                 tt = fc*dc - fs*ds;
853                                 fs = fs*dc + fc*ds;
854                                 fc = tt;
855                                 bd -= 2;
856                         }
857                 }
858
859                 if (n>1) {
860                         for (k=n2; k<len; k+=istep) {
861                                 t1 = data_n[k];
862                                 data_n[k] = data[k] - t1;
863                                 data[k] += t1;
864                         }
865                 }
866
867                 n = istep;
868                 a *= 0.5;
869         } while (n<len);
870
871         if (inverse) {
872                 fREAL sc = (fREAL)1 / (fREAL)len;
873                 for (k=0; k<len; ++k)
874                         data[k] *= sc;
875         }
876 }
877 //------------------------------------------------------------------------------
878 /* 2D Fast Hartley Transform, Mx/My -> log2 of width/height,
879  * nzp -> the row where zero pad data starts,
880  * inverse -> see above */
881 static void FHT2D(fREAL *data, unsigned int Mx, unsigned int My,
882                 unsigned int nzp, unsigned int inverse)
883 {
884         unsigned int i, j, Nx, Ny, maxy;
885         fREAL t;
886
887         Nx = 1 << Mx;
888         Ny = 1 << My;
889
890         // rows (forward transform skips 0 pad data)
891         maxy = inverse ? Ny : nzp;
892         for (j=0; j<maxy; ++j)
893                 FHT(&data[Nx*j], Mx, inverse);
894
895         // transpose data
896         if (Nx==Ny) {  // square
897                 for (j=0; j<Ny; ++j)
898                         for (i=j+1; i<Nx; ++i) {
899                                 unsigned int op = i + (j << Mx), np = j + (i << My);
900                                 t=data[op], data[op]=data[np], data[np]=t;
901                         }
902         }
903         else {  // rectangular
904                 unsigned int k, Nym = Ny-1, stm = 1 << (Mx + My);
905                 for (i=0; stm>0; i++) {
906                         #define PRED(k) (((k & Nym) << Mx) + (k >> My))
907                         for (j=PRED(i); j>i; j=PRED(j));
908                         if (j < i) continue;
909                         for (k=i, j=PRED(i); j!=i; k=j, j=PRED(j), stm--) {
910                                 t=data[j], data[j]=data[k], data[k]=t;
911                         }
912                         #undef PRED
913                         stm--;
914                 }
915         }
916         // swap Mx/My & Nx/Ny
917         i = Nx, Nx = Ny, Ny = i;
918         i = Mx, Mx = My, My = i;
919
920         // now columns == transposed rows
921         for (j=0; j<Ny; ++j)
922                 FHT(&data[Nx*j], Mx, inverse);
923
924         // finalize
925         for (j=0; j<=(Ny >> 1); j++) {
926                 unsigned int jm = (Ny - j) & (Ny-1);
927                 unsigned int ji = j << Mx;
928                 unsigned int jmi = jm << Mx;
929                 for (i=0; i<=(Nx >> 1); i++) {
930                         unsigned int im = (Nx - i) & (Nx-1);
931                         fREAL A = data[ji + i];
932                         fREAL B = data[jmi + i];
933                         fREAL C = data[ji + im];
934                         fREAL D = data[jmi + im];
935                         fREAL E = (fREAL)0.5*((A + D) - (B + C));
936                         data[ji + i] = A - E;
937                         data[jmi + i] = B + E;
938                         data[ji + im] = C + E;
939                         data[jmi + im] = D - E;
940                 }
941         }
942
943 }
944
945 //------------------------------------------------------------------------------
946
947 /* 2D convolution calc, d1 *= d2, M/N - > log2 of width/height */
948 static void fht_convolve(fREAL* d1, fREAL* d2, unsigned int M, unsigned int N)
949 {
950         fREAL a, b;
951         unsigned int i, j, k, L, mj, mL;
952         unsigned int m = 1 << M, n = 1 << N;
953         unsigned int m2 = 1 << (M-1), n2 = 1 << (N-1);
954         unsigned int mn2 = m << (N-1);
955
956         d1[0] *= d2[0];
957         d1[mn2] *= d2[mn2];
958         d1[m2] *= d2[m2];
959         d1[m2 + mn2] *= d2[m2 + mn2];
960         for (i=1; i<m2; i++) {
961                 k = m - i;
962                 a = d1[i]*d2[i] - d1[k]*d2[k];
963                 b = d1[k]*d2[i] + d1[i]*d2[k];
964                 d1[i] = (b + a)*(fREAL)0.5;
965                 d1[k] = (b - a)*(fREAL)0.5;
966                 a = d1[i + mn2]*d2[i + mn2] - d1[k + mn2]*d2[k + mn2];
967                 b = d1[k + mn2]*d2[i + mn2] + d1[i + mn2]*d2[k + mn2];
968                 d1[i + mn2] = (b + a)*(fREAL)0.5;
969                 d1[k + mn2] = (b - a)*(fREAL)0.5;
970         }
971         for (j=1; j<n2; j++) {
972                 L = n - j;
973                 mj = j << M;
974                 mL = L << M;
975                 a = d1[mj]*d2[mj] - d1[mL]*d2[mL];
976                 b = d1[mL]*d2[mj] + d1[mj]*d2[mL];
977                 d1[mj] = (b + a)*(fREAL)0.5;
978                 d1[mL] = (b - a)*(fREAL)0.5;
979                 a = d1[m2 + mj]*d2[m2 + mj] - d1[m2 + mL]*d2[m2 + mL];
980                 b = d1[m2 + mL]*d2[m2 + mj] + d1[m2 + mj]*d2[m2 + mL];
981                 d1[m2 + mj] = (b + a)*(fREAL)0.5;
982                 d1[m2 + mL] = (b - a)*(fREAL)0.5;
983         }
984         for (i=1; i<m2; i++) {
985                 k = m - i;
986                 for (j=1; j<n2; j++) {
987                         L = n - j;
988                         mj = j << M;
989                         mL = L << M;
990                         a = d1[i + mj]*d2[i + mj] - d1[k + mL]*d2[k + mL];
991                         b = d1[k + mL]*d2[i + mj] + d1[i + mj]*d2[k + mL];
992                         d1[i + mj] = (b + a)*(fREAL)0.5;
993                         d1[k + mL] = (b - a)*(fREAL)0.5;
994                         a = d1[i + mL]*d2[i + mL] - d1[k + mj]*d2[k + mj];
995                         b = d1[k + mj]*d2[i + mL] + d1[i + mL]*d2[k + mj];
996                         d1[i + mL] = (b + a)*(fREAL)0.5;
997                         d1[k + mj] = (b - a)*(fREAL)0.5;
998                 }
999         }
1000 }
1001
1002 //------------------------------------------------------------------------------
1003
1004 void convolve(CompBuf* dst, CompBuf* in1, CompBuf* in2)
1005 {
1006         fREAL *data1, *data2, *fp;
1007         unsigned int w2, h2, hw, hh, log2_w, log2_h;
1008         fRGB wt, *colp;
1009         int x, y, ch;
1010         int xbl, ybl, nxb, nyb, xbsz, ybsz;
1011         int in2done = FALSE;
1012
1013         CompBuf* rdst = alloc_compbuf(in1->x, in1->y, in1->type, 1);
1014
1015         // convolution result width & height
1016         w2 = 2*in2->x - 1;
1017         h2 = 2*in2->y - 1;
1018         // FFT pow2 required size & log2
1019         w2 = nextPow2(w2, &log2_w);
1020         h2 = nextPow2(h2, &log2_h);
1021
1022         // alloc space
1023         data1 = (fREAL*)MEM_callocN(3*w2*h2*sizeof(fREAL), "convolve_fast FHT data1");
1024         data2 = (fREAL*)MEM_callocN(w2*h2*sizeof(fREAL), "convolve_fast FHT data2");
1025
1026         // normalize convolutor
1027         wt[0] = wt[1] = wt[2] = 0.f;
1028         for (y=0; y<in2->y; y++) {
1029                 colp = (fRGB*)&in2->rect[y*in2->x*in2->type];
1030                 for (x=0; x<in2->x; x++)
1031                         add_v3_v3(wt, colp[x]);
1032         }
1033         if (wt[0] != 0.f) wt[0] = 1.f/wt[0];
1034         if (wt[1] != 0.f) wt[1] = 1.f/wt[1];
1035         if (wt[2] != 0.f) wt[2] = 1.f/wt[2];
1036         for (y=0; y<in2->y; y++) {
1037                 colp = (fRGB*)&in2->rect[y*in2->x*in2->type];
1038                 for (x=0; x<in2->x; x++)
1039                         mul_v3_v3(colp[x], wt);
1040         }
1041
1042         // copy image data, unpacking interleaved RGBA into separate channels
1043         // only need to calc data1 once
1044
1045         // block add-overlap
1046         hw = in2->x >> 1;
1047         hh = in2->y >> 1;
1048         xbsz = (w2 + 1) - in2->x;
1049         ybsz = (h2 + 1) - in2->y;
1050         nxb = in1->x / xbsz;
1051         if (in1->x % xbsz) nxb++;
1052         nyb = in1->y / ybsz;
1053         if (in1->y % ybsz) nyb++;
1054         for (ybl=0; ybl<nyb; ybl++) {
1055                 for (xbl=0; xbl<nxb; xbl++) {
1056
1057                         // each channel one by one
1058                         for (ch=0; ch<3; ch++) {
1059                                 fREAL* data1ch = &data1[ch*w2*h2];
1060
1061                                 // only need to calc fht data from in2 once, can re-use for every block
1062                                 if (!in2done) {
1063                                         // in2, channel ch -> data1
1064                                         for (y=0; y<in2->y; y++) {
1065                                                 fp = &data1ch[y*w2];
1066                                                 colp = (fRGB*)&in2->rect[y*in2->x*in2->type];
1067                                                 for (x=0; x<in2->x; x++)
1068                                                         fp[x] = colp[x][ch];
1069                                         }
1070                                 }
1071
1072                                 // in1, channel ch -> data2
1073                                 memset(data2, 0, w2*h2*sizeof(fREAL));
1074                                 for (y=0; y<ybsz; y++) {
1075                                         int yy = ybl*ybsz + y;
1076                                         if (yy >= in1->y) continue;
1077                                         fp = &data2[y*w2];
1078                                         colp = (fRGB*)&in1->rect[yy*in1->x*in1->type];
1079                                         for (x=0; x<xbsz; x++) {
1080                                                 int xx = xbl*xbsz + x;
1081                                                 if (xx >= in1->x) continue;
1082                                                 fp[x] = colp[xx][ch];
1083                                         }
1084                                 }
1085
1086                                 // forward FHT
1087                                 // zero pad data start is different for each == height+1
1088                                 if (!in2done) FHT2D(data1ch, log2_w, log2_h, in2->y+1, 0);
1089                                 FHT2D(data2, log2_w, log2_h, in2->y+1, 0);
1090
1091                                 // FHT2D transposed data, row/col now swapped
1092                                 // convolve & inverse FHT
1093                                 fht_convolve(data2, data1ch, log2_h, log2_w);
1094                                 FHT2D(data2, log2_h, log2_w, 0, 1);
1095                                 // data again transposed, so in order again
1096
1097                                 // overlap-add result
1098                                 for (y=0; y<(int)h2; y++) {
1099                                         const int yy = ybl*ybsz + y - hh;
1100                                         if ((yy < 0) || (yy >= in1->y)) continue;
1101                                         fp = &data2[y*w2];
1102                                         colp = (fRGB*)&rdst->rect[yy*in1->x*in1->type];
1103                                         for (x=0; x<(int)w2; x++) {
1104                                                 const int xx = xbl*xbsz + x - hw;
1105                                                 if ((xx < 0) || (xx >= in1->x)) continue;
1106                                                 colp[xx][ch] += fp[x];
1107                                         }
1108                                 }
1109
1110                         }
1111                         in2done = TRUE;
1112                 }
1113         }
1114
1115         MEM_freeN(data2);
1116         MEM_freeN(data1);
1117         memcpy(dst->rect, rdst->rect, sizeof(float)*dst->x*dst->y*dst->type);
1118         free_compbuf(rdst);
1119 }
1120
1121
1122 /*
1123  *
1124  * Utility functions qd_* should probably be intergrated better with other functions here.
1125  *
1126  */
1127 // sets fcol to pixelcolor at (x, y)
1128 void qd_getPixel(CompBuf* src, int x, int y, float* col)
1129 {
1130         if (src->rect_procedural) {
1131                 float bc[4];
1132                 src->rect_procedural(src, bc, (float)x/(float)src->xrad, (float)y/(float)src->yrad);
1133
1134                 switch (src->type) {
1135                         /* these fallthrough to get all the channels */
1136                         case CB_RGBA: col[3]=bc[3]; 
1137                         case CB_VEC3: col[2]=bc[2];
1138                         case CB_VEC2: col[1]=bc[1];
1139                         case CB_VAL: col[0]=bc[0];
1140                 }
1141         }
1142         else if ((x >= 0) && (x < src->x) && (y >= 0) && (y < src->y)) {
1143                 float* bc = &src->rect[(x + y*src->x)*src->type];
1144                 switch (src->type) {
1145                         /* these fallthrough to get all the channels */
1146                         case CB_RGBA: col[3]=bc[3]; 
1147                         case CB_VEC3: col[2]=bc[2];
1148                         case CB_VEC2: col[1]=bc[1];
1149                         case CB_VAL: col[0]=bc[0];
1150                 }
1151         }
1152         else {
1153                 switch (src->type) {
1154                         /* these fallthrough to get all the channels */
1155                         case CB_RGBA: col[3]=0.0; 
1156                         case CB_VEC3: col[2]=0.0; 
1157                         case CB_VEC2: col[1]=0.0; 
1158                         case CB_VAL: col[0]=0.0; 
1159                 }
1160         }
1161 }
1162
1163 // sets pixel (x, y) to color col
1164 void qd_setPixel(CompBuf* src, int x, int y, float* col)
1165 {
1166         if ((x >= 0) && (x < src->x) && (y >= 0) && (y < src->y)) {
1167                 float* bc = &src->rect[(x + y*src->x)*src->type];
1168                 switch (src->type) {
1169                         /* these fallthrough to get all the channels */
1170                         case CB_RGBA: bc[3]=col[3]; 
1171                         case CB_VEC3: bc[2]=col[2];
1172                         case CB_VEC2: bc[1]=col[1];
1173                         case CB_VAL: bc[0]=col[0];
1174                 }
1175         }
1176 }
1177
1178 // adds fcol to pixelcolor (x, y)
1179 void qd_addPixel(CompBuf* src, int x, int y, float* col)
1180 {
1181         if ((x >= 0) && (x < src->x) && (y >= 0) && (y < src->y)) {
1182                 float* bc = &src->rect[(x + y*src->x)*src->type];
1183                 bc[0] += col[0], bc[1] += col[1], bc[2] += col[2];
1184         }
1185 }
1186
1187 // multiplies pixel by factor value f
1188 void qd_multPixel(CompBuf* src, int x, int y, float f)
1189 {
1190         if ((x >= 0) && (x < src->x) && (y >= 0) && (y < src->y)) {
1191                 float* bc = &src->rect[(x + y*src->x)*src->type];
1192                 bc[0] *= f, bc[1] *= f, bc[2] *= f;
1193         }
1194 }
1195
1196 // bilinear interpolation with wraparound
1197 void qd_getPixelLerpWrap(CompBuf* src, float u, float v, float* col)
1198 {
1199         const float ufl = floor(u), vfl = floor(v);
1200         const int nx = (int)ufl % src->x, ny = (int)vfl % src->y;
1201         const int x1 = (nx < 0) ? (nx + src->x) : nx;
1202         const int y1 = (ny < 0) ? (ny + src->y) : ny;
1203         const int x2 = (x1 + 1) % src->x, y2 = (y1 + 1) % src->y;
1204         const float* c00 = &src->rect[(x1 + y1*src->x)*src->type];
1205         const float* c10 = &src->rect[(x2 + y1*src->x)*src->type];
1206         const float* c01 = &src->rect[(x1 + y2*src->x)*src->type];
1207         const float* c11 = &src->rect[(x2 + y2*src->x)*src->type];
1208         const float uf = u - ufl, vf = v - vfl;
1209         const float w00=(1.f-uf)*(1.f-vf), w10=uf*(1.f-vf), w01=(1.f-uf)*vf, w11=uf*vf;
1210         col[0] = w00*c00[0] + w10*c10[0] + w01*c01[0] + w11*c11[0];
1211         if (src->type != CB_VAL) {
1212                 col[1] = w00*c00[1] + w10*c10[1] + w01*c01[1] + w11*c11[1];
1213                 col[2] = w00*c00[2] + w10*c10[2] + w01*c01[2] + w11*c11[2];
1214                 col[3] = w00*c00[3] + w10*c10[3] + w01*c01[3] + w11*c11[3];
1215         }
1216 }
1217
1218 // as above, without wrap around
1219 void qd_getPixelLerp(CompBuf* src, float u, float v, float* col)
1220 {
1221         const float ufl = floor(u), vfl = floor(v);
1222         const int x1 = (int)ufl, y1 = (int)vfl;
1223         const int x2 = (int)ceil(u), y2 = (int)ceil(v);
1224         if ((x2 >= 0) && (y2 >= 0) && (x1 < src->x) && (y1 < src->y)) {
1225                 const float B[4] = {0, 0, 0, 0};
1226                 const int ox1 = (x1 < 0), oy1 = (y1 < 0), ox2 = (x2 >= src->x), oy2 = (y2 >= src->y);
1227                 const float* c00 = (ox1 || oy1) ? B : &src->rect[(x1 + y1*src->x)*src->type];
1228                 const float* c10 = (ox2 || oy1) ? B : &src->rect[(x2 + y1*src->x)*src->type];
1229                 const float* c01 = (ox1 || oy2) ? B : &src->rect[(x1 + y2*src->x)*src->type];
1230                 const float* c11 = (ox2 || oy2) ? B : &src->rect[(x2 + y2*src->x)*src->type];
1231                 const float uf = u - ufl, vf = v - vfl;
1232                 const float w00=(1.f-uf)*(1.f-vf), w10=uf*(1.f-vf), w01=(1.f-uf)*vf, w11=uf*vf;
1233                 col[0] = w00*c00[0] + w10*c10[0] + w01*c01[0] + w11*c11[0];
1234                 if (src->type != CB_VAL) {
1235                         col[1] = w00*c00[1] + w10*c10[1] + w01*c01[1] + w11*c11[1];
1236                         col[2] = w00*c00[2] + w10*c10[2] + w01*c01[2] + w11*c11[2];
1237                         col[3] = w00*c00[3] + w10*c10[3] + w01*c01[3] + w11*c11[3];
1238                 }
1239         }
1240         else col[0] = col[1] = col[2] = col[3] = 0.f;
1241 }
1242
1243 // as above, sampling only one channel
1244 void qd_getPixelLerpChan(CompBuf* src, float u, float v, int chan, float* out)
1245 {
1246         const float ufl = floor(u), vfl = floor(v);
1247         const int x1 = (int)ufl, y1 = (int)vfl;
1248         const int x2 = (int)ceil(u), y2 = (int)ceil(v);
1249         if (chan >= src->type) chan = 0;
1250         if ((x2 >= 0) && (y2 >= 0) && (x1 < src->x) && (y1 < src->y)) {
1251                 const float B[4] = {0, 0, 0, 0};
1252                 const int ox1 = (x1 < 0), oy1 = (y1 < 0), ox2 = (x2 >= src->x), oy2 = (y2 >= src->y);
1253                 const float* c00 = (ox1 || oy1) ? B : &src->rect[(x1 + y1*src->x)*src->type + chan];
1254                 const float* c10 = (ox2 || oy1) ? B : &src->rect[(x2 + y1*src->x)*src->type + chan];
1255                 const float* c01 = (ox1 || oy2) ? B : &src->rect[(x1 + y2*src->x)*src->type + chan];
1256                 const float* c11 = (ox2 || oy2) ? B : &src->rect[(x2 + y2*src->x)*src->type + chan];
1257                 const float uf = u - ufl, vf = v - vfl;
1258                 const float w00=(1.f-uf)*(1.f-vf), w10=uf*(1.f-vf), w01=(1.f-uf)*vf, w11=uf*vf;
1259                 out[0] = w00*c00[0] + w10*c10[0] + w01*c01[0] + w11*c11[0];
1260         }
1261         else *out = 0.f;
1262 }
1263
1264
1265 CompBuf* qd_downScaledCopy(CompBuf* src, int scale)
1266 {
1267         CompBuf* fbuf;
1268         if (scale <= 1)
1269                 fbuf = dupalloc_compbuf(src);
1270         else {
1271                 int nw = src->x/scale, nh = src->y/scale;
1272                 if ((2*(src->x % scale)) > scale) nw++;
1273                 if ((2*(src->y % scale)) > scale) nh++;
1274                 fbuf = alloc_compbuf(nw, nh, src->type, 1);
1275                 {
1276                         int x, y, xx, yy, sx, sy, mx, my;
1277                         float colsum[4] = {0.0f, 0.0f, 0.0f, 0.0f};
1278                         float fscale = 1.f/(float)(scale*scale);
1279                         for (y=0; y<nh; y++) {
1280                                 fRGB* fcolp = (fRGB*)&fbuf->rect[y*fbuf->x*fbuf->type];
1281                                 yy = y*scale;
1282                                 my = yy + scale;
1283                                 if (my > src->y) my = src->y;
1284                                 for (x=0; x<nw; x++) {
1285                                         xx = x*scale;
1286                                         mx = xx + scale;
1287                                         if (mx > src->x) mx = src->x;
1288                                         zero_v3(colsum);
1289                                         for (sy=yy; sy<my; sy++) {
1290                                                 fRGB* scolp = (fRGB*)&src->rect[sy*src->x*src->type];
1291                                                 for (sx=xx; sx<mx; sx++)
1292                                                         add_v3_v3(colsum, scolp[sx]);
1293                                         }
1294                                         mul_v3_fl(colsum, fscale);
1295                                         copy_v3_v3(fcolp[x], colsum);
1296                                 }
1297                         }
1298                 }
1299         }
1300         return fbuf;
1301 }
1302
1303 // fast g.blur, per channel
1304 // xy var. bits 1 & 2 ca be used to blur in x or y direction separately
1305 void IIR_gauss(CompBuf* src, float sigma, int chan, int xy)
1306 {
1307         double q, q2, sc, cf[4], tsM[9], tsu[3], tsv[3];
1308         double *X, *Y, *W;
1309         const unsigned int src_width = src->x;
1310         const unsigned int src_height = src->y;
1311         unsigned int i, x, y, sz;
1312
1313         // <0.5 not valid, though can have a possibly useful sort of sharpening effect
1314         if (sigma < 0.5f) return;
1315
1316         if ((xy < 1) || (xy > 3)) xy = 3;
1317
1318         // XXX The YVV macro defined below explicitly expects sources of at least 3x3 pixels,
1319         //     so just skiping blur along faulty direction if src's def is below that limit!
1320         if (src_width < 3) xy &= ~(int) 1;
1321         if (src_height < 3) xy &= ~(int) 2;
1322         if (xy < 1) return;
1323
1324         // see "Recursive Gabor Filtering" by Young/VanVliet
1325         // all factors here in double.prec. Required, because for single.prec it seems to blow up if sigma > ~200
1326         if (sigma >= 3.556f)
1327                 q = 0.9804f * (sigma - 3.556f) + 2.5091f;
1328         else     // sigma >= 0.5
1329                 q = (0.0561f * sigma + 0.5784f) * sigma - 0.2568f;
1330         q2 = q * q;
1331         sc = (1.1668 + q) * (3.203729649  + (2.21566 + q) * q);
1332         // no gabor filtering here, so no complex multiplies, just the regular coefs.
1333         // all negated here, so as not to have to recalc Triggs/Sdika matrix
1334         cf[1] = q * (5.788961737 + (6.76492 + 3.0 * q) * q) / sc;
1335         cf[2] = -q2 * (3.38246 + 3.0 * q) / sc;
1336         // 0 & 3 unchanged
1337         cf[3] = q2 * q / sc;
1338         cf[0] = 1.0 - cf[1] - cf[2] - cf[3];
1339
1340         // Triggs/Sdika border corrections,
1341         // it seems to work, not entirely sure if it is actually totally correct,
1342         // Besides J.M.Geusebroek's anigauss.c (see http://www.science.uva.nl/~mark),
1343         // found one other implementation by Cristoph Lampert,
1344         // but neither seem to be quite the same, result seems to be ok so far anyway.
1345         // Extra scale factor here to not have to do it in filter,
1346         // though maybe this had something to with the precision errors
1347         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]));
1348         tsM[0] = sc * (-cf[3] * cf[1] + 1.0 - cf[3] * cf[3] - cf[2]);
1349         tsM[1] = sc * ((cf[3] + cf[1]) * (cf[2] + cf[3] * cf[1]));
1350         tsM[2] = sc * (cf[3] * (cf[1] + cf[3] * cf[2]));
1351         tsM[3] = sc * (cf[1] + cf[3] * cf[2]);
1352         tsM[4] = sc * (-(cf[2] - 1.0) * (cf[2] + cf[3] * cf[1]));
1353         tsM[5] = sc * (-(cf[3] * cf[1] + cf[3] * cf[3] + cf[2] - 1.0) * cf[3]);
1354         tsM[6] = sc * (cf[3] * cf[1] + cf[2] + cf[1] * cf[1] - cf[2] * cf[2]);
1355         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]);
1356         tsM[8] = sc * (cf[3] * (cf[1] + cf[3] * cf[2]));
1357
1358 #define YVV(L)                                                                          \
1359 {                                                                                       \
1360         W[0] = cf[0] * X[0] + cf[1] * X[0] + cf[2] * X[0] + cf[3] * X[0];                   \
1361         W[1] = cf[0] * X[1] + cf[1] * W[0] + cf[2] * X[0] + cf[3] * X[0];                   \
1362         W[2] = cf[0] * X[2] + cf[1] * W[1] + cf[2] * W[0] + cf[3] * X[0];                   \
1363         for (i = 3; i < L; i++) {                                                           \
1364                 W[i] = cf[0] * X[i] + cf[1] * W[i - 1] + cf[2] * W[i - 2] + cf[3] * W[i - 3];   \
1365         }                                                                                   \
1366         tsu[0] = W[L - 1] - X[L - 1];                                                       \
1367         tsu[1] = W[L - 2] - X[L - 1];                                                       \
1368         tsu[2] = W[L - 3] - X[L - 1];                                                       \
1369         tsv[0] = tsM[0] * tsu[0] + tsM[1] * tsu[1] + tsM[2] * tsu[2] + X[L - 1];            \
1370         tsv[1] = tsM[3] * tsu[0] + tsM[4] * tsu[1] + tsM[5] * tsu[2] + X[L - 1];            \
1371         tsv[2] = tsM[6] * tsu[0] + tsM[7] * tsu[1] + tsM[8] * tsu[2] + X[L - 1];            \
1372         Y[L - 1] = cf[0] * W[L - 1] + cf[1] * tsv[0] + cf[2] * tsv[1] + cf[3] * tsv[2];     \
1373         Y[L - 2] = cf[0] * W[L - 2] + cf[1] * Y[L - 1] + cf[2] * tsv[0] + cf[3] * tsv[1];   \
1374         Y[L - 3] = cf[0] * W[L - 3] + cf[1] * Y[L - 2] + cf[2] * Y[L - 1] + cf[3] * tsv[0]; \
1375         /* 'i != UINT_MAX' is really 'i >= 0', but necessary for unsigned int wrapping */   \
1376         for (i = L - 4; i != UINT_MAX; i--) {                                               \
1377                 Y[i] = cf[0] * W[i] + cf[1] * Y[i + 1] + cf[2] * Y[i + 2] + cf[3] * Y[i + 3];   \
1378         }                                                                                   \
1379 } (void)0
1380
1381         // intermediate buffers
1382         sz = MAX2(src_width, src_height);
1383         X = MEM_callocN(sz * sizeof(double), "IIR_gauss X buf");
1384         Y = MEM_callocN(sz * sizeof(double), "IIR_gauss Y buf");
1385         W = MEM_callocN(sz * sizeof(double), "IIR_gauss W buf");
1386         if (xy & 1) {       // H
1387                 for (y = 0; y < src_height; ++y) {
1388                         const int yx = y * src_width;
1389                         for (x = 0; x < src_width; ++x)
1390                                 X[x] = src->rect[(x + yx) * src->type + chan];
1391                         YVV(src_width);
1392                         for (x = 0; x < src_width; ++x)
1393                                 src->rect[(x + yx) * src->type + chan] = Y[x];
1394                 }
1395         }
1396         if (xy & 2) {       // V
1397                 for (x = 0; x < src_width; ++x) {
1398                         for (y = 0; y < src_height; ++y)
1399                                 X[y] = src->rect[(x + y * src_width) * src->type + chan];
1400                         YVV(src_height);
1401                         for (y = 0; y < src_height; ++y)
1402                                 src->rect[(x + y * src_width) * src->type + chan] = Y[y];
1403                 }
1404         }
1405
1406         MEM_freeN(X);
1407         MEM_freeN(W);
1408         MEM_freeN(Y);
1409 #undef YVV
1410 }
1411
1412 #endif  /* WITH_COMPOSITOR_LEGACY */