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