style cleanup: blenkernel
[blender.git] / source / blender / blenkernel / intern / colortools.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) 2005 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/BL DUAL LICENSE BLOCK *****
26  */
27
28 /** \file blender/blenkernel/intern/colortools.c
29  *  \ingroup bke
30  */
31
32
33 #include <string.h>
34 #include <math.h>
35 #include <stdlib.h>
36 #include <float.h>
37
38 #include "MEM_guardedalloc.h"
39
40 #include "DNA_color_types.h"
41 #include "DNA_curve_types.h"
42
43 #include "BLI_blenlib.h"
44 #include "BLI_math.h"
45 #include "BLI_utildefines.h"
46
47 #include "BKE_colortools.h"
48 #include "BKE_curve.h"
49 #include "BKE_fcurve.h"
50
51
52 #include "IMB_imbuf.h"
53 #include "IMB_imbuf_types.h"
54
55 /* ********************************* color curve ********************* */
56
57 /* ***************** operations on full struct ************* */
58
59 CurveMapping *curvemapping_add(int tot, float minx, float miny, float maxx, float maxy)
60 {
61         CurveMapping *cumap;
62         int a;
63         float clipminx, clipminy, clipmaxx, clipmaxy;
64         
65         cumap = MEM_callocN(sizeof(CurveMapping), "new curvemap");
66         cumap->flag = CUMA_DO_CLIP;
67         if (tot == 4) cumap->cur = 3;   /* rhms, hack for 'col' curve? */
68         
69         clipminx = MIN2(minx, maxx);
70         clipminy = MIN2(miny, maxy);
71         clipmaxx = MAX2(minx, maxx);
72         clipmaxy = MAX2(miny, maxy);
73         
74         BLI_init_rctf(&cumap->curr, clipminx, clipmaxx, clipminy, clipmaxy);
75         cumap->clipr = cumap->curr;
76         
77         cumap->white[0] = cumap->white[1] = cumap->white[2] = 1.0f;
78         cumap->bwmul[0] = cumap->bwmul[1] = cumap->bwmul[2] = 1.0f;
79         
80         for (a = 0; a < tot; a++) {
81                 cumap->cm[a].flag = CUMA_EXTEND_EXTRAPOLATE;
82                 cumap->cm[a].totpoint = 2;
83                 cumap->cm[a].curve = MEM_callocN(2 * sizeof(CurveMapPoint), "curve points");
84
85                 cumap->cm[a].curve[0].x = minx;
86                 cumap->cm[a].curve[0].y = miny;
87                 cumap->cm[a].curve[1].x = maxx;
88                 cumap->cm[a].curve[1].y = maxy;
89         }       
90
91         cumap->changed_timestamp = 0;
92
93         return cumap;
94 }
95
96 void curvemapping_free(CurveMapping *cumap)
97 {
98         int a;
99         
100         if (cumap) {
101                 for (a = 0; a < CM_TOT; a++) {
102                         if (cumap->cm[a].curve) MEM_freeN(cumap->cm[a].curve);
103                         if (cumap->cm[a].table) MEM_freeN(cumap->cm[a].table);
104                         if (cumap->cm[a].premultable) MEM_freeN(cumap->cm[a].premultable);
105                 }
106                 MEM_freeN(cumap);
107         }
108 }
109
110 CurveMapping *curvemapping_copy(CurveMapping *cumap)
111 {
112         int a;
113         
114         if (cumap) {
115                 CurveMapping *cumapn = MEM_dupallocN(cumap);
116                 for (a = 0; a < CM_TOT; a++) {
117                         if (cumap->cm[a].curve) 
118                                 cumapn->cm[a].curve = MEM_dupallocN(cumap->cm[a].curve);
119                         if (cumap->cm[a].table) 
120                                 cumapn->cm[a].table = MEM_dupallocN(cumap->cm[a].table);
121                         if (cumap->cm[a].premultable) 
122                                 cumapn->cm[a].premultable = MEM_dupallocN(cumap->cm[a].premultable);
123                 }
124                 return cumapn;
125         }
126         return NULL;
127 }
128
129 void curvemapping_set_black_white(CurveMapping *cumap, const float black[3], const float white[3])
130 {
131         int a;
132         
133         if (white)
134                 copy_v3_v3(cumap->white, white);
135         if (black)
136                 copy_v3_v3(cumap->black, black);
137         
138         for (a = 0; a < 3; a++) {
139                 if (cumap->white[a] == cumap->black[a])
140                         cumap->bwmul[a] = 0.0f;
141                 else
142                         cumap->bwmul[a] = 1.0f / (cumap->white[a] - cumap->black[a]);
143         }       
144 }
145
146 /* ***************** operations on single curve ************* */
147 /* ********** NOTE: requires curvemapping_changed() call after ******** */
148
149 /* removes with flag set */
150 void curvemap_remove(CurveMap *cuma, int flag)
151 {
152         CurveMapPoint *cmp = MEM_mallocN((cuma->totpoint) * sizeof(CurveMapPoint), "curve points");
153         int a, b, removed = 0;
154         
155         /* well, lets keep the two outer points! */
156         cmp[0] = cuma->curve[0];
157         for (a = 1, b = 1; a < cuma->totpoint - 1; a++) {
158                 if (!(cuma->curve[a].flag & flag)) {
159                         cmp[b] = cuma->curve[a];
160                         b++;
161                 }
162                 else removed++;
163         }
164         cmp[b] = cuma->curve[a];
165         
166         MEM_freeN(cuma->curve);
167         cuma->curve = cmp;
168         cuma->totpoint -= removed;
169 }
170
171 void curvemap_insert(CurveMap *cuma, float x, float y)
172 {
173         CurveMapPoint *cmp = MEM_callocN((cuma->totpoint + 1) * sizeof(CurveMapPoint), "curve points");
174         int a, b, foundloc = 0;
175                 
176         /* insert fragments of the old one and the new point to the new curve */
177         cuma->totpoint++;
178         for (a = 0, b = 0; a < cuma->totpoint; a++) {
179                 if ((x < cuma->curve[a].x) && !foundloc) {
180                         cmp[a].x = x;
181                         cmp[a].y = y;
182                         cmp[a].flag = CUMA_SELECT;
183                         foundloc = 1;
184                 }
185                 else {
186                         cmp[a].x = cuma->curve[b].x;
187                         cmp[a].y = cuma->curve[b].y;
188                         cmp[a].flag = cuma->curve[b].flag;
189                         cmp[a].flag &= ~CUMA_SELECT; /* make sure old points don't remain selected */
190                         cmp[a].shorty = cuma->curve[b].shorty;
191                         b++;
192                 }
193         }
194
195         /* free old curve and replace it with new one */
196         MEM_freeN(cuma->curve);
197         cuma->curve = cmp;
198 }
199
200 void curvemap_reset(CurveMap *cuma, rctf *clipr, int preset, int slope)
201 {
202         if (cuma->curve)
203                 MEM_freeN(cuma->curve);
204
205         switch (preset) {
206                 case CURVE_PRESET_LINE: cuma->totpoint = 2; break;
207                 case CURVE_PRESET_SHARP: cuma->totpoint = 4; break;
208                 case CURVE_PRESET_SMOOTH: cuma->totpoint = 4; break;
209                 case CURVE_PRESET_MAX: cuma->totpoint = 2; break;
210                 case CURVE_PRESET_MID9: cuma->totpoint = 9; break;
211                 case CURVE_PRESET_ROUND: cuma->totpoint = 4; break;
212                 case CURVE_PRESET_ROOT: cuma->totpoint = 4; break;
213         }
214
215         cuma->curve = MEM_callocN(cuma->totpoint * sizeof(CurveMapPoint), "curve points");
216
217         switch (preset) {
218                 case CURVE_PRESET_LINE:
219                         cuma->curve[0].x = clipr->xmin;
220                         cuma->curve[0].y = clipr->ymax;
221                         cuma->curve[0].flag = 0;
222                         cuma->curve[1].x = clipr->xmax;
223                         cuma->curve[1].y = clipr->ymin;
224                         cuma->curve[1].flag = 0;
225                         break;
226                 case CURVE_PRESET_SHARP:
227                         cuma->curve[0].x = 0;
228                         cuma->curve[0].y = 1;
229                         cuma->curve[1].x = 0.25;
230                         cuma->curve[1].y = 0.50;
231                         cuma->curve[2].x = 0.75;
232                         cuma->curve[2].y = 0.04;
233                         cuma->curve[3].x = 1;
234                         cuma->curve[3].y = 0;
235                         break;
236                 case CURVE_PRESET_SMOOTH:
237                         cuma->curve[0].x = 0;
238                         cuma->curve[0].y = 1;
239                         cuma->curve[1].x = 0.25;
240                         cuma->curve[1].y = 0.94;
241                         cuma->curve[2].x = 0.75;
242                         cuma->curve[2].y = 0.06;
243                         cuma->curve[3].x = 1;
244                         cuma->curve[3].y = 0;
245                         break;
246                 case CURVE_PRESET_MAX:
247                         cuma->curve[0].x = 0;
248                         cuma->curve[0].y = 1;
249                         cuma->curve[1].x = 1;
250                         cuma->curve[1].y = 1;
251                         break;
252                 case CURVE_PRESET_MID9:
253                 {
254                         int i;
255                         for (i = 0; i < cuma->totpoint; i++) {
256                                 cuma->curve[i].x = i / ((float)cuma->totpoint - 1);
257                                 cuma->curve[i].y = 0.5;
258                         }
259                 }
260                 break;
261                 case CURVE_PRESET_ROUND:
262                         cuma->curve[0].x = 0;
263                         cuma->curve[0].y = 1;
264                         cuma->curve[1].x = 0.5;
265                         cuma->curve[1].y = 0.90;
266                         cuma->curve[2].x = 0.86;
267                         cuma->curve[2].y = 0.5;
268                         cuma->curve[3].x = 1;
269                         cuma->curve[3].y = 0;
270                         break;
271                 case CURVE_PRESET_ROOT:
272                         cuma->curve[0].x = 0;
273                         cuma->curve[0].y = 1;
274                         cuma->curve[1].x = 0.25;
275                         cuma->curve[1].y = 0.95;
276                         cuma->curve[2].x = 0.75;
277                         cuma->curve[2].y = 0.44;
278                         cuma->curve[3].x = 1;
279                         cuma->curve[3].y = 0;
280                         break;
281         }
282
283         /* mirror curve in x direction to have positive slope
284          * rather than default negative slope */
285         if (slope == CURVEMAP_SLOPE_POSITIVE) {
286                 int i, last = cuma->totpoint - 1;
287                 CurveMapPoint *newpoints = MEM_dupallocN(cuma->curve);
288                 
289                 for (i = 0; i < cuma->totpoint; i++) {
290                         newpoints[i].y = cuma->curve[last - i].y;
291                 }
292                 
293                 MEM_freeN(cuma->curve);
294                 cuma->curve = newpoints;
295         }
296         
297         if (cuma->table) {
298                 MEM_freeN(cuma->table);
299                 cuma->table = NULL;
300         }
301 }
302
303 /* if type==1: vector, else auto */
304 void curvemap_sethandle(CurveMap *cuma, int type)
305 {
306         int a;
307         
308         for (a = 0; a < cuma->totpoint; a++) {
309                 if (cuma->curve[a].flag & CUMA_SELECT) {
310                         if (type) cuma->curve[a].flag |= CUMA_VECTOR;
311                         else cuma->curve[a].flag &= ~CUMA_VECTOR;
312                 }
313         }
314 }
315
316 /* *********************** Making the tables and display ************** */
317
318 /* reduced copy of garbled calchandleNurb() code in curve.c */
319 static void calchandle_curvemap(BezTriple *bezt, BezTriple *prev, BezTriple *next, int UNUSED(mode))
320 {
321         float *p1, *p2, *p3, pt[3];
322         float len, len_a, len_b;
323         float dvec_a[2], dvec_b[2];
324
325         if (bezt->h1 == 0 && bezt->h2 == 0) {
326                 return;
327         }
328         
329         p2 = bezt->vec[1];
330         
331         if (prev == NULL) {
332                 p3 = next->vec[1];
333                 pt[0] = 2.0f * p2[0] - p3[0];
334                 pt[1] = 2.0f * p2[1] - p3[1];
335                 p1 = pt;
336         }
337         else {
338                 p1 = prev->vec[1];
339         }
340         
341         if (next == NULL) {
342                 p1 = prev->vec[1];
343                 pt[0] = 2.0f * p2[0] - p1[0];
344                 pt[1] = 2.0f * p2[1] - p1[1];
345                 p3 = pt;
346         }
347         else {
348                 p3 = next->vec[1];
349         }
350
351         sub_v2_v2v2(dvec_a, p2, p1);
352         sub_v2_v2v2(dvec_b, p3, p2);
353
354         len_a = len_v2(dvec_a);
355         len_b = len_v2(dvec_b);
356
357         if (len_a == 0.0f) len_a = 1.0f;
358         if (len_b == 0.0f) len_b = 1.0f;
359
360         if (bezt->h1 == HD_AUTO || bezt->h2 == HD_AUTO) { /* auto */
361                 float tvec[2];
362                 tvec[0] = dvec_b[0] / len_b + dvec_a[0] / len_a;
363                 tvec[1] = dvec_b[1] / len_b + dvec_a[1] / len_a;
364
365                 len = len_v2(tvec) * 2.5614f;
366                 if (len != 0.0f) {
367                         
368                         if (bezt->h1 == HD_AUTO) {
369                                 len_a /= len;
370                                 madd_v2_v2v2fl(p2 - 3, p2, tvec, -len_a);
371                         }
372                         if (bezt->h2 == HD_AUTO) {
373                                 len_b /= len;
374                                 madd_v2_v2v2fl(p2 + 3, p2, tvec,  len_b);
375                         }
376                 }
377         }
378
379         if (bezt->h1 == HD_VECT) {    /* vector */
380                 madd_v2_v2v2fl(p2 - 3, p2, dvec_a, -1.0f / 3.0f);
381         }
382         if (bezt->h2 == HD_VECT) {
383                 madd_v2_v2v2fl(p2 + 3, p2, dvec_b,  1.0f / 3.0f);
384         }
385 }
386
387 /* in X, out Y. 
388  * X is presumed to be outside first or last */
389 static float curvemap_calc_extend(CurveMap *cuma, float x, const float first[2], const float last[2])
390 {
391         if (x <= first[0]) {
392                 if ((cuma->flag & CUMA_EXTEND_EXTRAPOLATE) == 0) {
393                         /* no extrapolate */
394                         return first[1];
395                 }
396                 else {
397                         if (cuma->ext_in[0] == 0.0f)
398                                 return first[1] + cuma->ext_in[1] * 10000.0f;
399                         else
400                                 return first[1] + cuma->ext_in[1] * (x - first[0]) / cuma->ext_in[0];
401                 }
402         }
403         else if (x >= last[0]) {
404                 if ((cuma->flag & CUMA_EXTEND_EXTRAPOLATE) == 0) {
405                         /* no extrapolate */
406                         return last[1];
407                 }
408                 else {
409                         if (cuma->ext_out[0] == 0.0f)
410                                 return last[1] - cuma->ext_out[1] * 10000.0f;
411                         else
412                                 return last[1] + cuma->ext_out[1] * (x - last[0]) / cuma->ext_out[0];
413                 }
414         }
415         return 0.0f;
416 }
417
418 /* only creates a table for a single channel in CurveMapping */
419 static void curvemap_make_table(CurveMap *cuma, rctf *clipr)
420 {
421         CurveMapPoint *cmp = cuma->curve;
422         BezTriple *bezt;
423         float *fp, *allpoints, *lastpoint, curf, range;
424         int a, totpoint;
425         
426         if (cuma->curve == NULL) return;
427         
428         /* default rect also is table range */
429         cuma->mintable = clipr->xmin;
430         cuma->maxtable = clipr->xmax;
431         
432         /* hrmf... we now rely on blender ipo beziers, these are more advanced */
433         bezt = MEM_callocN(cuma->totpoint * sizeof(BezTriple), "beztarr");
434         
435         for (a = 0; a < cuma->totpoint; a++) {
436                 cuma->mintable = MIN2(cuma->mintable, cmp[a].x);
437                 cuma->maxtable = MAX2(cuma->maxtable, cmp[a].x);
438                 bezt[a].vec[1][0] = cmp[a].x;
439                 bezt[a].vec[1][1] = cmp[a].y;
440                 if (cmp[a].flag & CUMA_VECTOR)
441                         bezt[a].h1 = bezt[a].h2 = HD_VECT;
442                 else
443                         bezt[a].h1 = bezt[a].h2 = HD_AUTO;
444         }
445         
446         for (a = 0; a < cuma->totpoint; a++) {
447                 if (a == 0)
448                         calchandle_curvemap(bezt, NULL, bezt + 1, 0);
449                 else if (a == cuma->totpoint - 1)
450                         calchandle_curvemap(bezt + a, bezt + a - 1, NULL, 0);
451                 else
452                         calchandle_curvemap(bezt + a, bezt + a - 1, bezt + a + 1, 0);
453         }
454         
455         /* first and last handle need correction, instead of pointing to center of next/prev, 
456          * we let it point to the closest handle */
457         if (cuma->totpoint > 2) {
458                 float hlen, nlen, vec[3];
459                 
460                 if (bezt[0].h2 == HD_AUTO) {
461                         
462                         hlen = len_v3v3(bezt[0].vec[1], bezt[0].vec[2]); /* original handle length */
463                         /* clip handle point */
464                         copy_v3_v3(vec, bezt[1].vec[0]);
465                         if (vec[0] < bezt[0].vec[1][0])
466                                 vec[0] = bezt[0].vec[1][0];
467                         
468                         sub_v3_v3(vec, bezt[0].vec[1]);
469                         nlen = len_v3(vec);
470                         if (nlen > FLT_EPSILON) {
471                                 mul_v3_fl(vec, hlen / nlen);
472                                 add_v3_v3v3(bezt[0].vec[2], vec, bezt[0].vec[1]);
473                                 sub_v3_v3v3(bezt[0].vec[0], bezt[0].vec[1], vec);
474                         }
475                 }
476                 a = cuma->totpoint - 1;
477                 if (bezt[a].h2 == HD_AUTO) {
478                         
479                         hlen = len_v3v3(bezt[a].vec[1], bezt[a].vec[0]); /* original handle length */
480                         /* clip handle point */
481                         copy_v3_v3(vec, bezt[a - 1].vec[2]);
482                         if (vec[0] > bezt[a].vec[1][0])
483                                 vec[0] = bezt[a].vec[1][0];
484                         
485                         sub_v3_v3(vec, bezt[a].vec[1]);
486                         nlen = len_v3(vec);
487                         if (nlen > FLT_EPSILON) {
488                                 mul_v3_fl(vec, hlen / nlen);
489                                 add_v3_v3v3(bezt[a].vec[0], vec, bezt[a].vec[1]);
490                                 sub_v3_v3v3(bezt[a].vec[2], bezt[a].vec[1], vec);
491                         }
492                 }
493         }       
494         /* make the bezier curve */
495         if (cuma->table)
496                 MEM_freeN(cuma->table);
497         totpoint = (cuma->totpoint - 1) * CM_RESOL;
498         fp = allpoints = MEM_callocN(totpoint * 2 * sizeof(float), "table");
499         
500         for (a = 0; a < cuma->totpoint - 1; a++, fp += 2 * CM_RESOL) {
501                 correct_bezpart(bezt[a].vec[1], bezt[a].vec[2], bezt[a + 1].vec[0], bezt[a + 1].vec[1]);
502                 BKE_curve_forward_diff_bezier(bezt[a].vec[1][0], bezt[a].vec[2][0], bezt[a + 1].vec[0][0], bezt[a + 1].vec[1][0], fp, CM_RESOL - 1, 2 * sizeof(float));
503                 BKE_curve_forward_diff_bezier(bezt[a].vec[1][1], bezt[a].vec[2][1], bezt[a + 1].vec[0][1], bezt[a + 1].vec[1][1], fp + 1, CM_RESOL - 1, 2 * sizeof(float));
504         }
505         
506         /* store first and last handle for extrapolation, unit length */
507         cuma->ext_in[0] = bezt[0].vec[0][0] - bezt[0].vec[1][0];
508         cuma->ext_in[1] = bezt[0].vec[0][1] - bezt[0].vec[1][1];
509         range = sqrt(cuma->ext_in[0] * cuma->ext_in[0] + cuma->ext_in[1] * cuma->ext_in[1]);
510         cuma->ext_in[0] /= range;
511         cuma->ext_in[1] /= range;
512
513         a = cuma->totpoint - 1;
514         cuma->ext_out[0] = bezt[a].vec[1][0] - bezt[a].vec[2][0];
515         cuma->ext_out[1] = bezt[a].vec[1][1] - bezt[a].vec[2][1];
516         range = sqrt(cuma->ext_out[0] * cuma->ext_out[0] + cuma->ext_out[1] * cuma->ext_out[1]);
517         cuma->ext_out[0] /= range;
518         cuma->ext_out[1] /= range;
519         
520         /* cleanup */
521         MEM_freeN(bezt);
522
523         range = CM_TABLEDIV * (cuma->maxtable - cuma->mintable);
524         cuma->range = 1.0f / range;
525         
526         /* now make a table with CM_TABLE equal x distances */
527         fp = allpoints;
528         lastpoint = allpoints + 2 * (totpoint - 1);
529         cmp = MEM_callocN((CM_TABLE + 1) * sizeof(CurveMapPoint), "dist table");
530
531         for (a = 0; a <= CM_TABLE; a++) {
532                 curf = cuma->mintable + range * (float)a;
533                 cmp[a].x = curf;
534                 
535                 /* get the first x coordinate larger than curf */
536                 while (curf >= fp[0] && fp != lastpoint) {
537                         fp += 2;
538                 }
539                 if (fp == allpoints || (curf >= fp[0] && fp == lastpoint))
540                         cmp[a].y = curvemap_calc_extend(cuma, curf, allpoints, lastpoint);
541                 else {
542                         float fac1 = fp[0] - fp[-2];
543                         float fac2 = fp[0] - curf;
544                         if (fac1 > FLT_EPSILON)
545                                 fac1 = fac2 / fac1;
546                         else
547                                 fac1 = 0.0f;
548                         cmp[a].y = fac1 * fp[-1] + (1.0f - fac1) * fp[1];
549                 }
550         }
551         
552         MEM_freeN(allpoints);
553         cuma->table = cmp;
554 }
555
556 /* call when you do images etc, needs restore too. also verifies tables */
557 /* it uses a flag to prevent premul or free to happen twice */
558 void curvemapping_premultiply(CurveMapping *cumap, int restore)
559 {
560         int a;
561         
562         if (restore) {
563                 if (cumap->flag & CUMA_PREMULLED) {
564                         for (a = 0; a < 3; a++) {
565                                 MEM_freeN(cumap->cm[a].table);
566                                 cumap->cm[a].table = cumap->cm[a].premultable;
567                                 cumap->cm[a].premultable = NULL;
568                         }
569                         
570                         cumap->flag &= ~CUMA_PREMULLED;
571                 }
572         }
573         else {
574                 if ((cumap->flag & CUMA_PREMULLED) == 0) {
575                         /* verify and copy */
576                         for (a = 0; a < 3; a++) {
577                                 if (cumap->cm[a].table == NULL)
578                                         curvemap_make_table(cumap->cm + a, &cumap->clipr);
579                                 cumap->cm[a].premultable = cumap->cm[a].table;
580                                 cumap->cm[a].table = MEM_mallocN((CM_TABLE + 1) * sizeof(CurveMapPoint), "premul table");
581                                 memcpy(cumap->cm[a].table, cumap->cm[a].premultable, (CM_TABLE + 1) * sizeof(CurveMapPoint));
582                         }
583                         
584                         if (cumap->cm[3].table == NULL)
585                                 curvemap_make_table(cumap->cm + 3, &cumap->clipr);
586                 
587                         /* premul */
588                         for (a = 0; a < 3; a++) {
589                                 int b;
590                                 for (b = 0; b <= CM_TABLE; b++) {
591                                         cumap->cm[a].table[b].y = curvemap_evaluateF(cumap->cm + 3, cumap->cm[a].table[b].y);
592                                 }
593                         }
594                         
595                         cumap->flag |= CUMA_PREMULLED;
596                 }
597         }
598 }
599
600 static int sort_curvepoints(const void *a1, const void *a2)
601 {
602         const struct CurveMapPoint *x1 = a1, *x2 = a2;
603         
604         if (x1->x > x2->x) return 1;
605         else if (x1->x < x2->x) return -1;
606         return 0;
607 }
608
609 /* ************************ more CurveMapping calls *************** */
610
611 /* note; only does current curvemap! */
612 void curvemapping_changed(CurveMapping *cumap, int rem_doubles)
613 {
614         CurveMap *cuma = cumap->cm + cumap->cur;
615         CurveMapPoint *cmp = cuma->curve;
616         rctf *clipr = &cumap->clipr;
617         float thresh = 0.01f * (clipr->xmax - clipr->xmin);
618         float dx = 0.0f, dy = 0.0f;
619         int a;
620
621         cumap->changed_timestamp++;
622
623         /* clamp with clip */
624         if (cumap->flag & CUMA_DO_CLIP) {
625                 for (a = 0; a < cuma->totpoint; a++) {
626                         if (cmp[a].flag & CUMA_SELECT) {
627                                 if (cmp[a].x < clipr->xmin)
628                                         dx = MIN2(dx, cmp[a].x - clipr->xmin);
629                                 else if (cmp[a].x > clipr->xmax)
630                                         dx = MAX2(dx, cmp[a].x - clipr->xmax);
631                                 if (cmp[a].y < clipr->ymin)
632                                         dy = MIN2(dy, cmp[a].y - clipr->ymin);
633                                 else if (cmp[a].y > clipr->ymax)
634                                         dy = MAX2(dy, cmp[a].y - clipr->ymax);
635                         }
636                 }
637                 for (a = 0; a < cuma->totpoint; a++) {
638                         if (cmp[a].flag & CUMA_SELECT) {
639                                 cmp[a].x -= dx;
640                                 cmp[a].y -= dy;
641                         }
642                 }
643         }
644         
645         
646         qsort(cmp, cuma->totpoint, sizeof(CurveMapPoint), sort_curvepoints);
647         
648         /* remove doubles, threshold set on 1% of default range */
649         if (rem_doubles && cuma->totpoint > 2) {
650                 for (a = 0; a < cuma->totpoint - 1; a++) {
651                         dx = cmp[a].x - cmp[a + 1].x;
652                         dy = cmp[a].y - cmp[a + 1].y;
653                         if (sqrtf(dx * dx + dy * dy) < thresh) {
654                                 if (a == 0) {
655                                         cmp[a + 1].flag |= 2;
656                                         if (cmp[a + 1].flag & CUMA_SELECT)
657                                                 cmp[a].flag |= CUMA_SELECT;
658                                 }
659                                 else {
660                                         cmp[a].flag |= 2;
661                                         if (cmp[a].flag & CUMA_SELECT)
662                                                 cmp[a + 1].flag |= CUMA_SELECT;
663                                 }
664                                 break;  /* we assume 1 deletion per edit is ok */
665                         }
666                 }
667                 if (a != cuma->totpoint - 1)
668                         curvemap_remove(cuma, 2);
669         }       
670         curvemap_make_table(cuma, clipr);
671 }
672
673 /* table should be verified */
674 float curvemap_evaluateF(CurveMap *cuma, float value)
675 {
676         float fi;
677         int i;
678
679         /* index in table */
680         fi = (value - cuma->mintable) * cuma->range;
681         i = (int)fi;
682         
683         /* fi is table float index and should check against table range i.e. [0.0 CM_TABLE] */
684         if (fi < 0.0f || fi > CM_TABLE)
685                 return curvemap_calc_extend(cuma, value, &cuma->table[0].x, &cuma->table[CM_TABLE].x);
686         else {
687                 if (i < 0) return cuma->table[0].y;
688                 if (i >= CM_TABLE) return cuma->table[CM_TABLE].y;
689                 
690                 fi = fi - (float)i;
691                 return (1.0f - fi) * cuma->table[i].y + (fi) * cuma->table[i + 1].y;
692         }
693 }
694
695 /* works with curve 'cur' */
696 float curvemapping_evaluateF(CurveMapping *cumap, int cur, float value)
697 {
698         CurveMap *cuma = cumap->cm + cur;
699         
700         /* allocate or bail out */
701         if (cuma->table == NULL) {
702                 curvemap_make_table(cuma, &cumap->clipr);
703                 if (cuma->table == NULL)
704                         return 1.0f - value;
705         }
706         return curvemap_evaluateF(cuma, value);
707 }
708
709 /* vector case */
710 void curvemapping_evaluate3F(CurveMapping *cumap, float vecout[3], const float vecin[3])
711 {
712         vecout[0] = curvemapping_evaluateF(cumap, 0, vecin[0]);
713         vecout[1] = curvemapping_evaluateF(cumap, 1, vecin[1]);
714         vecout[2] = curvemapping_evaluateF(cumap, 2, vecin[2]);
715 }
716
717 /* RGB case, no black/white points, no premult */
718 void curvemapping_evaluateRGBF(CurveMapping *cumap, float vecout[3], const float vecin[3])
719 {
720         vecout[0] = curvemapping_evaluateF(cumap, 0, curvemapping_evaluateF(cumap, 3, vecin[0]));
721         vecout[1] = curvemapping_evaluateF(cumap, 1, curvemapping_evaluateF(cumap, 3, vecin[1]));
722         vecout[2] = curvemapping_evaluateF(cumap, 2, curvemapping_evaluateF(cumap, 3, vecin[2]));
723 }
724
725
726 /* RGB with black/white points and premult. tables are checked */
727 void curvemapping_evaluate_premulRGBF(CurveMapping *cumap, float vecout[3], const float vecin[3])
728 {
729         float fac;
730         
731         fac = (vecin[0] - cumap->black[0]) * cumap->bwmul[0];
732         vecout[0] = curvemap_evaluateF(cumap->cm, fac);
733         
734         fac = (vecin[1] - cumap->black[1]) * cumap->bwmul[1];
735         vecout[1] = curvemap_evaluateF(cumap->cm + 1, fac);
736         
737         fac = (vecin[2] - cumap->black[2]) * cumap->bwmul[2];
738         vecout[2] = curvemap_evaluateF(cumap->cm + 2, fac);
739 }
740
741
742 /* only used for image editor curves */
743 void curvemapping_do_ibuf(CurveMapping *cumap, ImBuf *ibuf)
744 {
745         ImBuf *tmpbuf;
746         int pixel;
747         float *pix_in;
748         float col[3];
749         int stride = 4;
750         float *pix_out;
751         
752         if (ibuf == NULL)
753                 return;
754         if (ibuf->rect_float == NULL)
755                 IMB_float_from_rect(ibuf);
756         else if (ibuf->rect == NULL)
757                 imb_addrectImBuf(ibuf);
758         
759         if (!ibuf->rect || !ibuf->rect_float)
760                 return;
761         
762         /* work on a temp buffer, so can color manage afterwards.
763          * No worse off memory wise than comp nodes */
764         tmpbuf = IMB_dupImBuf(ibuf);
765         
766         curvemapping_premultiply(cumap, 0);
767         
768         pix_in = ibuf->rect_float;
769         pix_out = tmpbuf->rect_float;
770
771         if (ibuf->channels)
772                 stride = ibuf->channels;
773         
774         for (pixel = ibuf->x * ibuf->y; pixel > 0; pixel--, pix_in += stride, pix_out += stride) {
775                 if (stride < 3) {
776                         col[0] = curvemap_evaluateF(cumap->cm, *pix_in);
777                         
778                         pix_out[1] = pix_out[2] = pix_out[3] = pix_out[0] = col[0];
779                 }
780                 else {
781                         curvemapping_evaluate_premulRGBF(cumap, col, pix_in);
782                         pix_out[0] = col[0];
783                         pix_out[1] = col[1];
784                         pix_out[2] = col[2];
785                         if (stride > 3)
786                                 pix_out[3] = pix_in[3];
787                         else
788                                 pix_out[3] = 1.f;
789                 }
790         }
791         
792         IMB_rect_from_float(tmpbuf);
793         SWAP(unsigned int *, tmpbuf->rect, ibuf->rect);
794         IMB_freeImBuf(tmpbuf);
795         
796         curvemapping_premultiply(cumap, 1);
797 }
798
799 int curvemapping_RGBA_does_something(CurveMapping *cumap)
800 {
801         int a;
802         
803         if (cumap->black[0] != 0.0f) return 1;
804         if (cumap->black[1] != 0.0f) return 1;
805         if (cumap->black[2] != 0.0f) return 1;
806         if (cumap->white[0] != 1.0f) return 1;
807         if (cumap->white[1] != 1.0f) return 1;
808         if (cumap->white[2] != 1.0f) return 1;
809         
810         for (a = 0; a < CM_TOT; a++) {
811                 if (cumap->cm[a].curve) {
812                         if (cumap->cm[a].totpoint != 2) return 1;
813                         
814                         if (cumap->cm[a].curve[0].x != 0.0f) return 1;
815                         if (cumap->cm[a].curve[0].y != 0.0f) return 1;
816                         if (cumap->cm[a].curve[1].x != 1.0f) return 1;
817                         if (cumap->cm[a].curve[1].y != 1.0f) return 1;
818                 }
819         }
820         return 0;
821 }
822
823 void curvemapping_initialize(CurveMapping *cumap)
824 {
825         int a;
826         
827         if (cumap == NULL) return;
828         
829         for (a = 0; a < CM_TOT; a++) {
830                 if (cumap->cm[a].table == NULL)
831                         curvemap_make_table(cumap->cm + a, &cumap->clipr);
832         }
833 }
834
835 void curvemapping_table_RGBA(CurveMapping *cumap, float **array, int *size)
836 {
837         int a;
838         
839         *size = CM_TABLE + 1;
840         *array = MEM_callocN(sizeof(float) * (*size) * 4, "CurveMapping");
841         curvemapping_initialize(cumap);
842
843         for (a = 0; a < *size; a++) {
844                 if (cumap->cm[0].table)
845                         (*array)[a * 4 + 0] = cumap->cm[0].table[a].y;
846                 if (cumap->cm[1].table)
847                         (*array)[a * 4 + 1] = cumap->cm[1].table[a].y;
848                 if (cumap->cm[2].table)
849                         (*array)[a * 4 + 2] = cumap->cm[2].table[a].y;
850                 if (cumap->cm[3].table)
851                         (*array)[a * 4 + 3] = cumap->cm[3].table[a].y;
852         }
853 }
854
855 /* ***************** Histogram **************** */
856
857 #define INV_255     (1.f / 255.f)
858
859 DO_INLINE int get_bin_float(float f)
860 {
861         int bin = (int)((f * 255.0f) + 0.5f);  /* 0.5 to prevent quantisation differences */
862
863         /* note: clamp integer instead of float to avoid problems with NaN */
864         CLAMP(bin, 0, 255);
865
866         return bin;
867 }
868
869 DO_INLINE void save_sample_line(Scopes *scopes, const int idx, const float fx, const float rgb[3], const float ycc[3])
870 {
871         float yuv[3];
872
873         /* vectorscope*/
874         rgb_to_yuv(rgb[0], rgb[1], rgb[2], &yuv[0], &yuv[1], &yuv[2]);
875         scopes->vecscope[idx + 0] = yuv[1];
876         scopes->vecscope[idx + 1] = yuv[2];
877
878         /* waveform */
879         switch (scopes->wavefrm_mode) {
880                 case SCOPES_WAVEFRM_RGB:
881                         scopes->waveform_1[idx + 0] = fx;
882                         scopes->waveform_1[idx + 1] = rgb[0];
883                         scopes->waveform_2[idx + 0] = fx;
884                         scopes->waveform_2[idx + 1] = rgb[1];
885                         scopes->waveform_3[idx + 0] = fx;
886                         scopes->waveform_3[idx + 1] = rgb[2];
887                         break;
888                 case SCOPES_WAVEFRM_LUMA:
889                         scopes->waveform_1[idx + 0] = fx;
890                         scopes->waveform_1[idx + 1] = ycc[0];
891                         break;
892                 case SCOPES_WAVEFRM_YCC_JPEG:
893                 case SCOPES_WAVEFRM_YCC_709:
894                 case SCOPES_WAVEFRM_YCC_601:
895                         scopes->waveform_1[idx + 0] = fx;
896                         scopes->waveform_1[idx + 1] = ycc[0];
897                         scopes->waveform_2[idx + 0] = fx;
898                         scopes->waveform_2[idx + 1] = ycc[1];
899                         scopes->waveform_3[idx + 0] = fx;
900                         scopes->waveform_3[idx + 1] = ycc[2];
901                         break;
902         }
903 }
904
905 void scopes_update(Scopes *scopes, ImBuf *ibuf, int use_color_management)
906 {
907         int x, y, c;
908         unsigned int n, nl;
909         double div, divl;
910         float *rf = NULL;
911         unsigned char *rc = NULL;
912         unsigned int *bin_r, *bin_g, *bin_b, *bin_lum;
913         int savedlines, saveline;
914         float rgb[3], ycc[3], luma;
915         int ycc_mode = -1;
916         const short is_float = (ibuf->rect_float != NULL);
917
918         if (ibuf->rect == NULL && ibuf->rect_float == NULL) return;
919
920         if (scopes->ok == 1) return;
921
922         if (scopes->hist.ymax == 0.f) scopes->hist.ymax = 1.f;
923
924         /* hmmmm */
925         if (!(ELEM(ibuf->channels, 3, 4))) return;
926
927         scopes->hist.channels = 3;
928         scopes->hist.x_resolution = 256;
929
930         switch (scopes->wavefrm_mode) {
931                 case SCOPES_WAVEFRM_RGB:
932                         ycc_mode = -1;
933                         break;
934                 case SCOPES_WAVEFRM_LUMA:
935                 case SCOPES_WAVEFRM_YCC_JPEG:
936                         ycc_mode = BLI_YCC_JFIF_0_255;
937                         break;
938                 case SCOPES_WAVEFRM_YCC_601:
939                         ycc_mode = BLI_YCC_ITU_BT601;
940                         break;
941                 case SCOPES_WAVEFRM_YCC_709:
942                         ycc_mode = BLI_YCC_ITU_BT709;
943                         break;
944         }
945
946         /* temp table to count pix value for histo */
947         bin_r = MEM_callocN(256 * sizeof(unsigned int), "temp historgram bins");
948         bin_g = MEM_callocN(256 * sizeof(unsigned int), "temp historgram bins");
949         bin_b = MEM_callocN(256 * sizeof(unsigned int), "temp historgram bins");
950         bin_lum = MEM_callocN(256 * sizeof(unsigned int), "temp historgram bins");
951
952         /* convert to number of lines with logarithmic scale */
953         scopes->sample_lines = (scopes->accuracy * 0.01f) * (scopes->accuracy * 0.01f) * ibuf->y;
954         
955         if (scopes->sample_full)
956                 scopes->sample_lines = ibuf->y;
957
958         /* scan the image */
959         savedlines = 0;
960         for (c = 0; c < 3; c++) {
961                 scopes->minmax[c][0] = 25500.0f;
962                 scopes->minmax[c][1] = -25500.0f;
963         }
964         
965         scopes->waveform_tot = ibuf->x * scopes->sample_lines;
966         
967         if (scopes->waveform_1)
968                 MEM_freeN(scopes->waveform_1);
969         if (scopes->waveform_2)
970                 MEM_freeN(scopes->waveform_2);
971         if (scopes->waveform_3)
972                 MEM_freeN(scopes->waveform_3);
973         if (scopes->vecscope)
974                 MEM_freeN(scopes->vecscope);
975         
976         scopes->waveform_1 = MEM_callocN(scopes->waveform_tot * 2 * sizeof(float), "waveform point channel 1");
977         scopes->waveform_2 = MEM_callocN(scopes->waveform_tot * 2 * sizeof(float), "waveform point channel 2");
978         scopes->waveform_3 = MEM_callocN(scopes->waveform_tot * 2 * sizeof(float), "waveform point channel 3");
979         scopes->vecscope = MEM_callocN(scopes->waveform_tot * 2 * sizeof(float), "vectorscope point channel");
980         
981         if (is_float)
982                 rf = ibuf->rect_float;
983         else
984                 rc = (unsigned char *)ibuf->rect;
985
986         for (y = 0; y < ibuf->y; y++) {
987                 if (savedlines < scopes->sample_lines && y >= ((savedlines) * ibuf->y) / (scopes->sample_lines + 1)) {
988                         saveline = 1;
989                 }
990                 else {
991                         saveline = 0;
992                 }
993                 for (x = 0; x < ibuf->x; x++) {
994
995                         if (is_float) {
996                                 if (use_color_management)
997                                         linearrgb_to_srgb_v3_v3(rgb, rf);
998                                 else
999                                         copy_v3_v3(rgb, rf);
1000                         }
1001                         else {
1002                                 for (c = 0; c < 3; c++)
1003                                         rgb[c] = rc[c] * INV_255;
1004                         }
1005
1006                         /* we still need luma for histogram */
1007                         luma = rgb_to_luma(rgb);
1008
1009                         /* check for min max */
1010                         if (ycc_mode == -1) {
1011                                 for (c = 0; c < 3; c++) {
1012                                         if (rgb[c] < scopes->minmax[c][0]) scopes->minmax[c][0] = rgb[c];
1013                                         if (rgb[c] > scopes->minmax[c][1]) scopes->minmax[c][1] = rgb[c];
1014                                 }
1015                         }
1016                         else {
1017                                 rgb_to_ycc(rgb[0], rgb[1], rgb[2], &ycc[0], &ycc[1], &ycc[2], ycc_mode);
1018                                 for (c = 0; c < 3; c++) {
1019                                         ycc[c] *= INV_255;
1020                                         if (ycc[c] < scopes->minmax[c][0]) scopes->minmax[c][0] = ycc[c];
1021                                         if (ycc[c] > scopes->minmax[c][1]) scopes->minmax[c][1] = ycc[c];
1022                                 }
1023                         }
1024                         /* increment count for histo*/
1025                         bin_r[get_bin_float(rgb[0])] += 1;
1026                         bin_g[get_bin_float(rgb[1])] += 1;
1027                         bin_b[get_bin_float(rgb[2])] += 1;
1028                         bin_lum[get_bin_float(luma)] += 1;
1029
1030                         /* save sample if needed */
1031                         if (saveline) {
1032                                 const float fx = (float)x / (float)ibuf->x;
1033                                 const int idx = 2 * (ibuf->x * savedlines + x);
1034                                 save_sample_line(scopes, idx, fx, rgb, ycc);
1035                         }
1036
1037                         rf += ibuf->channels;
1038                         rc += ibuf->channels;
1039                 }
1040                 if (saveline)
1041                         savedlines += 1;
1042         }
1043
1044         /* convert hist data to float (proportional to max count) */
1045         n = 0;
1046         nl = 0;
1047         for (x = 0; x < 256; x++) {
1048                 if (bin_r[x] > n)
1049                         n = bin_r[x];
1050                 if (bin_g[x] > n)
1051                         n = bin_g[x];
1052                 if (bin_b[x] > n)
1053                         n = bin_b[x];
1054                 if (bin_lum[x] > nl)
1055                         nl = bin_lum[x];
1056         }
1057         div = 1.0 / (double)n;
1058         divl = 1.0 / (double)nl;
1059         for (x = 0; x < 256; x++) {
1060                 scopes->hist.data_r[x] = bin_r[x] * div;
1061                 scopes->hist.data_g[x] = bin_g[x] * div;
1062                 scopes->hist.data_b[x] = bin_b[x] * div;
1063                 scopes->hist.data_luma[x] = bin_lum[x] * divl;
1064         }
1065         MEM_freeN(bin_r);
1066         MEM_freeN(bin_g);
1067         MEM_freeN(bin_b);
1068         MEM_freeN(bin_lum);
1069
1070         scopes->ok = 1;
1071 }
1072
1073 void scopes_free(Scopes *scopes)
1074 {
1075         if (scopes->waveform_1) {
1076                 MEM_freeN(scopes->waveform_1);
1077                 scopes->waveform_1 = NULL;
1078         }
1079         if (scopes->waveform_2) {
1080                 MEM_freeN(scopes->waveform_2);
1081                 scopes->waveform_2 = NULL;
1082         }
1083         if (scopes->waveform_3) {
1084                 MEM_freeN(scopes->waveform_3);
1085                 scopes->waveform_3 = NULL;
1086         }
1087         if (scopes->vecscope) {
1088                 MEM_freeN(scopes->vecscope);
1089                 scopes->vecscope = NULL;
1090         }
1091 }
1092
1093 void scopes_new(Scopes *scopes)
1094 {
1095         scopes->accuracy = 30.0;
1096         scopes->hist.mode = HISTO_MODE_RGB;
1097         scopes->wavefrm_alpha = 0.3;
1098         scopes->vecscope_alpha = 0.3;
1099         scopes->wavefrm_height = 100;
1100         scopes->vecscope_height = 100;
1101         scopes->hist.height = 100;
1102         scopes->ok = 0;
1103         scopes->waveform_1 = NULL;
1104         scopes->waveform_2 = NULL;
1105         scopes->waveform_3 = NULL;
1106         scopes->vecscope = NULL;
1107 }