diff --git a/rtengine/color.cc b/rtengine/color.cc index d55e002e4..4c95ee5b2 100644 --- a/rtengine/color.cc +++ b/rtengine/color.cc @@ -1,21 +1,21 @@ /* - * This file is part of RawTherapee. - * - * Copyright (c) 2004-2010 Gabor Horvath - * - * RawTherapee is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * RawTherapee is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with RawTherapee. If not, see . - */ +* This file is part of RawTherapee. +* +* Copyright (c) 2004-2010 Gabor Horvath +* +* RawTherapee is free software: you can redistribute it and/or modify +* it under the terms of the GNU General Public License as published by +* the Free Software Foundation, either version 3 of the License, or +* (at your option) any later version. +* +* RawTherapee is distributed in the hope that it will be useful, +* but WITHOUT ANY WARRANTY; without even the implied warranty of +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +* GNU General Public License for more details. +* +* You should have received a copy of the GNU General Public License +* along with RawTherapee. If not, see . +*/ #include "rt_math.h" #include "color.h" @@ -25,330 +25,362 @@ using namespace std; namespace rtengine { -#define eps_max 580.40756 //(MAXVAL* 216.0f/24389.0); -#define kappa 903.29630 //24389.0/27.0; + LUTf Color::cachef; + LUTf Color::gamma2curve = 0; -LUTf Color::cachef ; -LUTf Color::gamma2curve = 0; + LUTf Color::gammatab; + LUTf Color::igammatab_srgb; + LUTf Color::gammatab_srgb; -LUTf Color::gammatab; -LUTf Color::igammatab_srgb; -LUTf Color::gammatab_srgb; + // Wikipedia sRGB: Unlike most other RGB color spaces, the sRGB gamma cannot be expressed as a single numerical value. + // The overall gamma is approximately 2.2, consisting of a linear (gamma 1.0) section near black, and a non-linear section elsewhere involving a 2.4 exponent + // and a gamma (slope of log output versus log input) changing from 1.0 through about 2.3. + const double Color::sRGBGamma = 2.2; + const double Color::sRGBGammaCurve = 2.4; -// Wikipedia sRGB: Unlike most other RGB color spaces, the sRGB gamma cannot be expressed as a single numerical value. -// The overall gamma is approximately 2.2, consisting of a linear (gamma 1.0) section near black, and a non-linear section elsewhere involving a 2.4 exponent -// and a gamma (slope of log output versus log input) changing from 1.0 through about 2.3. -const double Color::sRGBGamma = 2.2; -const double Color::sRGBGammaCurve = 2.4; + const double Color::eps_max=580.40756; //(MAXVAL* 216.0f/24389.0); + const double Color::kappa=903.29630; //24389.0/27.0; -void Color::init () { + const float Color::D50x=0.96422; + const float Color::D50z=0.82521; + const double Color::u0=4.0*D50x/(D50x+15+3*D50z); + const double Color::v0=9.0/(D50x+15+3*D50z); - int maxindex = 65536; - cachef(maxindex,0/*LUT_CLIP_BELOW*/); + void Color::init () { - gamma2curve(maxindex,0); + int maxindex = 65536; + cachef(maxindex,0/*LUT_CLIP_BELOW*/); - for (int i=0; ieps_max) { - cachef[i] = 327.68*( exp(1.0/3.0 * log((double)i / MAXVAL) )); - } - else { - cachef[i] = 327.68*((kappa*i/MAXVAL+16.0)/116.0); - } - } + gamma2curve(maxindex,0); - for (int i=0; ieps_max) { + cachef[i] = 327.68*( exp(1.0/3.0 * log((double)i / MAXVAL) )); + } + else { + cachef[i] = 327.68*((kappa*i/MAXVAL+16.0)/116.0); + } + } - /*******************************************/ + for (int i=0; i 1 ) h -= 1; - } -} - -void Color::rgb2hsv (int r, int g, int b, float &h, float &s, float &v) { - - double var_R = r / 65535.0; - double var_G = g / 65535.0; - double var_B = b / 65535.0; - - double var_Min = min(var_R,var_G,var_B); - double var_Max = max(var_R,var_G,var_B); - double del_Max = var_Max - var_Min; - v = var_Max; - if (fabs(del_Max)<0.00001) { - h = 0; - s = 0; - } - else { - s = del_Max/var_Max; - - if ( var_R == var_Max ) h = (var_G - var_B)/del_Max; - else if ( var_G == var_Max ) h = 2.0 + (var_B - var_R)/del_Max; - else if ( var_B == var_Max ) h = 4.0 + (var_R - var_G)/del_Max; - h /= 6.0; - - if ( h < 0 ) h += 1; - if ( h > 1 ) h -= 1; - } -} - -void Color::hsv2rgb (float h, float s, float v, float &r, float &g, float &b) { - - float h1 = h*6; // sector 0 to 5 - int i = floor( h1 ); - float f = h1 - i; // fractional part of h - - float p = v * ( 1 - s ); - float q = v * ( 1 - s * f ); - float t = v * ( 1 - s * ( 1 - f ) ); - - float r1,g1,b1; - - if (i==0) {r1 = v; g1 = t; b1 = p;} - else if (i==1) {r1 = q; g1 = v; b1 = p;} - else if (i==2) {r1 = p; g1 = v; b1 = t;} - else if (i==3) {r1 = p; g1 = q; b1 = v;} - else if (i==4) {r1 = t; g1 = p; b1 = v;} - else if (i==5) {r1 = v; g1 = p; b1 = q;} - - r = ((r1)*65535.0); - g = ((g1)*65535.0); - b = ((b1)*65535.0); -} - - -// The same function but set float values instead if int -// Function copied for speed concerns -/* Not exactly the same as above ; this one return a result in the [0.0 ; 1.0] range -void Color::hsv2rgb (float h, float s, float v, float &r, float &g, float &b) { - - float h1 = h*6; // sector 0 to 5 - int i = floor( h1 ); - float f = h1 - i; // fractional part of h - - float p = v * ( 1 - s ); - float q = v * ( 1 - s * f ); - float t = v * ( 1 - s * ( 1 - f ) ); - - if (i==0) {r = v; g = t; b = p;} - else if (i==1) {r = q; g = v; b = p;} - else if (i==2) {r = p; g = v; b = t;} - else if (i==3) {r = p; g = q; b = v;} - else if (i==4) {r = t; g = p; b = v;} - else if (i==5) {r = v; g = p; b = q;} -} -*/ - -void Color::hsv2rgb (float h, float s, float v, int &r, int &g, int &b) { - - float h1 = h*6; // sector 0 to 5 - int i = floor( h1 ); - float f = h1 - i; // fractional part of h - - float p = v * ( 1 - s ); - float q = v * ( 1 - s * f ); - float t = v * ( 1 - s * ( 1 - f ) ); - - float r1,g1,b1; - - if (i==0) {r1 = v; g1 = t; b1 = p;} - else if (i==1) {r1 = q; g1 = v; b1 = p;} - else if (i==2) {r1 = p; g1 = v; b1 = t;} - else if (i==3) {r1 = p; g1 = q; b1 = v;} - else if (i==4) {r1 = t; g1 = p; b1 = v;} - else if (i==5) {r1 = v; g1 = p; b1 = q;} - - r = (int)( r1 * 65535); - g = (int)( g1 * 65535); - b = (int)( b1 * 65535); -} - -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -void Color::xyz2srgb (float x, float y, float z, float &r, float &g, float &b) { - - //Transform to output color. Standard sRGB is D65, but internal representation is D50 - //Note that it is only at this point that we should have need of clipping color data - - /*float x65 = d65_d50[0][0]*x + d65_d50[0][1]*y + d65_d50[0][2]*z ; - float y65 = d65_d50[1][0]*x + d65_d50[1][1]*y + d65_d50[1][2]*z ; - float z65 = d65_d50[2][0]*x + d65_d50[2][1]*y + d65_d50[2][2]*z ; - - r = sRGB_xyz[0][0]*x65 + sRGB_xyz[0][1]*y65 + sRGB_xyz[0][2]*z65; - g = sRGB_xyz[1][0]*x65 + sRGB_xyz[1][1]*y65 + sRGB_xyz[1][2]*z65; - b = sRGB_xyz[2][0]*x65 + sRGB_xyz[2][1]*y65 + sRGB_xyz[2][2]*z65;*/ - - /*r = sRGBd65_xyz[0][0]*x + sRGBd65_xyz[0][1]*y + sRGBd65_xyz[0][2]*z ; - g = sRGBd65_xyz[1][0]*x + sRGBd65_xyz[1][1]*y + sRGBd65_xyz[1][2]*z ; - b = sRGBd65_xyz[2][0]*x + sRGBd65_xyz[2][1]*y + sRGBd65_xyz[2][2]*z ;*/ - - r = ((sRGB_xyz[0][0]*x + sRGB_xyz[0][1]*y + sRGB_xyz[0][2]*z)) ; - g = ((sRGB_xyz[1][0]*x + sRGB_xyz[1][1]*y + sRGB_xyz[1][2]*z)) ; - b = ((sRGB_xyz[2][0]*x + sRGB_xyz[2][1]*y + sRGB_xyz[2][2]*z)) ; - -} - - -void Color::xyz2rgb (float x, float y, float z, float &r, float &g, float &b, float rgb_xyz[3][3]) { - - //Transform to output color. Standard sRGB is D65, but internal representation is D50 - //Note that it is only at this point that we should have need of clipping color data - - /*float x65 = d65_d50[0][0]*x + d65_d50[0][1]*y + d65_d50[0][2]*z ; - float y65 = d65_d50[1][0]*x + d65_d50[1][1]*y + d65_d50[1][2]*z ; - float z65 = d65_d50[2][0]*x + d65_d50[2][1]*y + d65_d50[2][2]*z ; - - r = sRGB_xyz[0][0]*x65 + sRGB_xyz[0][1]*y65 + sRGB_xyz[0][2]*z65; - g = sRGB_xyz[1][0]*x65 + sRGB_xyz[1][1]*y65 + sRGB_xyz[1][2]*z65; - b = sRGB_xyz[2][0]*x65 + sRGB_xyz[2][1]*y65 + sRGB_xyz[2][2]*z65;*/ - - /*r = sRGBd65_xyz[0][0]*x + sRGBd65_xyz[0][1]*y + sRGBd65_xyz[0][2]*z ; - g = sRGBd65_xyz[1][0]*x + sRGBd65_xyz[1][1]*y + sRGBd65_xyz[1][2]*z ; - b = sRGBd65_xyz[2][0]*x + sRGBd65_xyz[2][1]*y + sRGBd65_xyz[2][2]*z ;*/ - - - - r = ((rgb_xyz[0][0]*x + rgb_xyz[0][1]*y + rgb_xyz[0][2]*z)) ; - g = ((rgb_xyz[1][0]*x + rgb_xyz[1][1]*y + rgb_xyz[1][2]*z)) ; - b = ((rgb_xyz[2][0]*x + rgb_xyz[2][1]*y + rgb_xyz[2][2]*z)) ; - -} - -void Color::calcGamma (double pwr, double ts, int mode, int imax, double &gamma0, double &gamma1, double &gamma2, double &gamma3, double &gamma4, double &gamma5) { -//from Dcraw (D.Coffin) - int i; - double g[6], bnd[2]={0,0}; - - g[0] = pwr; - g[1] = ts; - g[2] = g[3] = g[4] = 0; - bnd[g[1] >= 1] = 1; - if (g[1] && (g[1]-1)*(g[0]-1) <= 0) { - for (i=0; i < 48; i++) { - g[2] = (bnd[0] + bnd[1])/2; - if (g[0]) - bnd[(pow(g[2]/g[1],-g[0]) - 1)/g[0] - 1/g[2] > -1] = g[2]; - else - bnd[g[2]/exp(1-1/g[2]) < g[1]] = g[2]; - } - g[3] = g[2] / g[1]; - if (g[0]) g[4] = g[2] * (1/g[0] - 1); - } - if (g[0]) - g[5] = 1 / (g[1]*SQR(g[3])/2 - g[4]*(1 - g[3]) + (1 - pow(g[3],1+g[0]))*(1 + g[4])/(1 + g[0])) - 1; - else - g[5] = 1 / (g[1]*SQR(g[3])/2 + 1 - g[2] - g[3] - g[2]*g[3]*(log(g[3]) - 1)) - 1; - if (!mode--) { - gamma0=g[0];gamma1=g[1];gamma2=g[2];gamma3=g[3];gamma4=g[4];gamma5=g[5]; - return; - } -} - -void Color::Lab2XYZ(float L, float a, float b, float &x, float &y, float &z) { - float fy = (0.00862069 * L) + 0.137932; // (L+16)/116 - float fx = (0.002 * a) + fy; - float fz = fy - (0.005 * b); - - x = 65535.0*f2xyz(fx)*D50x; - y = 65535.0*f2xyz(fy); - z = 65535.0*f2xyz(fz)*D50z; -} - -void Color::XYZ2Lab(float X, float Y, float Z, float &L, float &a, float &b) { - - float X1 = X/D50x; - float Z1 = Z/D50z; - - float fx = (X1<65535.0 ? cachef[X1] : (327.68*exp(log(X1/MAXVAL)/3.0 ))); - float fy = (Y<65535.0 ? cachef[Y] : (327.68*exp(log(Y/MAXVAL)/3.0 ))); - float fz = (Z1<65535.0 ? cachef[Z1] : (327.68*exp(log(Z1/MAXVAL)/3.0 ))); - - L = (116.0 * fy - 5242.88); //5242.88=16.0*327.68; - a = (500.0 * (fx - fy) ); - b = (200.0 * (fy - fz) ); - -} - -void Color::Lab2Yuv(float L, float a, float b, float &Y, float &u, float &v) { - float fy = (0.00862069 * L/327.68) + 0.137932; // (L+16)/116 - float fx = (0.002 * a/327.68) + fy; - float fz = fy - (0.005 * b/327.68); - - float X = 65535.0*f2xyz(fx)*D50x; - Y = 65535.0*f2xyz(fy); - float Z = 65535.0*f2xyz(fz)*D50z; - - u = 4.0*X/(X+15*Y+3*Z)-u0; - v = 9.0*Y/(X+15*Y+3*Z)-v0; - -} - -void Color::Yuv2Lab(float Yin, float u, float v, float &L, float &a, float &b, double wp[3][3]) { - - float u1 = u + u0; - float v1 = v + v0; - - float Y = Yin; - float X = (9*u1*Y)/(4*v1*D50x); - float Z = (12 - 3*u1 - 20*v1)*Y/(4*v1*D50z); - - gamutmap(X,Y,Z,wp); - - float fx = (X<65535.0 ? cachef[X] : (327.68*exp(log(X/MAXVAL)/3.0 ))); - float fy = (Y<65535.0 ? cachef[Y] : (327.68*exp(log(Y/MAXVAL)/3.0 ))); - float fz = (Z<65535.0 ? cachef[Z] : (327.68*exp(log(Z/MAXVAL)/3.0 ))); - - L = (116.0 * fy - 5242.88); //5242.88=16.0*327.68; - a = (500.0 * (fx - fy) ); - b = (200.0 * (fy - fz) ); - -} + for (int i=0; i<65536; i++) + gammatab_srgb[i] = (65535.0 * gamma2 (i/65535.0)); + for (int i=0; i<65536; i++) + igammatab_srgb[i] = (65535.0 * igamma2 (i/65535.0)); + for (int i=0; i<65536; i++) + gammatab[i] = (65535.0 * pow (i/65535.0, 0.454545)); + /*FILE* f = fopen ("c.txt", "wt"); + for (int i=0; i<256; i++) + fprintf (f, "%g %g\n", i/255.0, clower (i/255.0, 2.0, 1.0)); + fclose (f);*/ + + } + + void Color::cleanup () { + + } + + void Color::rgb2hsv(float r, float g, float b, float &h, float &s, float &v) { + double var_R = r / 65535.0; + double var_G = g / 65535.0; + double var_B = b / 65535.0; + + double var_Min = min(var_R,var_G,var_B); + double var_Max = max(var_R,var_G,var_B); + double del_Max = var_Max - var_Min; + v = var_Max; + if (del_Max<0.00001 && del_Max>-0.00001) { // no fabs, slow! + h = 0; + s = 0; + } + else { + s = del_Max/var_Max; + + if ( var_R == var_Max ) h = (var_G - var_B)/del_Max; + else if ( var_G == var_Max ) h = 2.0 + (var_B - var_R)/del_Max; + else if ( var_B == var_Max ) h = 4.0 + (var_R - var_G)/del_Max; + h /= 6.0; + + if ( h < 0 ) h += 1; + if ( h > 1 ) h -= 1; + } + } + + void Color::hsv2rgb (float h, float s, float v, float &r, float &g, float &b) { + + float h1 = h*6; // sector 0 to 5 + int i = (int)h1; // floor() is very slow, and h1 is always >0 + float f = h1 - i; // fractional part of h + + float p = v * ( 1 - s ); + float q = v * ( 1 - s * f ); + float t = v * ( 1 - s * ( 1 - f ) ); + + float r1,g1,b1; + + if (i==1) {r1 = q; g1 = v; b1 = p;} + else if (i==2) {r1 = p; g1 = v; b1 = t;} + else if (i==3) {r1 = p; g1 = q; b1 = v;} + else if (i==4) {r1 = t; g1 = p; b1 = v;} + else if (i==5) {r1 = v; g1 = p; b1 = q;} + else /*i==(0|6)*/ {r1 = v; g1 = t; b1 = p;} + + r = ((r1)*65535.0); + g = ((g1)*65535.0); + b = ((b1)*65535.0); + } + + // Function copied for speed concerns + // Not exactly the same as above ; this one return a result in the [0.0 ; 1.0] range + void Color::hsv2rgb01 (float h, float s, float v, float &r, float &g, float &b) { + float h1 = h*6; // sector 0 to 5 + int i = int(h1); + float f = h1 - i; // fractional part of h + + float p = v * ( 1 - s ); + float q = v * ( 1 - s * f ); + float t = v * ( 1 - s * ( 1 - f ) ); + + if (i==1) {r = q; g = v; b = p;} + else if (i==2) {r = p; g = v; b = t;} + else if (i==3) {r = p; g = q; b = v;} + else if (i==4) {r = t; g = p; b = v;} + else if (i==5) {r = v; g = p; b = q;} + else /*(i==0|6)*/ {r = v; g = t; b = p;} + } + + void Color::hsv2rgb (float h, float s, float v, int &r, int &g, int &b) { + + float h1 = h*6; // sector 0 to 5 + int i = floor( h1 ); + float f = h1 - i; // fractional part of h + + float p = v * ( 1 - s ); + float q = v * ( 1 - s * f ); + float t = v * ( 1 - s * ( 1 - f ) ); + + float r1,g1,b1; + + if (i==0) {r1 = v; g1 = t; b1 = p;} + else if (i==1) {r1 = q; g1 = v; b1 = p;} + else if (i==2) {r1 = p; g1 = v; b1 = t;} + else if (i==3) {r1 = p; g1 = q; b1 = v;} + else if (i==4) {r1 = t; g1 = p; b1 = v;} + else if (i==5) {r1 = v; g1 = p; b1 = q;} + + r = (int)( r1 * 65535); + g = (int)( g1 * 65535); + b = (int)( b1 * 65535); + } + + //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + void Color::xyz2srgb (float x, float y, float z, float &r, float &g, float &b) { + + //Transform to output color. Standard sRGB is D65, but internal representation is D50 + //Note that it is only at this point that we should have need of clipping color data + + /*float x65 = d65_d50[0][0]*x + d65_d50[0][1]*y + d65_d50[0][2]*z ; + float y65 = d65_d50[1][0]*x + d65_d50[1][1]*y + d65_d50[1][2]*z ; + float z65 = d65_d50[2][0]*x + d65_d50[2][1]*y + d65_d50[2][2]*z ; + + r = sRGB_xyz[0][0]*x65 + sRGB_xyz[0][1]*y65 + sRGB_xyz[0][2]*z65; + g = sRGB_xyz[1][0]*x65 + sRGB_xyz[1][1]*y65 + sRGB_xyz[1][2]*z65; + b = sRGB_xyz[2][0]*x65 + sRGB_xyz[2][1]*y65 + sRGB_xyz[2][2]*z65;*/ + + /*r = sRGBd65_xyz[0][0]*x + sRGBd65_xyz[0][1]*y + sRGBd65_xyz[0][2]*z ; + g = sRGBd65_xyz[1][0]*x + sRGBd65_xyz[1][1]*y + sRGBd65_xyz[1][2]*z ; + b = sRGBd65_xyz[2][0]*x + sRGBd65_xyz[2][1]*y + sRGBd65_xyz[2][2]*z ;*/ + + r = ((sRGB_xyz[0][0]*x + sRGB_xyz[0][1]*y + sRGB_xyz[0][2]*z)) ; + g = ((sRGB_xyz[1][0]*x + sRGB_xyz[1][1]*y + sRGB_xyz[1][2]*z)) ; + b = ((sRGB_xyz[2][0]*x + sRGB_xyz[2][1]*y + sRGB_xyz[2][2]*z)) ; + + } + + void Color::xyz2rgb (float x, float y, float z, float &r, float &g, float &b, double rgb_xyz[3][3]) { + //Transform to output color. Standard sRGB is D65, but internal representation is D50 + //Note that it is only at this point that we should have need of clipping color data + + /*float x65 = d65_d50[0][0]*x + d65_d50[0][1]*y + d65_d50[0][2]*z ; + float y65 = d65_d50[1][0]*x + d65_d50[1][1]*y + d65_d50[1][2]*z ; + float z65 = d65_d50[2][0]*x + d65_d50[2][1]*y + d65_d50[2][2]*z ; + + r = sRGB_xyz[0][0]*x65 + sRGB_xyz[0][1]*y65 + sRGB_xyz[0][2]*z65; + g = sRGB_xyz[1][0]*x65 + sRGB_xyz[1][1]*y65 + sRGB_xyz[1][2]*z65; + b = sRGB_xyz[2][0]*x65 + sRGB_xyz[2][1]*y65 + sRGB_xyz[2][2]*z65;*/ + + /*r = sRGBd65_xyz[0][0]*x + sRGBd65_xyz[0][1]*y + sRGBd65_xyz[0][2]*z ; + g = sRGBd65_xyz[1][0]*x + sRGBd65_xyz[1][1]*y + sRGBd65_xyz[1][2]*z ; + b = sRGBd65_xyz[2][0]*x + sRGBd65_xyz[2][1]*y + sRGBd65_xyz[2][2]*z ;*/ + + r = ((rgb_xyz[0][0]*x + rgb_xyz[0][1]*y + rgb_xyz[0][2]*z)) ; + g = ((rgb_xyz[1][0]*x + rgb_xyz[1][1]*y + rgb_xyz[1][2]*z)) ; + b = ((rgb_xyz[2][0]*x + rgb_xyz[2][1]*y + rgb_xyz[2][2]*z)) ; + } + + // same for float + void Color::xyz2rgb (float x, float y, float z, float &r, float &g, float &b, float rgb_xyz[3][3]) { + r = ((rgb_xyz[0][0]*x + rgb_xyz[0][1]*y + rgb_xyz[0][2]*z)) ; + g = ((rgb_xyz[1][0]*x + rgb_xyz[1][1]*y + rgb_xyz[1][2]*z)) ; + b = ((rgb_xyz[2][0]*x + rgb_xyz[2][1]*y + rgb_xyz[2][2]*z)) ; + } + + void Color::calcGamma (double pwr, double ts, int mode, int imax, double &gamma0, double &gamma1, double &gamma2, double &gamma3, double &gamma4, double &gamma5) { + //from Dcraw (D.Coffin) + int i; + double g[6], bnd[2]={0,0}; + + g[0] = pwr; + g[1] = ts; + g[2] = g[3] = g[4] = 0; + bnd[g[1] >= 1] = 1; + if (g[1] && (g[1]-1)*(g[0]-1) <= 0) { + for (i=0; i < 48; i++) { + g[2] = (bnd[0] + bnd[1])/2; + if (g[0]) + bnd[(pow(g[2]/g[1],-g[0]) - 1)/g[0] - 1/g[2] > -1] = g[2]; + else + bnd[g[2]/exp(1-1/g[2]) < g[1]] = g[2]; + } + g[3] = g[2] / g[1]; + if (g[0]) g[4] = g[2] * (1/g[0] - 1); + } + if (g[0]) + g[5] = 1 / (g[1]*SQR(g[3])/2 - g[4]*(1 - g[3]) + (1 - pow(g[3],1+g[0]))*(1 + g[4])/(1 + g[0])) - 1; + else + g[5] = 1 / (g[1]*SQR(g[3])/2 + 1 - g[2] - g[3] - g[2]*g[3]*(log(g[3]) - 1)) - 1; + if (!mode--) { + gamma0=g[0];gamma1=g[1];gamma2=g[2];gamma3=g[3];gamma4=g[4];gamma5=g[5]; + return; + } + } + + void Color::Lab2XYZ(float L, float a, float b, float &x, float &y, float &z) { + float fy = (0.00862069 * L) + 0.137932; // (L+16)/116 + float fx = (0.002 * a) + fy; + float fz = fy - (0.005 * b); + + x = 65535.0*f2xyz(fx)*D50x; + y = 65535.0*f2xyz(fy); + z = 65535.0*f2xyz(fz)*D50z; + } + + void Color::XYZ2Lab(float X, float Y, float Z, float &L, float &a, float &b) { + + float X1 = X/D50x; + float Z1 = Z/D50z; + + float fx = (X1<65535.0 ? cachef[X1] : (327.68*exp(log(X1/MAXVAL)/3.0 ))); + float fy = (Y<65535.0 ? cachef[Y] : (327.68*exp(log(Y/MAXVAL)/3.0 ))); + float fz = (Z1<65535.0 ? cachef[Z1] : (327.68*exp(log(Z1/MAXVAL)/3.0 ))); + + L = (116.0 * fy - 5242.88); //5242.88=16.0*327.68; + a = (500.0 * (fx - fy) ); + b = (200.0 * (fy - fz) ); + + } + + void Color::Lab2Yuv(float L, float a, float b, float &Y, float &u, float &v) { + float fy = (0.00862069 * L/327.68) + 0.137932; // (L+16)/116 + float fx = (0.002 * a/327.68) + fy; + float fz = fy - (0.005 * b/327.68); + + float X = 65535.0*f2xyz(fx)*D50x; + Y = 65535.0*f2xyz(fy); + float Z = 65535.0*f2xyz(fz)*D50z; + + u = 4.0*X/(X+15*Y+3*Z)-u0; + v = 9.0*Y/(X+15*Y+3*Z)-v0; + } + + void Color::Yuv2Lab(float Yin, float u, float v, float &L, float &a, float &b, double wp[3][3]) { + float u1 = u + u0; + float v1 = v + v0; + + float Y = Yin; + float X = (9*u1*Y)/(4*v1*D50x); + float Z = (12 - 3*u1 - 20*v1)*Y/(4*v1*D50z); + + gamutmap(X,Y,Z,wp); + + float fx = (X<65535.0 ? cachef[X] : (327.68*exp(log(X/MAXVAL)/3.0 ))); + float fy = (Y<65535.0 ? cachef[Y] : (327.68*exp(log(Y/MAXVAL)/3.0 ))); + float fz = (Z<65535.0 ? cachef[Z] : (327.68*exp(log(Z/MAXVAL)/3.0 ))); + + L = (116.0 * fy - 5242.88); //5242.88=16.0*327.68; + a = (500.0 * (fx - fy) ); + b = (200.0 * (fy - fz) ); + } + + double Color::f2xyz(double f) { + const double epsilonExpInv3 = 6.0/29.0; + const double kappaInv = 27.0/24389.0; // inverse of kappa + + return (f > epsilonExpInv3) ? f*f*f : (116 * f - 16) * kappaInv; + } + +/* Gamut mapping algorithm + copyright (c) 2010-2011 Emil Martinec + solutions to scaling u and v to XYZ paralleliped boundaries + some equations: + fu(X,Y,Z) = 4 X/(X + 15 Y + 3 Z); + fv(X,Y,Z) = 9 Y/(X + 15 Y + 3 Z); + + take the plane spanned by X=a*Xr+b*Xg+c*Xb etc with one of a,b,c equal to 0 or 1, + and itersect with the line u0+lam*u, or in other words solve + + u0+lam*u=fu(X,Y,Z) + v0+lam*v=fv(X,Y,Z) + + the value of lam is the scale factor that takes the color to the gamut boundary + + columns of the matrix p=xyz_rgb are RGB tristimulus primaries in XYZ + c is the color fixed on the boundary; and m=0 for c=0, m=1 for c=255 +*/ + + void Color::gamutmap(float &X, float &Y, float &Z, const double p[3][3]) + { + float u = 4*X/(X+15*Y+3*Z)-u0; + float v = 9*Y/(X+15*Y+3*Z)-v0; + + float lam[3][2]; + float lam_min = 1.0; + + for (int c=0; c<3; c++) + for (int m=0; m<2; m++) { + + int c1=(c+1)%3; + int c2=(c+2)%3; + + lam[c][m] = (-(p[0][c1]*p[1][c]*((-12 + 3*u0 + 20*v0)*Y + 4*m*65535*v0*p[2][c2])) + + p[0][c]*p[1][c1]*((-12 + 3*u0 + 20*v0)*Y + 4*m*65535*v0*p[2][c2]) - + 4*v0*p[0][c1]*(Y - m*65535*p[1][c2])*p[2][c] + 4*v0*p[0][c]*(Y - m*65535*p[1][c2])*p[2][c1] - + (4*m*65535*v0*p[0][c2] - 9*u0*Y)*(p[1][c1]*p[2][c] - p[1][c]*p[2][c1])); + + lam[c][m] /= (3*u*Y*(p[0][c1]*p[1][c] - p[1][c1]*(p[0][c] + 3*p[2][c]) + 3*p[1][c]*p[2][c1]) + + 4*v*(p[0][c1]*(5*Y*p[1][c] + m*65535*p[1][c]*p[2][c2] + Y*p[2][c] - m*65535*p[1][c2]*p[2][c]) - + p[0][c]*(5*Y*p[1][c1] + m*65535*p[1][c1]*p[2][c2] + Y*p[2][c1] - m*65535*p[1][c2]*p[2][c1]) + + m*65535*p[0][c2]*(p[1][c1]*p[2][c] - p[1][c]*p[2][c1]))); + + if (lam[c][m]0) lam_min=lam[c][m]; + + } + + u = u*lam_min + u0; + v = v*lam_min + v0; + + X = (9*u*Y)/(4*v); + Z = (12 - 3*u - 20*v)*Y/(4*v); + } } diff --git a/rtengine/color.h b/rtengine/color.h index a066a78a0..7b6a3b57a 100644 --- a/rtengine/color.h +++ b/rtengine/color.h @@ -23,11 +23,6 @@ #include #include "LUT.h" -#define D50x 0.96422 -#define D50z 0.82521 -#define u0 4.0*D50x/(D50x+15+3*D50z) -#define v0 9.0/(D50x+15+3*D50z) - namespace rtengine { class Color { @@ -35,6 +30,9 @@ class Color { public: const static double sRGBGamma; // standard average gamma const static double sRGBGammaCurve; // 2.4 in the curve + const static double eps_max, kappa; + const static float D50x, D50z; + const static double u0, v0; static LUTf cachef; static LUTf gamma2curve; @@ -50,15 +48,17 @@ public: static void cleanup (); static void rgb2hsv (float r, float g, float b, float &h, float &s, float &v); - static void rgb2hsv (int r, int g, int b, float &h, float &s, float &v); static void hsv2rgb (float h, float s, float v, float &r, float &g, float &b); static void hsv2rgb (float h, float s, float v, int &r, int &g, int &b); + static void hsv2rgb01 (float h, float s, float v, float &r, float &g, float &b); static void xyz2srgb (float x, float y, float z, float &r, float &g, float &b); - static void xyz2rgb (float x, float y, float z, float &r, float &g, float &b, float rgb_xyz[3][3]); + static void xyz2rgb (float x, float y, float z, float &r, float &g, float &b, double rgb_xyz[3][3]); + static void xyz2rgb (float x, float y, float z, float &r, float &g, float &b, float rgb_xyz[3][3]); static void Lab2XYZ(float L, float a, float b, float &x, float &y, float &z); static void XYZ2Lab(float X, float Y, float Z, float &L, float &a, float &b); static void Lab2Yuv(float L, float a, float b, float &Y, float &u, float &v); static void Yuv2Lab(float Y, float u, float v, float &L, float &a, float &b, double wp[3][3]); + static double f2xyz(double f); static void calcGamma (double pwr, double ts, int mode, int imax, double &gamma0, double &gamma1, double &gamma2, double &gamma3, double &gamma4,double &gamma5); // standard srgb gamma and its inverse diff --git a/rtengine/dcp.cc b/rtengine/dcp.cc index 60acca23f..a490b0090 100644 --- a/rtengine/dcp.cc +++ b/rtengine/dcp.cc @@ -271,12 +271,12 @@ void DCPProfile::Apply(Imagefloat *pImg, DCPLightType preferredProfile, Glib::us // if point is in negative area, just the matrix, but not the LUT if (newr>=0 && newg>=0 && newb>=0) { - ImProcFunctions::rgb2hsv(newr, newg, newb, h , s, v); + Color::rgb2hsv(newr, newg, newb, h , s, v); h*=6.f; // RT calculates in [0,1] if (useRawWhite) { // Retro-calculate what the point was like before RAW white came in - ImProcFunctions::rgb2hsv(newr/rawWhiteFac, newg/rawWhiteFac, newb/rawWhiteFac, hs, ss, vs); + Color::rgb2hsv(newr/rawWhiteFac, newg/rawWhiteFac, newb/rawWhiteFac, hs, ss, vs); hs*=6.f; // RT calculates in [0,1] } else { hs=h; ss=s; vs=v; @@ -422,7 +422,7 @@ void DCPProfile::Apply(Imagefloat *pImg, DCPLightType preferredProfile, Glib::us if (h < 0.0f) h += 6.0f; if (h >= 6.0f) h -= 6.0f; h/=6.f; - ImProcFunctions::hsv2rgb( h, s, v, newr, newg, newb); + Color::hsv2rgb( h, s, v, newr, newg, newb); } pImg->r[y][x] = m2Work[0][0]*newr + m2Work[0][1]*newg + m2Work[0][2]*newb; diff --git a/rtengine/dcrop.cc b/rtengine/dcrop.cc index 409358184..41f63d635 100644 --- a/rtengine/dcrop.cc +++ b/rtengine/dcrop.cc @@ -401,7 +401,7 @@ void Crop::fullUpdate () { parent->updaterThreadStart.lock (); if (parent->updaterRunning && parent->thread) { // Do NOT reset changes here, since in a long chain of events it will lead to chroma_scale not being updated, - // causing ImProcFunctions::lab2rgb to return a black image on some opens + // causing Color::lab2rgb to return a black image on some opens //parent->changeSinceLast = 0; parent->thread->join (); } diff --git a/rtengine/gamutbdy.cc b/rtengine/gamutbdy.cc index e815f3816..2e3ba23d9 100644 --- a/rtengine/gamutbdy.cc +++ b/rtengine/gamutbdy.cc @@ -1,97 +1 @@ -//////////////////////////////////////////////////////////////// -// -// Gamut mapping algorithm -// -// copyright (c) 2010-2011 Emil Martinec -// -// -// code dated: February 2, 2011 -// -// sRGBgamutbdy.cc is free software: you can redistribute it and/or modify -// it under the terms of the GNU General Public License as published by -// the Free Software Foundation, either version 3 of the License, or -// (at your option) any later version. -// -// This program is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -// GNU General Public License for more details. -// -// You should have received a copy of the GNU General Public License -// along with this program. If not, see . -// -//////////////////////////////////////////////////////////////// -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -#include "rt_math.h" - -#define u0 4.0*D50x/(D50x+15+3*D50z) -#define v0 9.0/(D50x+15+3*D50z) - -#define sgn(x) (x<0 ? -1 : 1) -#define sqrt2 (sqrt(2)) -#define sqrt3 (sqrt(3)) -#define sqrt6 (sqrt(6)) -#define K 20 -#define rmax (3*K + 1) -#define rctr (2*K) -#define cctr (2*K + 1) -#define maxindx (2*(3*K+1) + SQR(3*K+1)) - - - -// solutions to scaling u and v to XYZ paralleliped boundaries -/* some equations: - fu(X,Y,Z) = 4 X/(X + 15 Y + 3 Z); - fv(X,Y,Z) = 9 Y/(X + 15 Y + 3 Z); - - take the plane spanned by X=a*Xr+b*Xg+c*Xb etc with one of a,b,c equal to 0 or 1, - and itersect with the line u0+lam*u, or in other words solve - - u0+lam*u=fu(X,Y,Z) - v0+lam*v=fv(X,Y,Z) - - the value of lam is the scale factor that takes the color to the gamut boundary -*/ -// columns of the matrix p=xyz_rgb are RGB tristimulus primaries in XYZ -// c is the color fixed on the boundary; and m=0 for c=0, m=1 for c=255 - -void Color::gamutmap(float &X, float &Y, float &Z, const double p[3][3]) -{ - float u = 4*X/(X+15*Y+3*Z)-u0; - float v = 9*Y/(X+15*Y+3*Z)-v0; - - float lam[3][2]; - float lam_min = 1.0; - - for (int c=0; c<3; c++) - for (int m=0; m<2; m++) { - - int c1=(c+1)%3; - int c2=(c+2)%3; - - lam[c][m] = (-(p[0][c1]*p[1][c]*((-12 + 3*u0 + 20*v0)*Y + 4*m*65535*v0*p[2][c2])) + - p[0][c]*p[1][c1]*((-12 + 3*u0 + 20*v0)*Y + 4*m*65535*v0*p[2][c2]) - - 4*v0*p[0][c1]*(Y - m*65535*p[1][c2])*p[2][c] + 4*v0*p[0][c]*(Y - m*65535*p[1][c2])*p[2][c1] - - (4*m*65535*v0*p[0][c2] - 9*u0*Y)*(p[1][c1]*p[2][c] - p[1][c]*p[2][c1])); - - lam[c][m] /= (3*u*Y*(p[0][c1]*p[1][c] - p[1][c1]*(p[0][c] + 3*p[2][c]) + 3*p[1][c]*p[2][c1]) + - 4*v*(p[0][c1]*(5*Y*p[1][c] + m*65535*p[1][c]*p[2][c2] + Y*p[2][c] - m*65535*p[1][c2]*p[2][c]) - - p[0][c]*(5*Y*p[1][c1] + m*65535*p[1][c1]*p[2][c2] + Y*p[2][c1] - m*65535*p[1][c2]*p[2][c1]) + - m*65535*p[0][c2]*(p[1][c1]*p[2][c] - p[1][c]*p[2][c1]))); - - if (lam[c][m]0) lam_min=lam[c][m]; - - } - - u = u*lam_min + u0; - v = v*lam_min + v0; - - X = (9*u*Y)/(4*v); - Z = (12 - 3*u - 20*v)*Y/(4*v); - - -} - - -//};//namespace +// File will be deleted on commit \ No newline at end of file diff --git a/rtengine/improcfun.cc b/rtengine/improcfun.cc index 4b1f497d8..50ab3df84 100644 --- a/rtengine/improcfun.cc +++ b/rtengine/improcfun.cc @@ -56,13 +56,7 @@ namespace rtengine { #define CLIPC(a) ((a)>-32000?((a)<32000?(a):32000):-32000) #define CLIP2(a) ((a)0.0?((a)<65535.5?(a):65535.5):0.0) - -#define u0 4.0*D50x/(D50x+15+3*D50z) -#define v0 9.0/(D50x+15+3*D50z) - -#define eps_max 580.40756 //(MAXVAL* 216.0f/24389.0); -#define kappa 903.29630 //24389.0/27.0; - + extern const Settings* settings; @@ -81,7 +75,7 @@ void ImProcFunctions::initCache () { cachef[i] = 327.68*( exp(1.0/3.0 * log((double)i / MAXVAL) )); } else { - cachef[i] = 327.68*((kappa*i/MAXVAL+16.0)/116.0); + cachef[i] = 327.68*((Color::kappa*i/MAXVAL+16.0)/116.0); } } @@ -212,17 +206,17 @@ void ImProcFunctions::rgbProc (Imagefloat* working, LabImage* lab, LUTf & hltone double toxyz[3][3] = { { - ( wprof[0][0] / D50x), - ( wprof[0][1] / D50x), - ( wprof[0][2] / D50x) + ( wprof[0][0] / Color::D50x), + ( wprof[0][1] / Color::D50x), + ( wprof[0][2] / Color::D50x) },{ ( wprof[1][0] ), ( wprof[1][1] ), ( wprof[1][2] ) },{ - ( wprof[2][0] / D50z), - ( wprof[2][1] / D50z), - ( wprof[2][2] / D50z) + ( wprof[2][0] / Color::D50z), + ( wprof[2][1] / Color::D50z), + ( wprof[2][2] / Color::D50z) } }; @@ -412,9 +406,9 @@ void ImProcFunctions::rgbProc (Imagefloat* working, LabImage* lab, LUTf & hltone float fx = (0.002 * lab->a[i][j])/327.68 + fy; float fz = fy - (0.005 * lab->b[i][j])/327.68; - float x_ = 65535*Lab2xyz(fx)*D50x; + float x_ = 65535*Lab2xyz(fx)*Color::D50x; float y_ = 65535*Lab2xyz(fy); - float z_ = 65535*Lab2xyz(fz)*D50z; + float z_ = 65535*Lab2xyz(fz)*Color::D50z; int R,G,B; xyz2srgb(x_,y_,z_,R,G,B); @@ -927,227 +921,5 @@ fclose(f);*/ } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - - void ImProcFunctions::rgb2hsv (float r, float g, float b, float &h, float &s, float &v) { - - double var_R = r / 65535.0; - double var_G = g / 65535.0; - double var_B = b / 65535.0; - - double var_Min = min(var_R,var_G,var_B); - double var_Max = max(var_R,var_G,var_B); - double del_Max = var_Max - var_Min; - v = var_Max; - if (del_Max<0.00001 && del_Max>-0.00001) { // no fabs, slow! - h = 0; - s = 0; - } - else { - s = del_Max/var_Max; - - if ( var_R == var_Max ) h = (var_G - var_B)/del_Max; - else if ( var_G == var_Max ) h = 2.0 + (var_B - var_R)/del_Max; - else if ( var_B == var_Max ) h = 4.0 + (var_R - var_G)/del_Max; - h /= 6.0; - - if ( h < 0 ) h += 1; - if ( h > 1 ) h -= 1; - } - } - - //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - // - // WARNING: rtengine/color.cc have similar functions. I guess that color.cc was mistakenly committed from the denoise branch - // in changeset #2493dd00b9f3. Code cleanup about those methods will have to be made when the denoise branch will be merged - // - //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - - void ImProcFunctions::hsv2rgb (float h, float s, float v, float &r, float &g, float &b) { - - float h1 = h*6; // sector 0 to 5 - int i = (int)h1; // floor() is very slow, and h1 is always >0 - float f = h1 - i; // fractional part of h - - float p = v * ( 1 - s ); - float q = v * ( 1 - s * f ); - float t = v * ( 1 - s * ( 1 - f ) ); - - float r1,g1,b1; - if (i==1) {r1 = q; g1 = v; b1 = p;} - else if (i==2) {r1 = p; g1 = v; b1 = t;} - else if (i==3) {r1 = p; g1 = q; b1 = v;} - else if (i==4) {r1 = t; g1 = p; b1 = v;} - else if (i==5) {r1 = v; g1 = p; b1 = q;} - else /*i==(0|6)*/ {r1 = v; g1 = t; b1 = p;} - - r = ((r1)*65535.0); - g = ((g1)*65535.0); - b = ((b1)*65535.0); - } - - // Function copied for speed concerns - // Not exactly the same as above ; this one return a result in the [0.0 ; 1.0] range - void ImProcFunctions::hsv2rgb01 (float h, float s, float v, float &r, float &g, float &b) { - - float h1 = h*6; // sector 0 to 5 - int i = int(h1); - float f = h1 - i; // fractional part of h - - float p = v * ( 1 - s ); - float q = v * ( 1 - s * f ); - float t = v * ( 1 - s * ( 1 - f ) ); - - if (i==1) {r = q; g = v; b = p;} - else if (i==2) {r = p; g = v; b = t;} - else if (i==3) {r = p; g = q; b = v;} - else if (i==4) {r = t; g = p; b = v;} - else if (i==5) {r = v; g = p; b = q;} - else /*(i==0|6)*/ {r = v; g = t; b = p;} - } - - void ImProcFunctions::xyz2srgb (float x, float y, float z, float &r, float &g, float &b) { - - //Transform to output color. Standard sRGB is D65, but internal representation is D50 - //Note that it is only at this point that we should have need of clipping color data - - /*float x65 = d65_d50[0][0]*x + d65_d50[0][1]*y + d65_d50[0][2]*z ; - float y65 = d65_d50[1][0]*x + d65_d50[1][1]*y + d65_d50[1][2]*z ; - float z65 = d65_d50[2][0]*x + d65_d50[2][1]*y + d65_d50[2][2]*z ; - - r = sRGB_xyz[0][0]*x65 + sRGB_xyz[0][1]*y65 + sRGB_xyz[0][2]*z65; - g = sRGB_xyz[1][0]*x65 + sRGB_xyz[1][1]*y65 + sRGB_xyz[1][2]*z65; - b = sRGB_xyz[2][0]*x65 + sRGB_xyz[2][1]*y65 + sRGB_xyz[2][2]*z65;*/ - - /*r = sRGBd65_xyz[0][0]*x + sRGBd65_xyz[0][1]*y + sRGBd65_xyz[0][2]*z ; - g = sRGBd65_xyz[1][0]*x + sRGBd65_xyz[1][1]*y + sRGBd65_xyz[1][2]*z ; - b = sRGBd65_xyz[2][0]*x + sRGBd65_xyz[2][1]*y + sRGBd65_xyz[2][2]*z ;*/ - - r = ((sRGB_xyz[0][0]*x + sRGB_xyz[0][1]*y + sRGB_xyz[0][2]*z)) ; - g = ((sRGB_xyz[1][0]*x + sRGB_xyz[1][1]*y + sRGB_xyz[1][2]*z)) ; - b = ((sRGB_xyz[2][0]*x + sRGB_xyz[2][1]*y + sRGB_xyz[2][2]*z)) ; - - } - - - void ImProcFunctions::xyz2rgb (float x, float y, float z, float &r, float &g, float &b, double rgb_xyz[3][3]) { - - //Transform to output color. Standard sRGB is D65, but internal representation is D50 - //Note that it is only at this point that we should have need of clipping color data - - /*float x65 = d65_d50[0][0]*x + d65_d50[0][1]*y + d65_d50[0][2]*z ; - float y65 = d65_d50[1][0]*x + d65_d50[1][1]*y + d65_d50[1][2]*z ; - float z65 = d65_d50[2][0]*x + d65_d50[2][1]*y + d65_d50[2][2]*z ; - - r = sRGB_xyz[0][0]*x65 + sRGB_xyz[0][1]*y65 + sRGB_xyz[0][2]*z65; - g = sRGB_xyz[1][0]*x65 + sRGB_xyz[1][1]*y65 + sRGB_xyz[1][2]*z65; - b = sRGB_xyz[2][0]*x65 + sRGB_xyz[2][1]*y65 + sRGB_xyz[2][2]*z65;*/ - - /*r = sRGBd65_xyz[0][0]*x + sRGBd65_xyz[0][1]*y + sRGBd65_xyz[0][2]*z ; - g = sRGBd65_xyz[1][0]*x + sRGBd65_xyz[1][1]*y + sRGBd65_xyz[1][2]*z ; - b = sRGBd65_xyz[2][0]*x + sRGBd65_xyz[2][1]*y + sRGBd65_xyz[2][2]*z ;*/ - - - - r = ((rgb_xyz[0][0]*x + rgb_xyz[0][1]*y + rgb_xyz[0][2]*z)) ; - g = ((rgb_xyz[1][0]*x + rgb_xyz[1][1]*y + rgb_xyz[1][2]*z)) ; - b = ((rgb_xyz[2][0]*x + rgb_xyz[2][1]*y + rgb_xyz[2][2]*z)) ; - - } - -void ImProcFunctions::calcGamma (double pwr, double ts, int mode, int imax, double &gamma0, double &gamma1, double &gamma2, double &gamma3, double &gamma4, double &gamma5) { -{//from Dcraw (D.Coffin) - int i; - double g[6], bnd[2]={0,0}; - - g[0] = pwr; - g[1] = ts; - g[2] = g[3] = g[4] = 0; - bnd[g[1] >= 1] = 1; - if (g[1] && (g[1]-1)*(g[0]-1) <= 0) { - for (i=0; i < 48; i++) { - g[2] = (bnd[0] + bnd[1])/2; - if (g[0]) bnd[(pow(g[2]/g[1],-g[0]) - 1)/g[0] - 1/g[2] > -1] = g[2]; - else bnd[g[2]/exp(1-1/g[2]) < g[1]] = g[2]; - } - g[3] = g[2] / g[1]; - if (g[0]) g[4] = g[2] * (1/g[0] - 1); - } - if (g[0]) g[5] = 1 / (g[1]*SQR(g[3])/2 - g[4]*(1 - g[3]) + - (1 - pow(g[3],1+g[0]))*(1 + g[4])/(1 + g[0])) - 1; - else g[5] = 1 / (g[1]*SQR(g[3])/2 + 1 - - g[2] - g[3] - g[2]*g[3]*(log(g[3]) - 1)) - 1; - if (!mode--) { - gamma0=g[0];gamma1=g[1];gamma2=g[2];gamma3=g[3];gamma4=g[4];gamma5=g[5]; - return; - } } -} - - void ImProcFunctions::Lab2XYZ(float L, float a, float b, float &x, float &y, float &z) { - float fy = (0.00862069 * L) + 0.137932; // (L+16)/116 - float fx = (0.002 * a) + fy; - float fz = fy - (0.005 * b); - - x = 65535*f2xyz(fx)*D50x; - y = 65535*f2xyz(fy); - z = 65535*f2xyz(fz)*D50z; - } - - void ImProcFunctions::XYZ2Lab(float X, float Y, float Z, float &L, float &a, float &b) { - - float fx = (X<65535.0 ? cachef[X] : (327.68*exp(log(X/MAXVAL)/3.0 ))); - float fy = (Y<65535.0 ? cachef[Y] : (327.68*exp(log(Y/MAXVAL)/3.0 ))); - float fz = (Z<65535.0 ? cachef[Z] : (327.68*exp(log(Z/MAXVAL)/3.0 ))); - - L = (116.0 * fy - 5242.88); //5242.88=16.0*327.68; - a = (500.0 * (fx - fy) ); - b = (200.0 * (fy - fz) ); - - } - - void ImProcFunctions::Lab2Yuv(float L, float a, float b, float &Y, float &u, float &v) { - float fy = (0.00862069 * L/327.68) + 0.137932; // (L+16)/116 - float fx = (0.002 * a/327.68) + fy; - float fz = fy - (0.005 * b/327.68); - - float X = 65535.0*f2xyz(fx)*D50x; - Y = 65535.0*f2xyz(fy); - float Z = 65535.0*f2xyz(fz)*D50z; - - u = 4.0*X/(X+15*Y+3*Z)-u0; - v = 9.0*Y/(X+15*Y+3*Z)-v0; - - /*float u0 = 4*D50x/(D50x+15+3*D50z); - float v0 = 9/(D50x+15+3*D50z); - u -= u0; - v -= v0;*/ - - } - - void ImProcFunctions::Yuv2Lab(float Yin, float u, float v, float &L, float &a, float &b, double wp[3][3]) { - - float u1 = u + u0; - float v1 = v + v0; - - float Y = Yin; - float X = (9*u1*Y)/(4*v1*D50x); - float Z = (12 - 3*u1 - 20*v1)*Y/(4*v1*D50z); - - Color::gamutmap(X,Y,Z,wp); - - float fx = (X<65535.0 ? cachef[X] : (327.68*exp(log(X/MAXVAL)/3.0 ))); - float fy = (Y<65535.0 ? cachef[Y] : (327.68*exp(log(Y/MAXVAL)/3.0 ))); - float fz = (Z<65535.0 ? cachef[Z] : (327.68*exp(log(Z/MAXVAL)/3.0 ))); - - L = (116.0 * fy - 5242.88); //5242.88=16.0*327.68; - a = (500.0 * (fx - fy) ); - b = (200.0 * (fy - fz) ); - - } - -#include "gamutbdy.cc" -} - -#undef eps_max -#undef kappa diff --git a/rtengine/improcfun.h b/rtengine/improcfun.h index 2fe8a7a76..e62663eef 100644 --- a/rtengine/improcfun.h +++ b/rtengine/improcfun.h @@ -187,27 +187,6 @@ namespace rtengine { void getAutoExp (LUTu & histogram, int histcompr, double defgain, double clip, double& expcomp, int& bright, int& contr, int& black, int& hlcompr, int& hlcomprthresh); static double getAutoDistor (const Glib::ustring& fname, int thumb_size); double getTransformAutoFill (int oW, int oH, const LCPMapper *pLCPMap=NULL); - - static void rgb2hsv (float r, float g, float b, float &h, float &s, float &v); - static void hsv2rgb (float h, float s, float v, float &r, float &g, float &b); - static void hsv2rgb01 (float h, float s, float v, float &r, float &g, float &b); - void xyz2srgb (float x, float y, float z, float &r, float &g, float &b); - void xyz2rgb (float x, float y, float z, float &r, float &g, float &b, double rgb_xyz[3][3]); - void Lab2XYZ(float L, float a, float b, float &x, float &y, float &z); - void XYZ2Lab(float X, float Y, float Z, float &L, float &a, float &b); - void Lab2Yuv(float L, float a, float b, float &Y, float &u, float &v); - void Yuv2Lab(float Y, float u, float v, float &L, float &a, float &b, double wp[3][3]); - void calcGamma (double pwr, double ts, int mode, int imax, double &gamma0, double &gamma1, double &gamma2, double &gamma3, double &gamma4,double &gamma5); - - //void gamutmap(LabImage* ); - void gamutmap(float &X, float &Y, float &Z, const double p[3][3]); - - static inline float f2xyz(float f) { - const float epsilonExpInv3 = 6.0/29.0; - const float kappaInv = 27.0/24389.0; // inverse of kappa - - return (f > epsilonExpInv3) ? f*f*f : (116 * f - 16) * kappaInv; - } }; } #endif diff --git a/rtengine/iplab2rgb.cc b/rtengine/iplab2rgb.cc index fbdad662b..4e62cfe88 100644 --- a/rtengine/iplab2rgb.cc +++ b/rtengine/iplab2rgb.cc @@ -74,9 +74,9 @@ void ImProcFunctions::lab2monitorRgb (LabImage* lab, Image8* image) { fx = (0.002 * ra[j]) / 327.68 + fy; fz = fy - (0.005 * rb[j]) / 327.68; - x_ = f2xyz(fx)*D50x; + x_ = f2xyz(fx)*Color::D50x; y_ = f2xyz(fy); - z_ = f2xyz(fz)*D50z; + z_ = f2xyz(fz)*Color::D50z; buffer[iy++] = (unsigned short)CLIP(x_* MAXVAL+0.5); buffer[iy++] = (unsigned short)CLIP(y_* MAXVAL+0.5); @@ -108,9 +108,9 @@ void ImProcFunctions::lab2monitorRgb (LabImage* lab, Image8* image) { fx = (0.002 * ra[j]) / 327.68 + fy; fz = fy - (0.005 * rb[j]) / 327.68; - x_ = 65535.0 * f2xyz(fx)*D50x; + x_ = 65535.0 * f2xyz(fx)*Color::D50x; y_ = 65535.0 * f2xyz(fy); - z_ = 65535.0 * f2xyz(fz)*D50z; + z_ = 65535.0 * f2xyz(fz)*Color::D50z; xyz2srgb(x_,y_,z_,R,G,B); @@ -165,9 +165,9 @@ Image8* ImProcFunctions::lab2rgb (LabImage* lab, int cx, int cy, int cw, int ch, float fx = (0.002 * ra[j])/327.68 + fy; float fz = fy - (0.005 * rb[j])/327.68; - float x_ = 65535.0 * f2xyz(fx)*D50x; + float x_ = 65535.0 * f2xyz(fx)*Color::D50x; float y_ = 65535.0 * f2xyz(fy); - float z_ = 65535.0 * f2xyz(fz)*D50z; + float z_ = 65535.0 * f2xyz(fz)*Color::D50z; buffer[iy++] = CLIP((int)(x_+0.5)); buffer[iy++] = CLIP((int)(y_+0.5)); @@ -205,11 +205,11 @@ Image8* ImProcFunctions::lab2rgb (LabImage* lab, int cx, int cy, int cw, int ch, float fx = (0.002 * ra[j])/327.68 + fy; float fz = fy - (0.005 * rb[j])/327.68; - float x_ = 65535.0 * f2xyz(fx)*D50x; + float x_ = 65535.0 * f2xyz(fx)*Color::D50x; float y_ = 65535.0 * f2xyz(fy); - float z_ = 65535.0 * f2xyz(fz)*D50z; + float z_ = 65535.0 * f2xyz(fz)*Color::D50z; - xyz2rgb(x_,y_,z_,R,G,B,rgb_xyz); + Color::xyz2rgb(x_,y_,z_,R,G,B,rgb_xyz); image->data[ix++] = (int)gamma2curve[CLIP(R)] >> 8; image->data[ix++] = (int)gamma2curve[CLIP(G)] >> 8; @@ -249,9 +249,9 @@ Image16* ImProcFunctions::lab2rgb16 (LabImage* lab, int cx, int cy, int cw, int float fx = (0.002 * ra[j])/327.68 + fy; float fz = fy - (0.005 * rb[j])/327.68; - float x_ = 65535.0 * f2xyz(fx)*D50x; + float x_ = 65535.0 * f2xyz(fx)*Color::D50x; float y_ = 65535.0 * f2xyz(fy); - float z_ = 65535.0 * f2xyz(fz)*D50z; + float z_ = 65535.0 * f2xyz(fz)*Color::D50z; xa[j-cx] = CLIP((int)(x_+0.5)); ya[j-cx] = CLIP((int)(y_+0.5)); @@ -280,9 +280,9 @@ Image16* ImProcFunctions::lab2rgb16 (LabImage* lab, int cx, int cy, int cw, int float fx = (0.002 * ra[j])/327.68 + fy; float fz = fy - (0.005 * rb[j])/327.68; - float x_ = 65535.0 * f2xyz(fx)*D50x; + float x_ = 65535.0 * f2xyz(fx)*Color::D50x; float y_ = 65535.0 * f2xyz(fy); - float z_ = 65535.0 * f2xyz(fz)*D50z; + float z_ = 65535.0 * f2xyz(fz)*Color::D50z; xyz2srgb(x_,y_,z_,R,G,B); @@ -295,7 +295,6 @@ Image16* ImProcFunctions::lab2rgb16 (LabImage* lab, int cx, int cy, int cw, int return image; } - // for gamma options (BT709...sRGB linear...) Image16* ImProcFunctions::lab2rgb16b (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile, Glib::ustring profi, Glib::ustring gam, bool freegamma, double gampos, double slpos, double &ga0, double &ga1, double &ga2, double &ga3, double &ga4, double &ga5, double &ga6) { @@ -342,7 +341,7 @@ Image16* ImProcFunctions::lab2rgb16b (LabImage* lab, int cx, int cy, int cw, int else //free gamma selected { if(slpos==0) slpos=eps; - calcGamma(pwr, ts, mode, imax,g_a0,g_a1,g_a2,g_a3,g_a4,g_a5);// call to calcGamma with selected gamma and slope : return parameters for LCMS2 + Color::calcGamma(pwr, ts, mode, imax,g_a0,g_a1,g_a2,g_a3,g_a4,g_a5);// call to calcGamma with selected gamma and slope : return parameters for LCMS2 ga4=g_a3*ts; //printf("g_a0=%f g_a1=%f g_a2=%f g_a3=%f g_a4=%f\n", g_a0,g_a1,g_a2,g_a3,g_a4); ga0=gampos;ga1=1./(1.0+g_a4);ga2=g_a4/(1.0 + g_a4);ga3=1./slpos;ga5=0.0; diff --git a/rtengine/ipvibrance.cc b/rtengine/ipvibrance.cc index c3e773884..f830d785c 100644 --- a/rtengine/ipvibrance.cc +++ b/rtengine/ipvibrance.cc @@ -2145,10 +2145,10 @@ void ImProcFunctions::vibrance (LabImage* lab) { fx = (0.002 * aprov1) + fy; fz = fy - (0.005 * bprov1); - x_ = 65535.0 * f2xyz(fx)*D50x; + x_ = 65535.0 * f2xyz(fx)*Color::D50x; y_ = 65535.0 * f2xyz(fy); - z_ = 65535.0 * f2xyz(fz)*D50z; - xyz2rgb(x_,y_,z_,R,G,B,wip); + z_ = 65535.0 * f2xyz(fz)*Color::D50z; + Color::xyz2rgb(x_,y_,z_,R,G,B,wip); // gamut control before saturation to put Lab values in future gamut, but not RGB if (allwaysingamut) { // gamut control @@ -2323,10 +2323,10 @@ void ImProcFunctions::vibrance (LabImage* lab) { fxx = (0.002 * aprovn) + fyy; fzz = fyy - (0.005 * bprovn); - xx_ = 65535.0 * f2xyz(fxx)*D50x; + xx_ = 65535.0 * f2xyz(fxx)*Color::D50x; yy_ = 65535.0 * f2xyz(fyy); - zz_ = 65535.0 * f2xyz(fzz)*D50z; - xyz2rgb(xx_,yy_,zz_,RR,GG,BB,wip); + zz_ = 65535.0 * f2xyz(fzz)*Color::D50z; + Color::xyz2rgb(xx_,yy_,zz_,RR,GG,BB,wip); if(RR<0.0 || GG < 0.0 || BB < 0.0) { #ifdef _DEBUG diff --git a/rtengine/rawimagesource.cc b/rtengine/rawimagesource.cc index 0a74a662f..8a9c7b869 100644 --- a/rtengine/rawimagesource.cc +++ b/rtengine/rawimagesource.cc @@ -2210,8 +2210,8 @@ void RawImageSource::HLRecovery_CIELab (float* rin, float* gin, float* bin, floa double fz = fy - y + z; double fx = fy + x - y; - double zr = ImProcFunctions::f2xyz(fz); - double xr = ImProcFunctions::f2xyz(fx); + double zr = Color::f2xyz(fz); + double xr = Color::f2xyz(fx); x = xr*65535.0 ; y = yy; diff --git a/rtengine/rt_math.h b/rtengine/rt_math.h index ae63da9ba..f5adb525d 100644 --- a/rtengine/rt_math.h +++ b/rtengine/rt_math.h @@ -5,9 +5,6 @@ #include "rtengine.h" -#define D50x 0.96422 -#define D50z 0.82521 - namespace rtengine { static const int MAXVAL = 0xffff; diff --git a/rtgui/hsvequalizer.cc b/rtgui/hsvequalizer.cc index f5dfce6ea..b5e77f1b1 100644 --- a/rtgui/hsvequalizer.cc +++ b/rtgui/hsvequalizer.cc @@ -18,7 +18,6 @@ */ #include "hsvequalizer.h" -#include "../rtengine/improcfun.h" #include "../rtengine/color.h" using namespace rtengine; @@ -196,19 +195,19 @@ void HSVEqualizer::colorForValue (double valX, double valY) { h -= 1.0f; else if (h < 0.0f) h += 1.0f; - ImProcFunctions::hsv2rgb01(h, 0.5f, 0.5f, r, g, b); + Color::hsv2rgb01(h, 0.5f, 0.5f, r, g, b); red = double(r); green = double(g); blue = double(b); } else if (ce == sshape) { // Saturation = f(Hue) - ImProcFunctions::hsv2rgb01(float(valX), float(valY), 0.5f, r, g, b); + Color::hsv2rgb01(float(valX), float(valY), 0.5f, r, g, b); red = double(r); green = double(g); blue = double(b); } else if (ce == vshape) { // Value = f(Hue) - ImProcFunctions::hsv2rgb01(float(valX), 0.5f, float(valY), r, g, b); + Color::hsv2rgb01(float(valX), 0.5f, float(valY), r, g, b); red = double(r); green = double(g); blue = double(b); diff --git a/rtgui/navigator.cc b/rtgui/navigator.cc index 7864b8deb..c5aec09bb 100644 --- a/rtgui/navigator.cc +++ b/rtgui/navigator.cc @@ -24,6 +24,8 @@ #include "../rtengine/color.h" #include "../rtengine/rt_math.h" +using namespace rtengine; + Navigator::Navigator () { set_label (M("MAIN_MSG_NAVIGATOR")); @@ -90,11 +92,11 @@ void Navigator::pointerMoved (bool validPos, Glib::ustring profile, int x, int y R->set_text (Glib::ustring::compose (M("NAVIGATOR_R_VALUE"), r)); G->set_text (Glib::ustring::compose (M("NAVIGATOR_G_VALUE"), g)); B->set_text (Glib::ustring::compose (M("NAVIGATOR_B_VALUE"), b)); - int h, s, v; - rgb2hsv (r, g, b, h, s, v); - H->set_text (Glib::ustring::compose (M("NAVIGATOR_H_VALUE"), h)); - S->set_text (Glib::ustring::compose (M("NAVIGATOR_S_VALUE"), s)); - V->set_text (Glib::ustring::compose (M("NAVIGATOR_V_VALUE"), v)); + float h, s, v; + Color::rgb2hsv (r*0xffff/0xff, g*0xffff/0xff, b*0xffff/0xff, h, s, v); + H->set_text (Glib::ustring::compose (M("NAVIGATOR_H_VALUE"), (int)(h*(float)0xff))); + S->set_text (Glib::ustring::compose (M("NAVIGATOR_S_VALUE"), (int)(s*(float)0xff))); + V->set_text (Glib::ustring::compose (M("NAVIGATOR_V_VALUE"), (int)(v*(float)0xff))); int LAB_a, LAB_b, LAB_l; //rgb2lab (r, g, b, LAB_l, LAB_a, LAB_b); rgb2lab (profile, r, g, b, LAB_l, LAB_a, LAB_b); @@ -106,67 +108,12 @@ void Navigator::pointerMoved (bool validPos, Glib::ustring profile, int x, int y } } -void Navigator::rgb2hsv (int r, int g, int b, int &h, int &s, int &v) { - - volatile double var_R = r / 255.0; - volatile double var_G = g / 255.0; - volatile double var_B = b / 255.0; - - /*volatile double var_Min = min(var_R,var_G,var_B); - volatile double var_Max = max(var_R,var_G,var_B); - double del_Max = var_Max - var_Min; - double V = (var_Max + var_Min) / 2; - double H, S; - if (fabs(del_Max)<0.0000001) { - H = 0; - S = 0; - } - else { - if (V < 0.5) S = del_Max / (var_Max + var_Min); - else S = del_Max / (2 - var_Max - var_Min); - - double del_R = ( ( ( var_Max - var_R ) / 6.0 ) + ( del_Max / 2.0 ) ) / del_Max; - double del_G = ( ( ( var_Max - var_G ) / 6.0 ) + ( del_Max / 2.0 ) ) / del_Max; - double del_B = ( ( ( var_Max - var_B ) / 6.0 ) + ( del_Max / 2.0 ) ) / del_Max; - if ( var_R == var_Max ) H = del_B - del_G; - else if ( var_G == var_Max ) H = (1.0 / 3.0) + del_R - del_B; - else if ( var_B == var_Max ) H = (2.0 / 3.0) + del_G - del_R; - - if ( H < 0 ) H += 1; - if ( H > 1 ) H -= 1; - }*/ - - double var_Min = rtengine::min(var_R,var_G,var_B); - double var_Max = rtengine::max(var_R,var_G,var_B); - double del_Max = var_Max - var_Min; - double V = var_Max; - double H, S; - if (fabs(del_Max)<0.001) { - H = 0; - S = 0; - } - else { - S = del_Max/var_Max; - - if ( var_R == var_Max ) H = (var_G - var_B)/del_Max; - else if ( var_G == var_Max ) H = 2.0 + (var_B - var_R)/del_Max; - else if ( var_B == var_Max ) H = 4.0 + (var_R - var_G)/del_Max; - H /= 6.0; - - if ( H < 0 ) H += 1; - if ( H > 1 ) H -= 1; - } - - h = (int)(H*255.0); - s = (int)(S*255.0); - v = (int)(V*255.0); -} void Navigator::rgb2lab (Glib::ustring profile, int r, int g, int b, int &LAB_l, int &LAB_a, int &LAB_b) { double xyz_rgb[3][3]; - double ep=216.0/24389.0; - double ka=24389.0/27.0; + const double ep=216.0/24389.0; + const double ka=24389.0/27.0; double var_R = r / 255.0; double var_G = g / 255.0; @@ -222,9 +169,9 @@ void Navigator::rgb2lab (Glib::ustring profile, int r, int g, int b, int &LAB_l, } double varxx,varyy,varzz; - double var_X = ( xyz_rgb[0][0]*var_R + xyz_rgb[0][1]*var_G + xyz_rgb[0][2]*var_B ) / D50x; + double var_X = ( xyz_rgb[0][0]*var_R + xyz_rgb[0][1]*var_G + xyz_rgb[0][2]*var_B ) / Color::D50x; double var_Y = ( xyz_rgb[1][0]*var_R + xyz_rgb[1][1]*var_G + xyz_rgb[1][2]*var_B ) ; - double var_Z = ( xyz_rgb[2][0]*var_R + xyz_rgb[2][1]*var_G + xyz_rgb[2][2]*var_B ) / D50z; + double var_Z = ( xyz_rgb[2][0]*var_R + xyz_rgb[2][1]*var_G + xyz_rgb[2][2]*var_B ) / Color::D50z; varxx = var_X>ep?pow (var_X, 1.0/3.0):( ka * var_X + 16.0) / 116.0 ; varyy = var_Y>ep?pow (var_Y, 1.0/3.0):( ka * var_Y + 16.0) / 116.0 ; diff --git a/rtgui/navigator.h b/rtgui/navigator.h index 3d8a0c35c..207a9daf3 100644 --- a/rtgui/navigator.h +++ b/rtgui/navigator.h @@ -34,8 +34,6 @@ class Navigator : public Gtk::Frame, public PointerMotionListener { Gtk::Label *H, *S, *V; Gtk::Label *LAB_A, *LAB_B, *LAB_L; - void rgb2hsv (int r, int g, int b, int &h, int &s, int &v); - //void rgb2lab (int r, int g, int b, int &LAB_l, int &LAB_a, int &LAB_b); void rgb2lab (Glib::ustring profile, int r, int g, int b, int &LAB_l, int &LAB_a, int &LAB_b); void setInvalid (int fullWidth=-1, int fullHeight=-1);