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