diff --git a/rtengine/CMakeLists.txt b/rtengine/CMakeLists.txt index 6eff1db5a..0cda358e0 100644 --- a/rtengine/CMakeLists.txt +++ b/rtengine/CMakeLists.txt @@ -23,6 +23,7 @@ set (RTENGINESOURCEFILES safegtk.cc colortemp.cc curves.cc flatcurves.cc diagona klt/convolve.cc klt/error.cc klt/klt.cc klt/klt_util.cc klt/pnmio.cc klt/pyramid.cc klt/selectGoodFeatures.cc klt/storeFeatures.cc klt/trackFeatures.cc klt/writeFeatures.cc clutstore.cc + ciecam02.cc ) include_directories (BEFORE "${CMAKE_CURRENT_BINARY_DIR}") diff --git a/rtengine/LUT.h b/rtengine/LUT.h index 464c588fb..c60765cf9 100644 --- a/rtengine/LUT.h +++ b/rtengine/LUT.h @@ -239,7 +239,7 @@ public: #if defined( __SSE2__ ) && defined( __x86_64__ ) __m128 operator[](__m128 indexv ) const { - printf("don't use this operator. It's not ready for production"); +// printf("don't use this operator. It's not ready for production"); return _mm_setzero_ps(); // convert floats to ints diff --git a/rtengine/ciecam02.cc b/rtengine/ciecam02.cc new file mode 100644 index 000000000..43dd2bc3b --- /dev/null +++ b/rtengine/ciecam02.cc @@ -0,0 +1,868 @@ +/* + * 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 "ciecam02.h" +#include "rtengine.h" +#include "curves.h" +#include +#include +#include +#include "sleef.c" +#include "settings.h" + +#undef CLIPD +#define CLIPD(a) ((a)>0.0?((a)<1.0?(a):1.0):0.0) +#define MAXR(a,b) ((a) > (b) ? (a) : (b)) +#define pow_F(a,b) (xexpf(b*xlogf(a))) + +namespace rtengine +{ + +extern const Settings* settings; + +void Ciecam02::curvecolor(double satind, double satval, double &sres, double parsat) +{ + if (satind >=0.0) { + sres = (1.-(satind)/100.)*satval+(satind)/100.*(1.-SQR(SQR(1.-min(satval,1.0)))); + if (sres>parsat) { sres=parsat; } + if (sres<0.) { sres=0.; } + } else { + if (satind < -0.1) { + sres = satval*(1.+(satind)/100.); + } + } +} + +void Ciecam02::curvecolorfloat(float satind, float satval, float &sres, float parsat) +{ + if (satind > 0.f) { + if (satval >= 1.f) { // The calculation below goes wrong direction when satval > 1 + sres = satval; + } else { + sres = (1.f-(satind)/100.f)*satval+(satind)/100.f*(1.f-SQR(SQR(1.f-min(satval,1.0f)))); + } + if (sres>parsat) { + sres = max(parsat,satval); + } + } else if (satind < 0.f) { + sres = satval*(1.f+(satind)/100.f); + } else { // satind == 0 means we don't want to change the value at all + sres = satval; + } +} + +void Ciecam02::curveJ (double br, double contr, int db, LUTf & outCurve, LUTu & histogram ) +{ + LUTf dcurve(65536,0); + int skip=1; + + // check if brightness curve is needed + if (br>0.00001 || br<-0.00001) { + + std::vector brightcurvePoints; + brightcurvePoints.resize(9); + brightcurvePoints.at(0) = double(DCT_NURBS); + + brightcurvePoints.at(1) = 0.; // black point. Value in [0 ; 1] range + brightcurvePoints.at(2) = 0.; // black point. Value in [0 ; 1] range + + if (br>0) { + brightcurvePoints.at(3) = 0.1; // toe point + brightcurvePoints.at(4) = 0.1+br/150.0; //value at toe point + + brightcurvePoints.at(5) = 0.7; // shoulder point + brightcurvePoints.at(6) = min(1.0,0.7+br/300.0); //value at shoulder point + } else { + brightcurvePoints.at(3) = 0.1-br/150.0; // toe point + brightcurvePoints.at(4) = 0.1; // value at toe point + + brightcurvePoints.at(5) = min(1.0,0.7-br/300.0); // shoulder point + brightcurvePoints.at(6) = 0.7; // value at shoulder point + } + brightcurvePoints.at(7) = 1.; // white point + brightcurvePoints.at(8) = 1.; // value at white point + + DiagonalCurve* brightcurve = new DiagonalCurve (brightcurvePoints, CURVES_MIN_POLY_POINTS/skip); + + // Applying brightness curve + for (int i=0; i<32768; i++) { + + // change to [0,1] range + float val = (float)i / 32767.0; + + // apply brightness curve + val = brightcurve->getVal (val); + + // store result in a temporary array + dcurve[i] = CLIPD(val); + } + delete brightcurve; + } else { + // for (int i=0; i<32768; i++) { // L values range up to 32767, higher values are for highlight overflow + for (int i=0; i<(32768*db); i++) { // L values range up to 32767, higher values are for highlight overflow + + // set the identity curve in the temporary array + dcurve[i] = (float)i / (db*32768.0f); + } + } + + + if (contr>0.00001 || contr<-0.00001) { + + // compute mean luminance of the image with the curve applied + int sum = 0; + float avg = 0; + //float sqavg = 0; + for (int i=0; i<32768; i++) { + avg += dcurve[i] * histogram[i];//approximation for average : usage of L (lab) instead of J + sum += histogram[i]; + } + avg /= sum; + std::vector contrastcurvePoints; + contrastcurvePoints.resize(9); + contrastcurvePoints.at(0) = double(DCT_NURBS); + + contrastcurvePoints.at(1) = 0.; // black point. Value in [0 ; 1] range + contrastcurvePoints.at(2) = 0.; // black point. Value in [0 ; 1] range + + contrastcurvePoints.at(3) = avg-avg*(0.6-contr/250.0); // toe point + contrastcurvePoints.at(4) = avg-avg*(0.6+contr/250.0); // value at toe point + + contrastcurvePoints.at(5) = avg+(1-avg)*(0.6-contr/250.0); // shoulder point + contrastcurvePoints.at(6) = avg+(1-avg)*(0.6+contr/250.0); // value at shoulder point + + contrastcurvePoints.at(7) = 1.; // white point + contrastcurvePoints.at(8) = 1.; // value at white point + + DiagonalCurve* contrastcurve = new DiagonalCurve (contrastcurvePoints, CURVES_MIN_POLY_POINTS/skip); + + // apply contrast enhancement + for (int i=0; i<(32768*db); i++) { + dcurve[i] = contrastcurve->getVal (dcurve[i]); + } + + delete contrastcurve; + } + + // for (int i=0; i<32768; i++) outCurve[i] = 32768.0*dcurve[i]; + for (int i=0; i<(db*32768); i++) { outCurve[i] = db*32768.0*dcurve[i]; } +} + +void Ciecam02::curveJfloat (float br, float contr, int db, LUTf & outCurve, LUTu & histogram ) +{ + LUTf dcurve(65536,0); + int skip=1; + + // check if brightness curve is needed + if (br>0.00001f || br<-0.00001f) { + + std::vector brightcurvePoints; + brightcurvePoints.resize(9); + brightcurvePoints.at(0) = double(DCT_NURBS); + + brightcurvePoints.at(1) = 0.f; // black point. Value in [0 ; 1] range + brightcurvePoints.at(2) = 0.f; // black point. Value in [0 ; 1] range + + if (br>0) { + brightcurvePoints.at(3) = 0.1f; // toe point + brightcurvePoints.at(4) = 0.1f+br/150.0f; //value at toe point + + brightcurvePoints.at(5) = 0.7f; // shoulder point + brightcurvePoints.at(6) = min(1.0f,0.7f+br/300.0f); //value at shoulder point + } else { + brightcurvePoints.at(3) = 0.1f-br/150.0f; // toe point + brightcurvePoints.at(4) = 0.1f; // value at toe point + + brightcurvePoints.at(5) = min(1.0f,0.7f-br/300.0f); // shoulder point + brightcurvePoints.at(6) = 0.7f; // value at shoulder point + } + brightcurvePoints.at(7) = 1.f; // white point + brightcurvePoints.at(8) = 1.f; // value at white point + + DiagonalCurve* brightcurve = new DiagonalCurve (brightcurvePoints, CURVES_MIN_POLY_POINTS/skip); + + // Applying brightness curve + for (int i=0; i<32768; i++) { + + // change to [0,1] range + float val = (float)i / 32767.0f; + + // apply brightness curve + val = brightcurve->getVal (val); + + // store result in a temporary array + dcurve[i] = CLIPD(val); + } + delete brightcurve; + } else { + // for (int i=0; i<32768; i++) { // L values range up to 32767, higher values are for highlight overflow + for (int i=0; i<(32768*db); i++) { // L values range up to 32767, higher values are for highlight overflow + + // set the identity curve in the temporary array + dcurve[i] = (float)i / (db*32768.0f); + } + } + + + if (contr>0.00001f || contr<-0.00001f) { + + // compute mean luminance of the image with the curve applied + int sum = 0; + float avg = 0; + //float sqavg = 0; + for (int i=0; i<32768; i++) { + avg += dcurve[i] * histogram[i];//approximation for average : usage of L (lab) instead of J + sum += histogram[i]; + } + avg /= sum; + //printf("avg=%f\n",avg); + std::vector contrastcurvePoints; + contrastcurvePoints.resize(9); + contrastcurvePoints.at(0) = double(DCT_NURBS); + + contrastcurvePoints.at(1) = 0.f; // black point. Value in [0 ; 1] range + contrastcurvePoints.at(2) = 0.f; // black point. Value in [0 ; 1] range + + contrastcurvePoints.at(3) = avg-avg*(0.6f-contr/250.0f); // toe point + contrastcurvePoints.at(4) = avg-avg*(0.6f+contr/250.0f); // value at toe point + + contrastcurvePoints.at(5) = avg+(1-avg)*(0.6f-contr/250.0f); // shoulder point + contrastcurvePoints.at(6) = avg+(1-avg)*(0.6f+contr/250.0f); // value at shoulder point + + contrastcurvePoints.at(7) = 1.f; // white point + contrastcurvePoints.at(8) = 1.f; // value at white point + + DiagonalCurve* contrastcurve = new DiagonalCurve (contrastcurvePoints, CURVES_MIN_POLY_POINTS/skip); + + // apply contrast enhancement + for (int i=0; i<(32768*db); i++) { + dcurve[i] = contrastcurve->getVal (dcurve[i]); + } + + delete contrastcurve; + } + + // for (int i=0; i<32768; i++) outCurve[i] = 32768.0*dcurve[i]; + for (int i=0; i<(db*32768); i++) { outCurve[i] = db*32768.0f*dcurve[i]; } +} + +/** + * Copyright (c) 2003 Billy Biggs + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + * + */ + +double Ciecam02::d_factor( double f, double la ) +{ + return f * (1.0 - ((1.0 / 3.6) * exp((-la - 42.0) / 92.0))); +} + +float Ciecam02::d_factorfloat( float f, float la ) +{ + return f * (1.0f - ((1.0f / 3.6f) * xexpf((-la - 42.0f) / 92.0f))); +} + +double Ciecam02::calculate_fl_from_la_ciecam02( double la ) +{ + double la5 = la * 5.0; + double k = 1.0 / (la5 + 1.0); + + /* Calculate k^4. */ + k = k * k; + k = k * k; + + return (0.2 * k * la5) + (0.1 * (1.0 - k) * (1.0 - k) * pow(la5, 1.0 / 3.0)); +} + +float Ciecam02::calculate_fl_from_la_ciecam02float( float la ) +{ + float la5 = la * 5.0f; + float k = 1.0f / (la5 + 1.0f); + + /* Calculate k^4. */ + k = k * k; + k = k * k; + + return (0.2f * k * la5) + (0.1f * (1.0f - k) * (1.0f - k) * pow_F(la5, 1.0f / 3.0f)); +} + +double Ciecam02::achromatic_response_to_white( double x, double y, double z, double d, double fl, double nbb, int gamu ) +{ + double r, g, b; + double rc, gc, bc; + double rp, gp, bp; + double rpa, gpa, bpa; + gamu=1; + xyz_to_cat02( r, g, b, x, y, z, gamu ); + + rc = r * (((y * d) / r) + (1.0 - d)); + gc = g * (((y * d) / g) + (1.0 - d)); + bc = b * (((y * d) / b) + (1.0 - d)); + + cat02_to_hpe( rp, gp, bp, rc, gc, bc, gamu ); + if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + rp=MAXR(rp,0.0); + gp=MAXR(gp,0.0); + bp=MAXR(bp,0.0); + } + + rpa = nonlinear_adaptation( rp, fl ); + gpa = nonlinear_adaptation( gp, fl ); + bpa = nonlinear_adaptation( bp, fl ); + + return ((2.0 * rpa) + gpa + ((1.0 / 20.0) * bpa) - 0.305) * nbb; +} + +float Ciecam02::achromatic_response_to_whitefloat( float x, float y, float z, float d, float fl, float nbb, int gamu ) +{ + float r, g, b; + float rc, gc, bc; + float rp, gp, bp; + float rpa, gpa, bpa; + gamu=1; + xyz_to_cat02float( r, g, b, x, y, z, gamu ); + + rc = r * (((y * d) / r) + (1.0f - d)); + gc = g * (((y * d) / g) + (1.0f - d)); + bc = b * (((y * d) / b) + (1.0f - d)); + + cat02_to_hpefloat( rp, gp, bp, rc, gc, bc, gamu ); + if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + rp=MAXR(rp,0.0f); + gp=MAXR(gp,0.0f); + bp=MAXR(bp,0.0f); + } + + rpa = nonlinear_adaptationfloat( rp, fl ); + gpa = nonlinear_adaptationfloat( gp, fl ); + bpa = nonlinear_adaptationfloat( bp, fl ); + + return ((2.0f * rpa) + gpa + ((1.0f / 20.0f) * bpa) - 0.305f) * nbb; +} + +void Ciecam02::xyz_to_cat02( double &r, double &g, double &b, double x, double y, double z, int gamu ) +{ + gamu=1; + if (gamu==0) { + r = ( 0.7328 * x) + (0.4296 * y) - (0.1624 * z); + g = (-0.7036 * x) + (1.6975 * y) + (0.0061 * z); + b = ( 0.0030 * x) + (0.0136 * y) + (0.9834 * z); + } else if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + //r = ( 0.7328 * x) + (0.4296 * y) - (0.1624 * z); + //g = (-0.7036 * x) + (1.6975 * y) + (0.0061 * z); + //b = ( 0.0000 * x) + (0.0000 * y) + (1.0000 * z); + r = ( 1.007245 * x) + (0.011136* y) - (0.018381 * z);//Changjun Li + g = (-0.318061 * x) + (1.314589 * y) + (0.003471 * z); + b = ( 0.0000 * x) + (0.0000 * y) + (1.0000 * z); + } +} + +void Ciecam02::xyz_to_cat02float( float &r, float &g, float &b, float x, float y, float z, int gamu ) +{ + gamu=1; + if (gamu==0) { + r = ( 0.7328f * x) + (0.4296f * y) - (0.1624f * z); + g = (-0.7036f * x) + (1.6975f * y) + (0.0061f * z); + b = ( 0.0030f * x) + (0.0136f * y) + (0.9834f * z); + } else if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + //r = ( 0.7328 * x) + (0.4296 * y) - (0.1624 * z); + //g = (-0.7036 * x) + (1.6975 * y) + (0.0061 * z); + //b = ( 0.0000 * x) + (0.0000 * y) + (1.0000 * z); + r = ( 1.007245f * x) + (0.011136f* y) - (0.018381f * z);//Changjun Li + g = (-0.318061f * x) + (1.314589f * y) + (0.003471f * z); + b = ( 0.0000f * x) + (0.0000f * y) + (1.0000f * z); + } +} + +void Ciecam02::cat02_to_xyz( double &x, double &y, double &z, double r, double g, double b, int gamu ) +{ + gamu=1; + if (gamu==0) { + x = ( 1.096124 * r) - (0.278869 * g) + (0.182745 * b); + y = ( 0.454369 * r) + (0.473533 * g) + (0.072098 * b); + z = (-0.009628 * r) - (0.005698 * g) + (1.015326 * b); + } else if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + //x = ( 1.0978566 * r) - (0.277843 * g) + (0.179987 * b); + //y = ( 0.455053 * r) + (0.473938 * g) + (0.0710096* b); + //z = ( 0.000000 * r) - (0.000000 * g) + (1.000000 * b); + x = ( 0.99015849 * r) - (0.00838772* g) + (0.018229217 * b);//Changjun Li + y = ( 0.239565979 * r) + (0.758664642 * g) + (0.001770137* b); + z = ( 0.000000 * r) - (0.000000 * g) + (1.000000 * b); + } +} + +void Ciecam02::cat02_to_xyzfloat( float &x, float &y, float &z, float r, float g, float b, int gamu ) +{ + gamu=1; + if (gamu==0) { + x = ( 1.096124f * r) - (0.278869f * g) + (0.182745f * b); + y = ( 0.454369f * r) + (0.473533f * g) + (0.072098f * b); + z = (-0.009628f * r) - (0.005698f * g) + (1.015326f * b); + } else if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + //x = ( 1.0978566 * r) - (0.277843 * g) + (0.179987 * b); + //y = ( 0.455053 * r) + (0.473938 * g) + (0.0710096* b); + //z = ( 0.000000 * r) - (0.000000 * g) + (1.000000 * b); + x = ( 0.99015849f * r) - (0.00838772f* g) + (0.018229217f * b);//Changjun Li + y = ( 0.239565979f * r) + (0.758664642f * g) + (0.001770137f* b); + z = ( 0.000000f * r) - (0.000000f * g) + (1.000000f * b); + } +} + + +void Ciecam02::hpe_to_xyz( double &x, double &y, double &z, double r, double g, double b ) +{ + x = (1.910197 * r) - (1.112124 * g) + (0.201908 * b); + y = (0.370950 * r) + (0.629054 * g) - (0.000008 * b); + z = b; +} + +void Ciecam02::hpe_to_xyzfloat( float &x, float &y, float &z, float r, float g, float b ) +{ + x = (1.910197f * r) - (1.112124f * g) + (0.201908f * b); + y = (0.370950f * r) + (0.629054f * g) - (0.000008f * b); + z = b; +} + +void Ciecam02::cat02_to_hpe( double &rh, double &gh, double &bh, double r, double g, double b, int gamu ) +{ + gamu=1; + if (gamu==0) { + rh = ( 0.7409792 * r) + (0.2180250 * g) + (0.0410058 * b); + gh = ( 0.2853532 * r) + (0.6242014 * g) + (0.0904454 * b); + bh = (-0.0096280 * r) - (0.0056980 * g) + (1.0153260 * b); + } else if (gamu==1) { //Changjun Li + rh = ( 0.550930835 * r) + (0.519435987* g) - ( 0.070356303* b); + gh = ( 0.055954056 * r) + (0.89973132 * g) + (0.044315524 * b); + bh = (0.0 * r) - (0.0* g) + (1.0 * b); + } +} + +void Ciecam02::cat02_to_hpefloat( float &rh, float &gh, float &bh, float r, float g, float b, int gamu ) +{ + gamu=1; + if (gamu==0) { + rh = ( 0.7409792f * r) + (0.2180250f * g) + (0.0410058f * b); + gh = ( 0.2853532f * r) + (0.6242014f * g) + (0.0904454f * b); + bh = (-0.0096280f * r) - (0.0056980f * g) + (1.0153260f * b); + } else if (gamu==1) { //Changjun Li + rh = ( 0.550930835f * r) + (0.519435987f* g) - ( 0.070356303f* b); + gh = ( 0.055954056f * r) + (0.89973132f * g) + (0.044315524f * b); + bh = (0.0f * r) - (0.0f* g) + (1.0f * b); + } +} + +void Ciecam02::Aab_to_rgb( double &r, double &g, double &b, double A, double aa, double bb, double nbb ) +{ + double x = (A / nbb) + 0.305; + + /* c1 c2 c3 */ + r = (0.32787 * x) + (0.32145 * aa) + (0.20527 * bb); + /* c1 c4 c5 */ + g = (0.32787 * x) - (0.63507 * aa) - (0.18603 * bb); + /* c1 c6 c7 */ + b = (0.32787 * x) - (0.15681 * aa) - (4.49038 * bb); +} +void Ciecam02::Aab_to_rgbfloat( float &r, float &g, float &b, float A, float aa, float bb, float nbb ) +{ + float x = (A / nbb) + 0.305f; + + /* c1 c2 c3 */ + r = (0.32787f * x) + (0.32145f * aa) + (0.20527f * bb); + /* c1 c4 c5 */ + g = (0.32787f * x) - (0.63507f * aa) - (0.18603f * bb); + /* c1 c6 c7 */ + b = (0.32787f * x) - (0.15681f * aa) - (4.49038f * bb); +} + +void Ciecam02::calculate_ab( double &aa, double &bb, double h, double e, double t, double nbb, double a ) +{ + double hrad = (h * M_PI) / 180.0; + double sinh = sin( hrad ); + double cosh = cos( hrad ); + double x = (a / nbb) + 0.305; + double p3 = 21.0/20.0; + if ( fabs( sinh ) >= fabs( cosh ) ) { + bb = ((0.32787 * x) * (2.0 + p3)) / + ((e / (t * sinh)) - + // ((0.32145 - 0.63507 - (p3 * 0.15681)) * (cosh / sinh)) - + // (0.20527 - 0.18603 - (p3 * 4.49038))); + ((-0.31362 - (p3 * 0.15681)) * (cosh / sinh)) - + (0.01924 - (p3 * 4.49038))); + + aa = (bb * cosh) / sinh; + } else { + aa = ((0.32787 * x) * (2.0 + p3)) / + ((e / (t * cosh)) - + // (0.32145 - 0.63507 - (p3 * 0.15681)) - + // ((0.20527 - 0.18603 - (p3 * 4.49038)) * (sinh / cosh))); + (-0.31362 - (p3 * 0.15681)) - + ((0.01924 - (p3 * 4.49038)) * (sinh / cosh))); + + bb = (aa * sinh) / cosh; + } +} + +void Ciecam02::calculate_abfloat( float &aa, float &bb, float h, float e, float t, float nbb, float a ) +{ + float2 sincosval = xsincosf((h * M_PI) / 180.0f); + float sinh = sincosval.x; + float cosh = sincosval.y; + float x = (a / nbb) + 0.305f; + float p3 = 1.05f; + bool swapValues = fabs( sinh ) > fabs( cosh ); + if (swapValues) { + std::swap(sinh,cosh); + } + + float div = ((e / (t * cosh)) - (-0.31362f - (p3 * 0.15681f)) - ((0.01924f - (p3 * 4.49038f)) * (sinh / cosh))); + // for large values of t the above calculation can change its sign which results in a hue shift of 180 degree + // so we have to check the sign to avoid this shift. + // Additionally it seems useful to limit the minimum value of div + // I limited it, but I'm sure the actual limit is not the best one + + if (signf(div) != signf(cosh) || fabsf(div) <= fabsf(cosh) * 2.f) { + div = cosh * 2.f; + } + + aa = ((0.32787f * x) * (2.0f + p3)) / div; + bb = (aa * sinh) / cosh; + + if (swapValues) { + std::swap(aa,bb); + } +} + +void Ciecam02::initcam1(double gamu, double yb, double pilotd, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, + double &cz, double &aw, double &wh, double &pfl, double &fl, double &c) +{ + n = yb / yw; + if (pilotd==2.0) { d = d_factor( f, la ); } + else { d=pilotd; } + fl = calculate_fl_from_la_ciecam02( la ); + nbb = ncb = 0.725 * pow( 1.0 / n, 0.2 ); + cz = 1.48 + sqrt( n ); + aw = achromatic_response_to_white( xw, yw, zw, d, fl, nbb, gamu ); + wh =( 4.0 / c ) * ( aw + 4.0 ) * pow( fl, 0.25 ); + pfl = pow( fl, 0.25 ); +#ifdef _DEBUG + if (settings->verbose) { printf("Source double d=%f aw=%f fl=%f wh=%f\n",d,aw,fl,wh); } +#endif +} + +void Ciecam02::initcam1float(float gamu, float yb, float pilotd, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, + float &cz, float &aw, float &wh, float &pfl, float &fl, float &c) +{ + n = yb / yw; + if (pilotd==2.0) { d = d_factorfloat( f, la ); } + else { d=pilotd; } + fl = calculate_fl_from_la_ciecam02float( la ); + nbb = ncb = 0.725f * pow_F( 1.0f / n, 0.2f ); + cz = 1.48f + sqrt( n ); + aw = achromatic_response_to_whitefloat( xw, yw, zw, d, fl, nbb, gamu ); + wh =( 4.0f / c ) * ( aw + 4.0f ) * pow_F( fl, 0.25f ); + pfl = pow_F( fl, 0.25f ); +#ifdef _DEBUG + if (settings->verbose) { printf("Source float d=%f aw=%f fl=%f wh=%f c=%f awc=%f\n",d,aw,fl,wh,c,(4.f/c)*(aw+4.f)); } +#endif +} + +void Ciecam02::initcam2(double gamu, double yb, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, + double &cz, double &aw, double &fl) +{ + n = yb / yw; + d = d_factor( f, la ); + fl = calculate_fl_from_la_ciecam02( la ); + nbb = ncb = 0.725 * pow( 1.0 / n, 0.2 ); + cz = 1.48 + sqrt( n ); + aw = achromatic_response_to_white( xw, yw, zw, d, fl, nbb, gamu ); +#ifdef _DEBUG + if (settings->verbose) { printf("Viewing double d=%f aw=%f fl=%f n=%f\n",d,aw,fl,n); } +#endif +} + +void Ciecam02::initcam2float(float gamu, float yb, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, + float &cz, float &aw, float &fl) +{ + n = yb / yw; + d = d_factorfloat( f, la ); + fl = calculate_fl_from_la_ciecam02float( la ); + nbb = ncb = 0.725f * pow_F( 1.0f / n, 0.2f ); + cz = 1.48f + sqrt( n ); + aw = achromatic_response_to_whitefloat( xw, yw, zw, d, fl, nbb, gamu ); +#ifdef _DEBUG + if (settings->verbose) { printf("Viewing float d=%f aw=%f fl=%f n=%f\n",d,aw,fl,n); } +#endif +} + +void Ciecam02::xyz2jchqms_ciecam02( double &J, double &C, double &h, double &Q, double &M, double &s,double &aw, double &fl, double &wh, + double x, double y, double z, double xw, double yw, double zw, + double yb, double la, double f, double c, double nc, double pilotd, int gamu , double n, double nbb, double ncb, double pfl, double cz, double d) +{ + double r, g, b; + double rw, gw, bw; + double rc, gc, bc; + double rp, gp, bp; + double rpa, gpa, bpa; + double a, ca, cb; + double e, t; + double myh; + gamu=1; + xyz_to_cat02( r, g, b, x, y, z, gamu ); + xyz_to_cat02( rw, gw, bw, xw, yw, zw, gamu ); + rc = r * (((yw * d) / rw) + (1.0 - d)); + gc = g * (((yw * d) / gw) + (1.0 - d)); + bc = b * (((yw * d) / bw) + (1.0 - d)); + + cat02_to_hpe( rp, gp, bp, rc, gc, bc, gamu ); + if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + rp=MAXR(rp,0.0); + gp=MAXR(gp,0.0); + bp=MAXR(bp,0.0); + } + rpa = nonlinear_adaptation( rp, fl ); + gpa = nonlinear_adaptation( gp, fl ); + bpa = nonlinear_adaptation( bp, fl ); + + ca = rpa - ((12.0 * gpa) / 11.0) + (bpa / 11.0); + cb = (1.0 / 9.0) * (rpa + gpa - (2.0 * bpa)); + + myh = (180.0 / M_PI) * atan2( cb, ca ); + if ( myh < 0.0 ) { myh += 360.0; } + //we can also calculate H, if necessary...but it's using time...for what usage ? + /*double temp; + if(myh<20.14) { + temp = ((myh + 122.47)/1.2) + ((20.14 - myh)/0.8); + H = 300 + (100*((myh + 122.47)/1.2)) / temp; + } + else if(myh < 90.0) { + temp = ((myh - 20.14)/0.8) + ((90.00 - myh)/0.7); + H = (100*((myh - 20.14)/0.8)) / temp; + } + else if (myh < 164.25) { + temp = ((myh - 90.00)/0.7) + ((164.25 - myh)/1.0); + H = 100 + ((100*((myh - 90.00)/0.7)) / temp); + } + else if (myh < 237.53) { + temp = ((myh - 164.25)/1.0) + ((237.53 - myh)/1.2); + H = 200 + ((100*((myh - 164.25)/1.0)) / temp); + } + else { + temp = ((myh - 237.53)/1.2) + ((360 - myh + 20.14)/0.8); + H = 300 + ((100*((myh - 237.53)/1.2)) / temp); + } + */ + a = ((2.0 * rpa) + gpa + ((1.0 / 20.0) * bpa) - 0.305) * nbb; + if (gamu==1) { a=MAXR(a,0.0); }//gamut correction M.H.Brill S.Susstrunk + J = 100.0 * pow( a / aw, c * cz ); + + e = ((12500.0 / 13.0) * nc * ncb) * (cos( ((myh * M_PI) / 180.0) + 2.0 ) + 3.8); + t = (e * sqrt( (ca * ca) + (cb * cb) )) / (rpa + gpa + ((21.0 / 20.0) * bpa)); + + C = pow( t, 0.9 ) * sqrt( J / 100.0 ) + * pow( 1.64 - pow( 0.29, n ), 0.73 ); + + Q = wh * sqrt( J / 100.0 ); + M = C * pfl; + s = 100.0 * sqrt( M / Q ); + h = myh; +} + +void Ciecam02::xyz2jchqms_ciecam02float( float &J, float &C, float &h, float &Q, float &M, float &s,float &aw, float &fl, float &wh, + float x, float y, float z, float xw, float yw, float zw, + float yb, float la, float f, float c, float nc, float pilotd, int gamu, float pow1, float nbb, float ncb, float pfl, float cz, float d) + +{ + float r, g, b; + float rw, gw, bw; + float rc, gc, bc; + float rp, gp, bp; + float rpa, gpa, bpa; + float a, ca, cb; + float e, t; + float myh; + gamu=1; + xyz_to_cat02float( r, g, b, x, y, z, gamu ); + xyz_to_cat02float( rw, gw, bw, xw, yw, zw, gamu ); + rc = r * (((yw * d) / rw) + (1.0 - d)); + gc = g * (((yw * d) / gw) + (1.0 - d)); + bc = b * (((yw * d) / bw) + (1.0 - d)); + + cat02_to_hpefloat( rp, gp, bp, rc, gc, bc, gamu ); + if (gamu==1) { //gamut correction M.H.Brill S.Susstrunk + rp=MAXR(rp,0.0f); + gp=MAXR(gp,0.0f); + bp=MAXR(bp,0.0f); + } + rpa = nonlinear_adaptationfloat( rp, fl ); + gpa = nonlinear_adaptationfloat( gp, fl ); + bpa = nonlinear_adaptationfloat( bp, fl ); + + ca = rpa - ((12.0f * gpa) - bpa) / 11.0f; + cb = (0.11111111f) * (rpa + gpa - (2.0f * bpa)); + + myh = xatan2f( cb, ca ); + if ( myh < 0.0f ) { + myh += (2.f * M_PI); + } + + a = ((2.0f * rpa) + gpa + (0.05f * bpa) - 0.305f) * nbb; + if (gamu==1) { + a=MAXR(a,0.0f); //gamut correction M.H.Brill S.Susstrunk + } + + J = pow_F( a / aw, c * cz * 0.5f); + + e = ((961.53846f) * nc * ncb) * (xcosf( myh + 2.0f ) + 3.8f); + t = (e * sqrtf( (ca * ca) + (cb * cb) )) / (rpa + gpa + (1.05f * bpa)); + + C = pow_F( t, 0.9f ) * J * pow1; + + Q = wh * J; + J *= J * 100.0f; + M = C * pfl; + Q = (Q == 0.f ? 0.0001f : Q); // avoid division by zero + s = 100.0f * sqrtf( M / Q ); + h = (myh * 180.f) / (float)M_PI; +} + +void Ciecam02::jch2xyz_ciecam02( double &x, double &y, double &z, double J, double C, double h, + double xw, double yw, double zw, double yb, double la, + double f, double c, double nc , int gamu, double n, double nbb, double ncb, double fl, double cz, double d, double aw ) +{ + double r, g, b; + double rc, gc, bc; + double rp, gp, bp; + double rpa, gpa, bpa; + double rw, gw, bw; + double a, ca, cb; + double e, t; + gamu=1; + xyz_to_cat02( rw, gw, bw, xw, yw, zw, gamu ); + e = ((12500.0 / 13.0) * nc * ncb) * (cos( ((h * M_PI) / 180.0) + 2.0 ) + 3.8); + a = pow( J / 100.0, 1.0 / (c * cz) ) * aw; + t = pow( C / (sqrt( J / 100) * pow( 1.64 - pow( 0.29, n ), 0.73 )), 10.0 / 9.0 ); + + calculate_ab( ca, cb, h, e, t, nbb, a ); + Aab_to_rgb( rpa, gpa, bpa, a, ca, cb, nbb ); + + rp = inverse_nonlinear_adaptation( rpa, fl ); + gp = inverse_nonlinear_adaptation( gpa, fl ); + bp = inverse_nonlinear_adaptation( bpa, fl ); + + hpe_to_xyz( x, y, z, rp, gp, bp ); + xyz_to_cat02( rc, gc, bc, x, y, z, gamu ); + + r = rc / (((yw * d) / rw) + (1.0 - d)); + g = gc / (((yw * d) / gw) + (1.0 - d)); + b = bc / (((yw * d) / bw) + (1.0 - d)); + + cat02_to_xyz( x, y, z, r, g, b, gamu ); +} + +void Ciecam02::jch2xyz_ciecam02float( float &x, float &y, float &z, float J, float C, float h, + float xw, float yw, float zw, float yb, float la, + float f, float c, float nc , int gamu, float pow1, float nbb, float ncb, float fl, float cz, float d, float aw) +{ + float r, g, b; + float rc, gc, bc; + float rp, gp, bp; + float rpa, gpa, bpa; + float rw, gw, bw; + float a, ca, cb; + float e, t; + gamu=1; + xyz_to_cat02float( rw, gw, bw, xw, yw, zw, gamu ); + e = ((961.53846f) * nc * ncb) * (xcosf( ((h * M_PI) / 180.0f) + 2.0f ) + 3.8f); + a = pow_F( J / 100.0f, 1.0f / (c * cz) ) * aw; + t = pow_F( 10.f * C / (sqrtf( J ) * pow1), 1.1111111f ); + + calculate_abfloat( ca, cb, h, e, t, nbb, a ); + Aab_to_rgbfloat( rpa, gpa, bpa, a, ca, cb, nbb ); + + rp = inverse_nonlinear_adaptationfloat( rpa, fl ); + gp = inverse_nonlinear_adaptationfloat( gpa, fl ); + bp = inverse_nonlinear_adaptationfloat( bpa, fl ); + + hpe_to_xyzfloat( x, y, z, rp, gp, bp ); + xyz_to_cat02float( rc, gc, bc, x, y, z, gamu ); + + r = rc / (((yw * d) / rw) + (1.0f - d)); + g = gc / (((yw * d) / gw) + (1.0f - d)); + b = bc / (((yw * d) / bw) + (1.0f - d)); + + cat02_to_xyzfloat( x, y, z, r, g, b, gamu ); +} + +double Ciecam02::nonlinear_adaptation( double c, double fl ) +{ + double p; + if (c<0.0) { p = pow( (-1.0*fl * c) / 100.0, 0.42 ); return ((-1.0*400.0 * p) / (27.13 + p)) + 0.1;} + else {p = pow( (fl * c) / 100.0, 0.42 ); return ((400.0 * p) / (27.13 + p)) + 0.1;} +} + +float Ciecam02::nonlinear_adaptationfloat( float c, float fl ) +{ + float p; + if (c<0.0f) { p = pow_F( (-1.0f*fl * c) / 100.0f, 0.42f ); return ((-1.0f*400.0f * p) / (27.13f + p)) + 0.1f;} + else {p = pow_F( (fl * c) / 100.0f, 0.42f ); return ((400.0f * p) / (27.13f + p)) + 0.1f;} +} + +double Ciecam02::inverse_nonlinear_adaptation( double c, double fl ) +{ + int c1; + if (c-0.1 < 0.0) { c1=-1; } + else { c1=1; } + return c1*(100.0 / fl) * pow( (27.13 * fabs( c - 0.1 )) / (400.0 - fabs( c - 0.1 )), 1.0 / 0.42 ); +} + +float Ciecam02::inverse_nonlinear_adaptationfloat( float c, float fl ) +{ + c -= 0.1f; + if (c < 0.f) { + fl *= -1.f; + if (c < -399.99f) { // avoid nan values + c = -399.99f; + } + } else if (c > 399.99f) { // avoid nan values + c = 399.99f; + } + return (100.0f / fl) * pow_F( (27.13f * fabsf( c )) / (400.0f - fabsf( c )), 2.38095238f ); +} + +//end CIECAM Billy Bigg + +} diff --git a/rtengine/ciecam02.h b/rtengine/ciecam02.h new file mode 100644 index 000000000..a078a1136 --- /dev/null +++ b/rtengine/ciecam02.h @@ -0,0 +1,110 @@ +/* + * 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 . + */ +#ifndef _CIECAM02_ +#define _CIECAM02_ +#include +#include "LUT.h" + +namespace rtengine +{ + +class Ciecam02 +{ +private: + static double d_factor( double f, double la ); + static float d_factorfloat( float f, float la ); + static double calculate_fl_from_la_ciecam02( double la ); + static float calculate_fl_from_la_ciecam02float( float la ); + static double achromatic_response_to_white( double x, double y, double z, double d, double fl, double nbb, int gamu ); + static float achromatic_response_to_whitefloat( float x, float y, float z, float d, float fl, float nbb, int gamu ); + static void xyz_to_cat02 ( double &r, double &g, double &b, double x, double y, double z, int gamu ); + static void cat02_to_hpe ( double &rh, double &gh, double &bh, double r, double g, double b, int gamu ); + static void cat02_to_xyz ( double &x, double &y, double &z, double r, double g, double b, int gamu ); + static void hpe_to_xyz ( double &x, double &y, double &z, double r, double g, double b ); + + static void xyz_to_cat02float ( float &r, float &g, float &b, float x, float y, float z, int gamu ); + static void cat02_to_hpefloat ( float &rh, float &gh, float &bh, float r, float g, float b, int gamu ); + static void cat02_to_xyzfloat ( float &x, float &y, float &z, float r, float g, float b, int gamu ); + static void hpe_to_xyzfloat ( float &x, float &y, float &z, float r, float g, float b ); + + static void Aab_to_rgb( double &r, double &g, double &b, double A, double aa, double bb, double nbb ); + static void Aab_to_rgbfloat( float &r, float &g, float &b, float A, float aa, float bb, float nbb ); + static void calculate_ab( double &aa, double &bb, double h, double e, double t, double nbb, double a ); + static void calculate_abfloat( float &aa, float &bb, float h, float e, float t, float nbb, float a ); + + static double nonlinear_adaptation( double c, double fl ); + static float nonlinear_adaptationfloat( float c, float fl ); + static double inverse_nonlinear_adaptation( double c, double fl ); + static float inverse_nonlinear_adaptationfloat( float c, float fl ); + +public: + Ciecam02 () {} + static void curvecolor(double satind, double satval, double &sres, double parsat); + static void curvecolorfloat(float satind, float satval, float &sres, float parsat); + static void curveJ (double br, double contr, int db, LUTf & outCurve , LUTu & histogram ) ; + static void curveJfloat (float br, float contr, int db, LUTf & outCurve , LUTu & histogram ) ; + + /** + * Inverse transform from CIECAM02 JCh to XYZ. + */ + static void jch2xyz_ciecam02( double &x, double &y, double &z, + double J, double C, double h, + double xw, double yw, double zw, + double yb, double la, + double f, double c, double nc, int gamu, double n, double nbb, double ncb, double fl, double cz, double d, double aw); + + static void jch2xyz_ciecam02float( float &x, float &y, float &z, + float J, float C, float h, + float xw, float yw, float zw, + float yb, float la, + float f, float c, float nc,int gamu,float n, float nbb, float ncb, float fl, float cz, float d, float aw ); + + /** + * Forward transform from XYZ to CIECAM02 JCh. + */ + static void initcam1(double gamu, double yb, double pilotd, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, + double &cz, double &aw, double &wh, double &pfl, double &fl, double &c); + + static void initcam2(double gamu, double yb, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, + double &cz, double &aw, double &fl); + + static void initcam1float(float gamu, float yb, float pilotd, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, + float &cz, float &aw, float &wh, float &pfl, float &fl, float &c); + + static void initcam2float(float gamu, float yb, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, + float &cz, float &aw, float &fl); + + static void xyz2jchqms_ciecam02( double &J, double &C, double &h, + double &Q, double &M, double &s,double &aw, double &fl, double &wh, + double x, double y, double z, + double xw, double yw, double zw, + double yb, double la, + double f, double c, double nc, double pilotd,int gamu , double n, double nbb, double ncb, double pfl, double cz, double d ); + + static void xyz2jchqms_ciecam02float( float &J, float &C, float &h, + float &Q, float &M, float &s,float &aw, float &fl, float &wh, + float x, float y, float z, + float xw, float yw, float zw, + float yb, float la, + float f, float c, float nc, float pilotd, int gamu, float n, float nbb, float ncb, float pfl, float cz, float d ); + + +}; +} +#endif diff --git a/rtengine/colortemp.cc b/rtengine/colortemp.cc index 3901345d9..db473ec70 100644 --- a/rtengine/colortemp.cc +++ b/rtengine/colortemp.cc @@ -17,16 +17,12 @@ * along with RawTherapee. If not, see . */ #include "colortemp.h" -#include "iccmatrices.h" #include "rtengine.h" -#include "improcfun.h" -#include "curves.h" #include #include #include -#include "mytime.h" #include "sleef.c" -#include "../rtgui/options.h" +#include "settings.h" #undef CLIPD #define CLIPD(a) ((a)>0.0?((a)<1.0?(a):1.0):0.0) @@ -34,7 +30,6 @@ #define MAXR(a,b) ((a) > (b) ? (a) : (b)) namespace rtengine { - using namespace procparams; extern const Settings* settings; @@ -914,221 +909,6 @@ int ColorTemp::XYZtoCorColorTemp(double x0, double y0, double z0, double &temp) return 0; /* success */ } -void ColorTemp::curvecolor(double satind, double satval, double &sres, double parsat) { - if (satind >=0.0) { - sres = (1.-(satind)/100.)*satval+(satind)/100.*(1.-SQR(SQR(1.-min(satval,1.0)))); - if (sres>parsat) sres=parsat; - if (sres<0.) sres=0.; - } else { - if (satind < -0.1) - sres = satval*(1.+(satind)/100.); - } -} -void ColorTemp::curvecolorfloat(float satind, float satval, float &sres, float parsat) { - if (satind >=0.0f) { - sres = (1.f-(satind)/100.f)*satval+(satind)/100.f*(1.f-SQR(SQR(1.f-min(satval,1.0f)))); - if (sres>parsat) sres=parsat; - if (sres<0.f) sres=0.f; - } else { - sres = satval*(1.f+(satind)/100.f); - } -} - - -void ColorTemp::curveJ (double br, double contr, int db, LUTf & outCurve, LUTu & histogram ) { - LUTf dcurve(65536,0); - int skip=1; - - // check if brightness curve is needed - if (br>0.00001 || br<-0.00001) { - - std::vector brightcurvePoints; - brightcurvePoints.resize(9); - brightcurvePoints.at(0) = double(DCT_NURBS); - - brightcurvePoints.at(1) = 0.; // black point. Value in [0 ; 1] range - brightcurvePoints.at(2) = 0.; // black point. Value in [0 ; 1] range - - if (br>0) { - brightcurvePoints.at(3) = 0.1; // toe point - brightcurvePoints.at(4) = 0.1+br/150.0; //value at toe point - - brightcurvePoints.at(5) = 0.7; // shoulder point - brightcurvePoints.at(6) = min(1.0,0.7+br/300.0); //value at shoulder point - } else { - brightcurvePoints.at(3) = 0.1-br/150.0; // toe point - brightcurvePoints.at(4) = 0.1; // value at toe point - - brightcurvePoints.at(5) = min(1.0,0.7-br/300.0); // shoulder point - brightcurvePoints.at(6) = 0.7; // value at shoulder point - } - brightcurvePoints.at(7) = 1.; // white point - brightcurvePoints.at(8) = 1.; // value at white point - - DiagonalCurve* brightcurve = new DiagonalCurve (brightcurvePoints, CURVES_MIN_POLY_POINTS/skip); - - // Applying brightness curve - for (int i=0; i<32768; i++) { - - // change to [0,1] range - float val = (float)i / 32767.0; - - // apply brightness curve - val = brightcurve->getVal (val); - - // store result in a temporary array - dcurve[i] = CLIPD(val); - } - delete brightcurve; - } - else { - // for (int i=0; i<32768; i++) { // L values range up to 32767, higher values are for highlight overflow - for (int i=0; i<(32768*db); i++) { // L values range up to 32767, higher values are for highlight overflow - - // set the identity curve in the temporary array - dcurve[i] = (float)i / (db*32768.0f); - } - } - - - if (contr>0.00001 || contr<-0.00001) { - - // compute mean luminance of the image with the curve applied - int sum = 0; - float avg = 0; - //float sqavg = 0; - for (int i=0; i<32768; i++) { - avg += dcurve[i] * histogram[i];//approximation for average : usage of L (lab) instead of J - sum += histogram[i]; - } - avg /= sum; - std::vector contrastcurvePoints; - contrastcurvePoints.resize(9); - contrastcurvePoints.at(0) = double(DCT_NURBS); - - contrastcurvePoints.at(1) = 0.; // black point. Value in [0 ; 1] range - contrastcurvePoints.at(2) = 0.; // black point. Value in [0 ; 1] range - - contrastcurvePoints.at(3) = avg-avg*(0.6-contr/250.0); // toe point - contrastcurvePoints.at(4) = avg-avg*(0.6+contr/250.0); // value at toe point - - contrastcurvePoints.at(5) = avg+(1-avg)*(0.6-contr/250.0); // shoulder point - contrastcurvePoints.at(6) = avg+(1-avg)*(0.6+contr/250.0); // value at shoulder point - - contrastcurvePoints.at(7) = 1.; // white point - contrastcurvePoints.at(8) = 1.; // value at white point - - DiagonalCurve* contrastcurve = new DiagonalCurve (contrastcurvePoints, CURVES_MIN_POLY_POINTS/skip); - - // apply contrast enhancement - for (int i=0; i<(32768*db); i++) { - dcurve[i] = contrastcurve->getVal (dcurve[i]); - } - - delete contrastcurve; - } - - // for (int i=0; i<32768; i++) outCurve[i] = 32768.0*dcurve[i]; - for (int i=0; i<(db*32768); i++) outCurve[i] = db*32768.0*dcurve[i]; -} -void ColorTemp::curveJfloat (float br, float contr, int db, LUTf & outCurve, LUTu & histogram ) { - LUTf dcurve(65536,0); - int skip=1; - - // check if brightness curve is needed - if (br>0.00001f || br<-0.00001f) { - - std::vector brightcurvePoints; - brightcurvePoints.resize(9); - brightcurvePoints.at(0) = double(DCT_NURBS); - - brightcurvePoints.at(1) = 0.f; // black point. Value in [0 ; 1] range - brightcurvePoints.at(2) = 0.f; // black point. Value in [0 ; 1] range - - if (br>0) { - brightcurvePoints.at(3) = 0.1f; // toe point - brightcurvePoints.at(4) = 0.1f+br/150.0f; //value at toe point - - brightcurvePoints.at(5) = 0.7f; // shoulder point - brightcurvePoints.at(6) = min(1.0f,0.7f+br/300.0f); //value at shoulder point - } else { - brightcurvePoints.at(3) = 0.1f-br/150.0f; // toe point - brightcurvePoints.at(4) = 0.1f; // value at toe point - - brightcurvePoints.at(5) = min(1.0f,0.7f-br/300.0f); // shoulder point - brightcurvePoints.at(6) = 0.7f; // value at shoulder point - } - brightcurvePoints.at(7) = 1.f; // white point - brightcurvePoints.at(8) = 1.f; // value at white point - - DiagonalCurve* brightcurve = new DiagonalCurve (brightcurvePoints, CURVES_MIN_POLY_POINTS/skip); - - // Applying brightness curve - for (int i=0; i<32768; i++) { - - // change to [0,1] range - float val = (float)i / 32767.0f; - - // apply brightness curve - val = brightcurve->getVal (val); - - // store result in a temporary array - dcurve[i] = CLIPD(val); - } - delete brightcurve; - } - else { - // for (int i=0; i<32768; i++) { // L values range up to 32767, higher values are for highlight overflow - for (int i=0; i<(32768*db); i++) { // L values range up to 32767, higher values are for highlight overflow - - // set the identity curve in the temporary array - dcurve[i] = (float)i / (db*32768.0f); - } - } - - - if (contr>0.00001f || contr<-0.00001f) { - - // compute mean luminance of the image with the curve applied - int sum = 0; - float avg = 0; - //float sqavg = 0; - for (int i=0; i<32768; i++) { - avg += dcurve[i] * histogram[i];//approximation for average : usage of L (lab) instead of J - sum += histogram[i]; - } - avg /= sum; - //printf("avg=%f\n",avg); - std::vector contrastcurvePoints; - contrastcurvePoints.resize(9); - contrastcurvePoints.at(0) = double(DCT_NURBS); - - contrastcurvePoints.at(1) = 0.f; // black point. Value in [0 ; 1] range - contrastcurvePoints.at(2) = 0.f; // black point. Value in [0 ; 1] range - - contrastcurvePoints.at(3) = avg-avg*(0.6f-contr/250.0f); // toe point - contrastcurvePoints.at(4) = avg-avg*(0.6f+contr/250.0f); // value at toe point - - contrastcurvePoints.at(5) = avg+(1-avg)*(0.6f-contr/250.0f); // shoulder point - contrastcurvePoints.at(6) = avg+(1-avg)*(0.6f+contr/250.0f); // value at shoulder point - - contrastcurvePoints.at(7) = 1.f; // white point - contrastcurvePoints.at(8) = 1.f; // value at white point - - DiagonalCurve* contrastcurve = new DiagonalCurve (contrastcurvePoints, CURVES_MIN_POLY_POINTS/skip); - - // apply contrast enhancement - for (int i=0; i<(32768*db); i++) { - dcurve[i] = contrastcurve->getVal (dcurve[i]); - } - - delete contrastcurve; - } - - // for (int i=0; i<32768; i++) outCurve[i] = 32768.0*dcurve[i]; - for (int i=0; i<(db*32768); i++) outCurve[i] = db*32768.0f*dcurve[i]; -} - void ColorTemp::cieCAT02(double Xw, double Yw, double Zw, double &CAM02BB00,double &CAM02BB01,double &CAM02BB02,double &CAM02BB10,double &CAM02BB11,double &CAM02BB12,double &CAM02BB20,double &CAM02BB21,double &CAM02BB22, double adap ) { // CIECAT02 - J.Desmis January 2012 review September 2012 @@ -1706,512 +1486,6 @@ void ColorTemp::temp2mul (double temp, double green, double equal, double& rmul, } } } -/* J.Desmis october 2012 - CIECAM02 : -I have notably modified the code of Billy Biggs to adapt to RT -The code is update from last works of M.Fairchild -and also gamut correction M.H.Brill S.Susstrunk -*/ - -/** - * Copyright (c) 2003 Billy Biggs - * - * Permission is hereby granted, free of charge, to any person obtaining - * a copy of this software and associated documentation files (the - * "Software"), to deal in the Software without restriction, including - * without limitation the rights to use, copy, modify, merge, publish, - * distribute, sublicense, and/or sell copies of the Software, and to - * permit persons to whom the Software is furnished to do so, subject to - * the following conditions: - * - * The above copyright notice and this permission notice shall be - * included in all copies or substantial portions of the Software. - * - * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, - * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF - * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND - * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS - * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN - * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN - * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE - * SOFTWARE. - * - */ -void ColorTemp::xyz_to_cat02( double &r, double &g, double &b, double x, double y, double z, int gamu ) -{ -gamu=1; -if(gamu==0){ - r = ( 0.7328 * x) + (0.4296 * y) - (0.1624 * z); - g = (-0.7036 * x) + (1.6975 * y) + (0.0061 * z); - b = ( 0.0030 * x) + (0.0136 * y) + (0.9834 * z); - } -else if (gamu==1) {//gamut correction M.H.Brill S.Susstrunk - //r = ( 0.7328 * x) + (0.4296 * y) - (0.1624 * z); - //g = (-0.7036 * x) + (1.6975 * y) + (0.0061 * z); - //b = ( 0.0000 * x) + (0.0000 * y) + (1.0000 * z); - r = ( 1.007245 * x) + (0.011136* y) - (0.018381 * z);//Changjun Li - g = (-0.318061 * x) + (1.314589 * y) + (0.003471 * z); - b = ( 0.0000 * x) + (0.0000 * y) + (1.0000 * z); - -} -} - - void ColorTemp::xyz_to_cat02float( float &r, float &g, float &b, float x, float y, float z, int gamu ) -{ -gamu=1; -if(gamu==0){ - r = ( 0.7328f * x) + (0.4296f * y) - (0.1624f * z); - g = (-0.7036f * x) + (1.6975f * y) + (0.0061f * z); - b = ( 0.0030f * x) + (0.0136f * y) + (0.9834f * z); - } -else if (gamu==1) {//gamut correction M.H.Brill S.Susstrunk - //r = ( 0.7328 * x) + (0.4296 * y) - (0.1624 * z); - //g = (-0.7036 * x) + (1.6975 * y) + (0.0061 * z); - //b = ( 0.0000 * x) + (0.0000 * y) + (1.0000 * z); - r = ( 1.007245f * x) + (0.011136f* y) - (0.018381f * z);//Changjun Li - g = (-0.318061f * x) + (1.314589f * y) + (0.003471f * z); - b = ( 0.0000f * x) + (0.0000f * y) + (1.0000f * z); - -} -} - - - void ColorTemp::cat02_to_xyz( double &x, double &y, double &z, double r, double g, double b, int gamu ) -{ -gamu=1; -if(gamu==0) { - x = ( 1.096124 * r) - (0.278869 * g) + (0.182745 * b); - y = ( 0.454369 * r) + (0.473533 * g) + (0.072098 * b); - z = (-0.009628 * r) - (0.005698 * g) + (1.015326 * b); - } -else if(gamu==1) {//gamut correction M.H.Brill S.Susstrunk - //x = ( 1.0978566 * r) - (0.277843 * g) + (0.179987 * b); - //y = ( 0.455053 * r) + (0.473938 * g) + (0.0710096* b); - //z = ( 0.000000 * r) - (0.000000 * g) + (1.000000 * b); - x = ( 0.99015849 * r) - (0.00838772* g) + (0.018229217 * b);//Changjun Li - y = ( 0.239565979 * r) + (0.758664642 * g) + (0.001770137* b); - z = ( 0.000000 * r) - (0.000000 * g) + (1.000000 * b); - -} -} - - void ColorTemp::cat02_to_xyzfloat( float &x, float &y, float &z, float r, float g, float b, int gamu ) -{ -gamu=1; -if(gamu==0) { - x = ( 1.096124f * r) - (0.278869f * g) + (0.182745f * b); - y = ( 0.454369f * r) + (0.473533f * g) + (0.072098f * b); - z = (-0.009628f * r) - (0.005698f * g) + (1.015326f * b); - } -else if(gamu==1) {//gamut correction M.H.Brill S.Susstrunk - //x = ( 1.0978566 * r) - (0.277843 * g) + (0.179987 * b); - //y = ( 0.455053 * r) + (0.473938 * g) + (0.0710096* b); - //z = ( 0.000000 * r) - (0.000000 * g) + (1.000000 * b); - x = ( 0.99015849f * r) - (0.00838772f* g) + (0.018229217f * b);//Changjun Li - y = ( 0.239565979f * r) + (0.758664642f * g) + (0.001770137f* b); - z = ( 0.000000f * r) - (0.000000f * g) + (1.000000f * b); - -} -} - - -void ColorTemp::hpe_to_xyz( double &x, double &y, double &z, double r, double g, double b ) -{ - x = (1.910197 * r) - (1.112124 * g) + (0.201908 * b); - y = (0.370950 * r) + (0.629054 * g) - (0.000008 * b); - z = b; - -} - - void ColorTemp::hpe_to_xyzfloat( float &x, float &y, float &z, float r, float g, float b ) -{ - x = (1.910197f * r) - (1.112124f * g) + (0.201908f * b); - y = (0.370950f * r) + (0.629054f * g) - (0.000008f * b); - z = b; - -} - - - -void ColorTemp::cat02_to_hpe( double &rh, double &gh, double &bh, double r, double g, double b, int gamu ) -{ gamu=1; - if(gamu==0){ - rh = ( 0.7409792 * r) + (0.2180250 * g) + (0.0410058 * b); - gh = ( 0.2853532 * r) + (0.6242014 * g) + (0.0904454 * b); - bh = (-0.0096280 * r) - (0.0056980 * g) + (1.0153260 * b); - } - else if (gamu==1) {//Changjun Li - rh = ( 0.550930835 * r) + (0.519435987* g) - ( 0.070356303* b); - gh = ( 0.055954056 * r) + (0.89973132 * g) + (0.044315524 * b); - bh = (0.0 * r) - (0.0* g) + (1.0 * b); - } -} -void ColorTemp::cat02_to_hpefloat( float &rh, float &gh, float &bh, float r, float g, float b, int gamu ) -{ gamu=1; - if(gamu==0){ - rh = ( 0.7409792f * r) + (0.2180250f * g) + (0.0410058f * b); - gh = ( 0.2853532f * r) + (0.6242014f * g) + (0.0904454f * b); - bh = (-0.0096280f * r) - (0.0056980f * g) + (1.0153260f * b); - } - else if (gamu==1) {//Changjun Li - rh = ( 0.550930835f * r) + (0.519435987f* g) - ( 0.070356303f* b); - gh = ( 0.055954056f * r) + (0.89973132f * g) + (0.044315524f * b); - bh = (0.0f * r) - (0.0f* g) + (1.0f * b); - } -} - - - void ColorTemp::Aab_to_rgb( double &r, double &g, double &b, double A, double aa, double bb, double nbb ) -{ - double x = (A / nbb) + 0.305; - - /* c1 c2 c3 */ - r = (0.32787 * x) + (0.32145 * aa) + (0.20527 * bb); - /* c1 c4 c5 */ - g = (0.32787 * x) - (0.63507 * aa) - (0.18603 * bb); - /* c1 c6 c7 */ - b = (0.32787 * x) - (0.15681 * aa) - (4.49038 * bb); -} - void ColorTemp::Aab_to_rgbfloat( float &r, float &g, float &b, float A, float aa, float bb, float nbb ) -{ - float x = (A / nbb) + 0.305f; - - /* c1 c2 c3 */ - r = (0.32787f * x) + (0.32145f * aa) + (0.20527f * bb); - /* c1 c4 c5 */ - g = (0.32787f * x) - (0.63507f * aa) - (0.18603f * bb); - /* c1 c6 c7 */ - b = (0.32787f * x) - (0.15681f * aa) - (4.49038f * bb); -} - -void ColorTemp::calculate_ab( double &aa, double &bb, double h, double e, double t, double nbb, double a ) -{ - double hrad = (h * M_PI) / 180.0; - double sinh = sin( hrad ); - double cosh = cos( hrad ); - double x = (a / nbb) + 0.305; - double p3 = 21.0/20.0; - if( fabs( sinh ) >= fabs( cosh ) ) { - bb = ((0.32787 * x) * (2.0 + p3)) / - ((e / (t * sinh)) - - // ((0.32145 - 0.63507 - (p3 * 0.15681)) * (cosh / sinh)) - - // (0.20527 - 0.18603 - (p3 * 4.49038))); - ((-0.31362 - (p3 * 0.15681)) * (cosh / sinh)) - - (0.01924 - (p3 * 4.49038))); - - aa = (bb * cosh) / sinh; - } else { - aa = ((0.32787 * x) * (2.0 + p3)) / - ((e / (t * cosh)) - - // (0.32145 - 0.63507 - (p3 * 0.15681)) - - // ((0.20527 - 0.18603 - (p3 * 4.49038)) * (sinh / cosh))); - (-0.31362 - (p3 * 0.15681)) - - ((0.01924 - (p3 * 4.49038)) * (sinh / cosh))); - - bb = (aa * sinh) / cosh; - } -} -void ColorTemp::calculate_abfloat( float &aa, float &bb, float h, float e, float t, float nbb, float a ) -{ - float2 sincosval = xsincosf((h * M_PI) / 180.0f); - float sinh = sincosval.x; - float cosh = sincosval.y; - float x = (a / nbb) + 0.305f; -// float p3 = 21.0f/20.0f; - float p3 = 1.05f; - if( fabs( sinh ) >= fabs( cosh ) ) { - bb = ((0.32787f * x) * (2.0f + p3)) / - ((e / (t * sinh)) - - // ((0.32145 - 0.63507 - (p3 * 0.15681)) * (cosh / sinh)) - - // (0.20527 - 0.18603 - (p3 * 4.49038))); - ((-0.31362f - (p3 * 0.15681f)) * (cosh / sinh)) - - (0.01924f - (p3 * 4.49038f))); - - aa = (bb * cosh) / sinh; - } else { - aa = ((0.32787f * x) * (2.0f + p3)) / - ((e / (t * cosh)) - - // (0.32145 - 0.63507 - (p3 * 0.15681)) - - // ((0.20527 - 0.18603 - (p3 * 4.49038)) * (sinh / cosh))); - (-0.31362f - (p3 * 0.15681f)) - - ((0.01924f - (p3 * 4.49038f)) * (sinh / cosh))); - - bb = (aa * sinh) / cosh; - } -} - - -void ColorTemp::initcam1(double gamu, double yb, double pilotd, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, - double &cz, double &aw, double &wh, double &pfl, double &fl, double &c) -{ - n = yb / yw; - if(pilotd==2.0) d = d_factor( f, la );else d=pilotd; - fl = calculate_fl_from_la_ciecam02( la ); - nbb = ncb = 0.725 * pow( 1.0 / n, 0.2 ); - cz = 1.48 + sqrt( n ); - aw = achromatic_response_to_white( xw, yw, zw, d, fl, nbb, gamu ); - wh =( 4.0 / c ) * ( aw + 4.0 ) * pow( fl, 0.25 ); - pfl = pow( fl, 0.25 ); -#ifdef _DEBUG - if (settings->verbose) printf("Source double d=%f aw=%f fl=%f wh=%f\n",d,aw,fl,wh); -#endif - -} -void ColorTemp::initcam1float(float gamu, float yb, float pilotd, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, - float &cz, float &aw, float &wh, float &pfl, float &fl, float &c) -{ - n = yb / yw; - if(pilotd==2.0) d = d_factorfloat( f, la );else d=pilotd; - fl = calculate_fl_from_la_ciecam02float( la ); - nbb = ncb = 0.725f * pow_F( 1.0f / n, 0.2f ); - cz = 1.48f + sqrt( n ); - aw = achromatic_response_to_whitefloat( xw, yw, zw, d, fl, nbb, gamu ); - wh =( 4.0f / c ) * ( aw + 4.0f ) * pow_F( fl, 0.25f ); - pfl = pow_F( fl, 0.25f ); -#ifdef _DEBUG - if (settings->verbose) printf("Source float d=%f aw=%f fl=%f wh=%f c=%f awc=%f\n",d,aw,fl,wh,c,(4.f/c)*(aw+4.f)); -#endif - -} - -void ColorTemp::initcam2(double gamu, double yb, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, - double &cz, double &aw, double &fl) -{ - n = yb / yw; - d = d_factor( f, la ); - fl = calculate_fl_from_la_ciecam02( la ); - nbb = ncb = 0.725 * pow( 1.0 / n, 0.2 ); - cz = 1.48 + sqrt( n ); - aw = achromatic_response_to_white( xw, yw, zw, d, fl, nbb, gamu ); -#ifdef _DEBUG - if (settings->verbose) printf("Viewing double d=%f aw=%f fl=%f n=%f\n",d,aw,fl,n); -#endif -} - -void ColorTemp::initcam2float(float gamu, float yb, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, - float &cz, float &aw, float &fl) -{ - n = yb / yw; - d = d_factorfloat( f, la ); - fl = calculate_fl_from_la_ciecam02float( la ); - nbb = ncb = 0.725f * pow_F( 1.0f / n, 0.2f ); - cz = 1.48f + sqrt( n ); - aw = achromatic_response_to_whitefloat( xw, yw, zw, d, fl, nbb, gamu ); -#ifdef _DEBUG - if (settings->verbose) printf("Viewing float d=%f aw=%f fl=%f n=%f\n",d,aw,fl,n); -#endif - -} - - - -void ColorTemp::xyz2jchqms_ciecam02( double &J, double &C, double &h, double &Q, double &M, double &s,double &aw, double &fl, double &wh, - double x, double y, double z, double xw, double yw, double zw, - double yb, double la, double f, double c, double nc, double pilotd, int gamu , double n, double nbb, double ncb, double pfl, double cz, double d) -{ - double r, g, b; - double rw, gw, bw; - double rc, gc, bc; - double rp, gp, bp; - double rpa, gpa, bpa; - double a, ca, cb; - double e, t; - double myh, myj, myc, myq, mym, mys; - gamu=1; - xyz_to_cat02( r, g, b, x, y, z, gamu ); - xyz_to_cat02( rw, gw, bw, xw, yw, zw, gamu ); - rc = r * (((yw * d) / rw) + (1.0 - d)); - gc = g * (((yw * d) / gw) + (1.0 - d)); - bc = b * (((yw * d) / bw) + (1.0 - d)); - - ColorTemp::cat02_to_hpe( rp, gp, bp, rc, gc, bc, gamu ); - if(gamu==1){//gamut correction M.H.Brill S.Susstrunk - rp=MAXR(rp,0.0); - gp=MAXR(gp,0.0); - bp=MAXR(bp,0.0); - } - rpa = nonlinear_adaptation( rp, fl ); - gpa = nonlinear_adaptation( gp, fl ); - bpa = nonlinear_adaptation( bp, fl ); - - ca = rpa - ((12.0 * gpa) / 11.0) + (bpa / 11.0); - cb = (1.0 / 9.0) * (rpa + gpa - (2.0 * bpa)); - - myh = (180.0 / M_PI) * atan2( cb, ca ); - if( myh < 0.0 ) myh += 360.0; - //we can also calculate H, if necessary...but it's using time...for what usage ? - /*double temp; - if(myh<20.14) { - temp = ((myh + 122.47)/1.2) + ((20.14 - myh)/0.8); - H = 300 + (100*((myh + 122.47)/1.2)) / temp; - } - else if(myh < 90.0) { - temp = ((myh - 20.14)/0.8) + ((90.00 - myh)/0.7); - H = (100*((myh - 20.14)/0.8)) / temp; - } - else if (myh < 164.25) { - temp = ((myh - 90.00)/0.7) + ((164.25 - myh)/1.0); - H = 100 + ((100*((myh - 90.00)/0.7)) / temp); - } - else if (myh < 237.53) { - temp = ((myh - 164.25)/1.0) + ((237.53 - myh)/1.2); - H = 200 + ((100*((myh - 164.25)/1.0)) / temp); - } - else { - temp = ((myh - 237.53)/1.2) + ((360 - myh + 20.14)/0.8); - H = 300 + ((100*((myh - 237.53)/1.2)) / temp); - } - */ - a = ((2.0 * rpa) + gpa + ((1.0 / 20.0) * bpa) - 0.305) * nbb; - if (gamu==1) a=MAXR(a,0.0);//gamut correction M.H.Brill S.Susstrunk - J = 100.0 * pow( a / aw, c * cz ); - - e = ((12500.0 / 13.0) * nc * ncb) * (cos( ((myh * M_PI) / 180.0) + 2.0 ) + 3.8); - t = (e * sqrt( (ca * ca) + (cb * cb) )) / (rpa + gpa + ((21.0 / 20.0) * bpa)); - - C = pow( t, 0.9 ) * sqrt( J / 100.0 ) - * pow( 1.64 - pow( 0.29, n ), 0.73 ); - - Q = wh * sqrt( J / 100.0 ); - M = C * pfl; - s = 100.0 * sqrt( M / Q ); - h = myh; - -} - - -void ColorTemp::xyz2jchqms_ciecam02float( float &J, float &C, float &h, float &Q, float &M, float &s,float &aw, float &fl, float &wh, - float x, float y, float z, float xw, float yw, float zw, - float yb, float la, float f, float c, float nc, float pilotd, int gamu, float pow1, float nbb, float ncb, float pfl, float cz, float d) - - { - float r, g, b; - float rw, gw, bw; - float rc, gc, bc; - float rp, gp, bp; - float rpa, gpa, bpa; - float a, ca, cb; - float e, t; - float myh, myj, myc, myq, mym, mys; - gamu=1; - xyz_to_cat02float( r, g, b, x, y, z, gamu ); - xyz_to_cat02float( rw, gw, bw, xw, yw, zw, gamu ); - rc = r * (((yw * d) / rw) + (1.0 - d)); - gc = g * (((yw * d) / gw) + (1.0 - d)); - bc = b * (((yw * d) / bw) + (1.0 - d)); - - ColorTemp::cat02_to_hpefloat( rp, gp, bp, rc, gc, bc, gamu ); - if(gamu==1) {//gamut correction M.H.Brill S.Susstrunk - rp=MAXR(rp,0.0f); - gp=MAXR(gp,0.0f); - bp=MAXR(bp,0.0f); - } - rpa = nonlinear_adaptationfloat( rp, fl ); - gpa = nonlinear_adaptationfloat( gp, fl ); - bpa = nonlinear_adaptationfloat( bp, fl ); - - ca = rpa - ((12.0f * gpa) - bpa) / 11.0f; - cb = (0.11111111f) * (rpa + gpa - (2.0f * bpa)); - - myh = xatan2f( cb, ca ); - if( myh < 0.0f ) - myh += (2.f * M_PI); - - a = ((2.0f * rpa) + gpa + (0.05f * bpa) - 0.305f) * nbb; - if (gamu==1) - a=MAXR(a,0.0f);//gamut correction M.H.Brill S.Susstrunk - - J = pow_F( a / aw, c * cz * 0.5f); - - e = ((961.53846f) * nc * ncb) * (xcosf( myh + 2.0f ) + 3.8f); - t = (e * sqrtf( (ca * ca) + (cb * cb) )) / (rpa + gpa + (1.05f * bpa)); - - C = pow_F( t, 0.9f ) * J * pow1; - - Q = wh * J; - J *= J * 100.0f; - M = C * pfl; - Q = (Q == 0.f ? 0.0001f : Q); // avoid division by zero - s = 100.0f * sqrtf( M / Q ); - h = (myh * 180.f) / (float)M_PI; - -} - - -void ColorTemp::jch2xyz_ciecam02( double &x, double &y, double &z, double J, double C, double h, - double xw, double yw, double zw, double yb, double la, - double f, double c, double nc , int gamu, double n, double nbb, double ncb, double fl, double cz, double d, double aw ) -{ - double r, g, b; - double rc, gc, bc; - double rp, gp, bp; - double rpa, gpa, bpa; - double rw, gw, bw; - double a, ca, cb; - double e, t; - gamu=1; - ColorTemp::xyz_to_cat02( rw, gw, bw, xw, yw, zw, gamu ); - e = ((12500.0 / 13.0) * nc * ncb) * (cos( ((h * M_PI) / 180.0) + 2.0 ) + 3.8); - a = pow( J / 100.0, 1.0 / (c * cz) ) * aw; - t = pow( C / (sqrt( J / 100) * pow( 1.64 - pow( 0.29, n ), 0.73 )), 10.0 / 9.0 ); - - calculate_ab( ca, cb, h, e, t, nbb, a ); - Aab_to_rgb( rpa, gpa, bpa, a, ca, cb, nbb ); - - rp = inverse_nonlinear_adaptation( rpa, fl ); - gp = inverse_nonlinear_adaptation( gpa, fl ); - bp = inverse_nonlinear_adaptation( bpa, fl ); - - ColorTemp::hpe_to_xyz( x, y, z, rp, gp, bp ); - ColorTemp::xyz_to_cat02( rc, gc, bc, x, y, z, gamu ); - - r = rc / (((yw * d) / rw) + (1.0 - d)); - g = gc / (((yw * d) / gw) + (1.0 - d)); - b = bc / (((yw * d) / bw) + (1.0 - d)); - - ColorTemp::cat02_to_xyz( x, y, z, r, g, b, gamu ); -} - -void ColorTemp::jch2xyz_ciecam02float( float &x, float &y, float &z, float J, float C, float h, - float xw, float yw, float zw, float yb, float la, - float f, float c, float nc , int gamu, float pow1, float nbb, float ncb, float fl, float cz, float d, float aw) - -{ - float r, g, b; - float rc, gc, bc; - float rp, gp, bp; - float rpa, gpa, bpa; - float rw, gw, bw; - float a, ca, cb; - float e, t; - gamu=1; - ColorTemp::xyz_to_cat02float( rw, gw, bw, xw, yw, zw, gamu ); - e = ((961.53846f) * nc * ncb) * (xcosf( ((h * M_PI) / 180.0f) + 2.0f ) + 3.8f); - a = pow_F( J / 100.0f, 1.0f / (c * cz) ) * aw; - t = pow_F( 10.f * C / (sqrtf( J ) * pow1), 1.1111111f ); - - calculate_abfloat( ca, cb, h, e, t, nbb, a ); - Aab_to_rgbfloat( rpa, gpa, bpa, a, ca, cb, nbb ); - - rp = inverse_nonlinear_adaptationfloat( rpa, fl ); - gp = inverse_nonlinear_adaptationfloat( gpa, fl ); - bp = inverse_nonlinear_adaptationfloat( bpa, fl ); - - ColorTemp::hpe_to_xyzfloat( x, y, z, rp, gp, bp ); - ColorTemp::xyz_to_cat02float( rc, gc, bc, x, y, z, gamu ); - - r = rc / (((yw * d) / rw) + (1.0f - d)); - g = gc / (((yw * d) / gw) + (1.0f - d)); - b = bc / (((yw * d) / bw) + (1.0f - d)); - - ColorTemp::cat02_to_xyzfloat( x, y, z, r, g, b, gamu ); -} - - -//end CIECAM Billy Bigg - - - - /* Calculate Planck's radiation diff --git a/rtengine/colortemp.h b/rtengine/colortemp.h index 9d92cf60c..7b1686ed8 100644 --- a/rtengine/colortemp.h +++ b/rtengine/colortemp.h @@ -18,19 +18,13 @@ */ #ifndef _COLORTEMP_ #define _COLORTEMP_ -#include "procparams.h" -#include +#include #include -#include "sleef.c" -#include "LUT.h" -#define MAXR(a,b) ((a) > (b) ? (a) : (b)) + #define pow_F(a,b) (xexpf(b*xlogf(a))) - namespace rtengine { -using namespace procparams; - #define MINTEMP 1500 #define MAXTEMP 60000 @@ -78,136 +72,6 @@ class ColorTemp { //static void CAT02 (Imagefloat* baseImg, const ProcParams* params); //static void ciecam_02 (LabImage* lab, const ProcParams* params); - static double d_factor( double f, double la ) { - return f * (1.0 - ((1.0 / 3.6) * exp((-la - 42.0) / 92.0))); - } - static float d_factorfloat( float f, float la ) { - return f * (1.0f - ((1.0f / 3.6f) * xexpf((-la - 42.0f) / 92.0f))); - } - - static double calculate_fl_from_la_ciecam02( double la ) { - double la5 = la * 5.0; - double k = 1.0 / (la5 + 1.0); - - /* Calculate k^4. */ - k = k * k; - k = k * k; - - return (0.2 * k * la5) + (0.1 * (1.0 - k) * (1.0 - k) * pow(la5, 1.0 / 3.0)); - } - static float calculate_fl_from_la_ciecam02float( float la ) { - float la5 = la * 5.0f; - float k = 1.0f / (la5 + 1.0f); - - /* Calculate k^4. */ - k = k * k; - k = k * k; - - return (0.2f * k * la5) + (0.1f * (1.0f - k) * (1.0f - k) * pow_F(la5, 1.0f / 3.0f)); - } - - static double achromatic_response_to_white( double x, double y, double z, double d, double fl, double nbb, int gamu ) { - double r, g, b; - double rc, gc, bc; - double rp, gp, bp; - double rpa, gpa, bpa; - gamu=1; - xyz_to_cat02( r, g, b, x, y, z, gamu ); - - rc = r * (((y * d) / r) + (1.0 - d)); - gc = g * (((y * d) / g) + (1.0 - d)); - bc = b * (((y * d) / b) + (1.0 - d)); - - cat02_to_hpe( rp, gp, bp, rc, gc, bc, gamu ); - if(gamu==1){//gamut correction M.H.Brill S.Susstrunk - rp=MAXR(rp,0.0); - gp=MAXR(gp,0.0); - bp=MAXR(bp,0.0); - } - - rpa = nonlinear_adaptation( rp, fl ); - gpa = nonlinear_adaptation( gp, fl ); - bpa = nonlinear_adaptation( bp, fl ); - - return ((2.0 * rpa) + gpa + ((1.0 / 20.0) * bpa) - 0.305) * nbb; - } - - static float achromatic_response_to_whitefloat( float x, float y, float z, float d, float fl, float nbb, int gamu ) { - float r, g, b; - float rc, gc, bc; - float rp, gp, bp; - float rpa, gpa, bpa; - gamu=1; - xyz_to_cat02float( r, g, b, x, y, z, gamu ); - - rc = r * (((y * d) / r) + (1.0f - d)); - gc = g * (((y * d) / g) + (1.0f - d)); - bc = b * (((y * d) / b) + (1.0f - d)); - - cat02_to_hpefloat( rp, gp, bp, rc, gc, bc, gamu ); - if(gamu==1){//gamut correction M.H.Brill S.Susstrunk - rp=MAXR(rp,0.0f); - gp=MAXR(gp,0.0f); - bp=MAXR(bp,0.0f); - } - - rpa = nonlinear_adaptationfloat( rp, fl ); - gpa = nonlinear_adaptationfloat( gp, fl ); - bpa = nonlinear_adaptationfloat( bp, fl ); - - return ((2.0f * rpa) + gpa + ((1.0f / 20.0f) * bpa) - 0.305f) * nbb; - } - - static void xyz_to_cat02 ( double &r, double &g, double &b, double x, double y, double z, int gamu ); - static void cat02_to_hpe ( double &rh, double &gh, double &bh, double r, double g, double b, int gamu ); - static void cat02_to_xyz ( double &x, double &y, double &z, double r, double g, double b, int gamu ); - static void hpe_to_xyz ( double &x, double &y, double &z, double r, double g, double b ); - - static void xyz_to_cat02float ( float &r, float &g, float &b, float x, float y, float z, int gamu ); - static void cat02_to_hpefloat ( float &rh, float &gh, float &bh, float r, float g, float b, int gamu ); - static void cat02_to_xyzfloat ( float &x, float &y, float &z, float r, float g, float b, int gamu ); - static void hpe_to_xyzfloat ( float &x, float &y, float &z, float r, float g, float b ); - - static void Aab_to_rgb( double &r, double &g, double &b, double A, double aa, double bb, double nbb ); - static void Aab_to_rgbfloat( float &r, float &g, float &b, float A, float aa, float bb, float nbb ); - static void calculate_ab( double &aa, double &bb, double h, double e, double t, double nbb, double a ); - static void calculate_abfloat( float &aa, float &bb, float h, float e, float t, float nbb, float a ); - - - static double nonlinear_adaptation( double c, double fl ) { - double p; - if(c<0.0){ p = pow( (-1.0*fl * c) / 100.0, 0.42 );return ((-1.0*400.0 * p) / (27.13 + p)) + 0.1;} - else {p = pow( (fl * c) / 100.0, 0.42 ); return ((400.0 * p) / (27.13 + p)) + 0.1;} - } - static float nonlinear_adaptationfloat( float c, float fl ) { - float p; - if(c<0.0f){ p = pow_F( (-1.0f*fl * c) / 100.0f, 0.42f );return ((-1.0f*400.0f * p) / (27.13f + p)) + 0.1f;} - else {p = pow_F( (fl * c) / 100.0f, 0.42f ); return ((400.0f * p) / (27.13f + p)) + 0.1f;} - } - - static double inverse_nonlinear_adaptation( double c, double fl ) { - int c1; - if(c-0.1 < 0.0) c1=-1; else c1=1; - return c1*(100.0 / fl) * pow( (27.13 * fabs( c - 0.1 )) / (400.0 - fabs( c - 0.1 )), 1.0 / 0.42 ); - } - static float inverse_nonlinear_adaptationfloat( float c, float fl ) { - c -= 0.1f; - if(c < 0.f) { - fl *= -1.f; - if(c < -399.99f) // avoid nan values - c = -399.99f; - } else if(c > 399.99f) { // avoid nan values - c = 399.99f; - } - return (100.0f / fl) * pow_F( (27.13f * fabsf( c )) / (400.0f - fabsf( c )), 2.38095238f ); - } - - - static void curvecolor(double satind, double satval, double &sres, double parsat); - static void curvecolorfloat(float satind, float satval, float &sres, float parsat); - static void curveJ (double br, double contr, int db, LUTf & outCurve , LUTu & histogram ) ; - static void curveJfloat (float br, float contr, int db, LUTf & outCurve , LUTu & histogram ) ; - bool operator== (const ColorTemp& other) { return fabs(temp-other.temp)<1e-10 && fabs(green-other.green)<1e-10; } bool operator!= (const ColorTemp& other) { return !(*this==other); } @@ -308,54 +172,7 @@ class ColorTemp { static void spectrum_to_color_xyz_daylight (const double* spec_color, double _m1, double _m2, double &xx, double &yy, double &zz); static void spectrum_to_color_xyz_blackbody (const double* spec_color, double _temp, double &xx, double &yy, double &zz); static void spectrum_to_color_xyz_preset (const double* spec_color, const double* spec_intens, double &xx, double &yy, double &zz); - - - -/** - * Inverse transform from CIECAM02 JCh to XYZ. - */ - static void jch2xyz_ciecam02( double &x, double &y, double &z, - double J, double C, double h, - double xw, double yw, double zw, - double yb, double la, - double f, double c, double nc, int gamu, double n, double nbb, double ncb, double fl, double cz, double d, double aw); - - static void jch2xyz_ciecam02float( float &x, float &y, float &z, - float J, float C, float h, - float xw, float yw, float zw, - float yb, float la, - float f, float c, float nc,int gamu,float n, float nbb, float ncb, float fl, float cz, float d, float aw ); - -/** - * Forward transform from XYZ to CIECAM02 JCh. - */ - static void initcam1(double gamu, double yb, double pilotd, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, - double &cz, double &aw, double &wh, double &pfl, double &fl, double &c); - - static void initcam2(double gamu, double yb, double f, double la, double xw, double yw, double zw, double &n, double &d, double &nbb, double &ncb, - double &cz, double &aw, double &fl); - - static void initcam1float(float gamu, float yb, float pilotd, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, - float &cz, float &aw, float &wh, float &pfl, float &fl, float &c); - - static void initcam2float(float gamu, float yb, float f, float la, float xw, float yw, float zw, float &n, float &d, float &nbb, float &ncb, - float &cz, float &aw, float &fl); - - static void xyz2jchqms_ciecam02( double &J, double &C, double &h, - double &Q, double &M, double &s,double &aw, double &fl, double &wh, - double x, double y, double z, - double xw, double yw, double zw, - double yb, double la, - double f, double c, double nc, double pilotd,int gamu , double n, double nbb, double ncb, double pfl, double cz, double d ); - - static void xyz2jchqms_ciecam02float( float &J, float &C, float &h, - float &Q, float &M, float &s,float &aw, float &fl, float &wh, - float x, float y, float z, - float xw, float yw, float zw, - float yb, float la, - float f, float c, float nc, float pilotd, int gamu, float n, float nbb, float ncb, float pfl, float cz, float d ); - - + }; } #endif diff --git a/rtengine/dcrop.cc b/rtengine/dcrop.cc index 13fef46f5..c0ac11fbd 100644 --- a/rtengine/dcrop.cc +++ b/rtengine/dcrop.cc @@ -22,7 +22,6 @@ #include "mytime.h" #include "refreshmap.h" #include "rt_math.h" -#include "colortemp.h" // "ceil" rounding #define SKIPS(a,b) ((a) / (b) + ((a) % (b) > 0)) @@ -213,8 +212,8 @@ void Crop::update (int todo) { int CenterPreview_Y=trafy+(trafh*skip)/2; int minimuX=20000; int minimuY=20000; - int poscenterX; - int poscenterY; + int poscenterX=0; + int poscenterY=0; for(int cc=0;ccipf.dirpyrequalizer (labnCrop, skip); // parent->ipf.Lanczoslab (labnCrop,labnCrop , 1.f/skip); } - TMatrix wprof = iccStore->workingSpaceMatrix (params.icm.working); - double wp[3][3] = { - {wprof[0][0],wprof[0][1],wprof[0][2]}, - {wprof[1][0],wprof[1][1],wprof[1][2]}, - {wprof[2][0],wprof[2][1],wprof[2][2]}}; - TMatrix wiprof = iccStore->workingSpaceInverseMatrix (params.icm.working); - double wip[3][3] = { - {wiprof[0][0],wiprof[0][1],wiprof[0][2]}, - {wiprof[1][0],wiprof[1][1],wiprof[1][2]}, - {wiprof[2][0],wiprof[2][1],wiprof[2][2]} - }; + - - - - int kall=0; + int kall=0; int minwin=min(labnCrop->W,labnCrop->H); int maxlevelcrop=10; // if(cp.mul[9]!=0)maxlevelcrop=10; diff --git a/rtengine/iimage.h b/rtengine/iimage.h index 812b260dd..64a30f771 100644 --- a/rtengine/iimage.h +++ b/rtengine/iimage.h @@ -19,11 +19,8 @@ #ifndef _IIMAGE_ #define _IIMAGE_ -#include -#include #include #include -#include "../rtgui/threadutils.h" #include "rt_math.h" #include "alignedbuffer.h" #include "imagedimensions.h" @@ -31,7 +28,6 @@ #include "coord2d.h" #include "procparams.h" #include "color.h" -#include "colortemp.h" #define TR_NONE 0 #define TR_R90 1 diff --git a/rtengine/imageio.h b/rtengine/imageio.h index 1ccc6a5ee..980d86f92 100644 --- a/rtengine/imageio.h +++ b/rtengine/imageio.h @@ -35,7 +35,7 @@ #include "../rtexif/rtexif.h" #include "imagedimensions.h" #include "iimage.h" -#include "../rtgui/threadutils.h" +#include "colortemp.h" namespace rtengine { diff --git a/rtengine/improcfun.cc b/rtengine/improcfun.cc index 36c434b7e..27eca0f82 100644 --- a/rtengine/improcfun.cc +++ b/rtengine/improcfun.cc @@ -39,6 +39,7 @@ #include "EdgePreservingDecomposition.h" #include "improccoordinator.h" #include "clutstore.h" +#include "ciecam02.h" #ifdef _OPENMP #include @@ -413,14 +414,14 @@ if(params->colorappearance.enabled) { CAMBrightCurveJ(65536,0); CAMBrightCurveJ.dirty = false; } - ColorTemp::curveJ (jli, contra, 1, CAMBrightCurveJ, hist16J);//lightness and contrast J + Ciecam02::curveJ (jli, contra, 1, CAMBrightCurveJ, hist16J);//lightness and contrast J } if (needQ) { if (!CAMBrightCurveQ) { CAMBrightCurveQ(65536,0); CAMBrightCurveQ.dirty = false; } - ColorTemp::curveJ (qbri, qcontra, 1, CAMBrightCurveQ, hist16Q);//brightness and contrast Q + Ciecam02::curveJ (qbri, qcontra, 1, CAMBrightCurveQ, hist16Q);//brightness and contrast Q } } if(settings->viewinggreySc==0){//auto @@ -450,9 +451,9 @@ if(settings->viewinggreySc==1) yb=18.0; if(params->colorappearance.wbmodel=="RawT") {xw1=96.46;yw1=100.0;zw1=82.445;xw2=xwd;yw2=ywd;zw2=zwd;} //use RT WB; CAT 02 is used for output device (see prefreneces) else if(params->colorappearance.wbmodel=="RawTCAT02") {xw1=xw;yw1=yw;zw1=zw;xw2=xwd;yw2=ywd;zw2=zwd;} // Settings RT WB are used for CAT02 => mix , CAT02 is use for output device (screen: D50 D65, projector: lamp, LED) see preferences double cz,wh, pfl; - ColorTemp::initcam1(gamu, yb, pilot, f, la,xw, yw, zw, n, d, nbb, ncb,cz, aw, wh, pfl, fl, c); + Ciecam02::initcam1(gamu, yb, pilot, f, la,xw, yw, zw, n, d, nbb, ncb,cz, aw, wh, pfl, fl, c); double nj,dj,nbbj,ncbj,czj,awj,flj; - ColorTemp::initcam2(gamu,yb2, f2, la2, xw2, yw2, zw2, nj, dj, nbbj, ncbj,czj, awj, flj); + Ciecam02::initcam2(gamu,yb2, f2, la2, xw2, yw2, zw2, nj, dj, nbbj, ncbj,czj, awj, flj); @@ -497,7 +498,7 @@ if(settings->viewinggreySc==1) yb=18.0; y=(double)y1/655.35; z=(double)z1/655.35; //process source==> normal - ColorTemp::xyz2jchqms_ciecam02( J, C, h, + Ciecam02::xyz2jchqms_ciecam02( J, C, h, Q, M, s, aw, fl, wh, x, y, z, xw1, yw1, zw1, @@ -524,7 +525,7 @@ if(settings->viewinggreySc==1) yb=18.0; double parsat=1.5; parsat=1.5;//parsat=1.5 =>saturation ; 1.8 => chroma ; 2.5 => colorfullness (personal evaluation) if(schr==-100.0) schr=-99.8; - ColorTemp::curvecolor(schr, Sp , sres, parsat); + Ciecam02::curvecolor(schr, Sp , sres, parsat); double coe=pow(fl,0.25); float dred=100.f;// in C mode float protect_red=80.0f; // in C mode @@ -550,7 +551,7 @@ if(settings->viewinggreySc==1) yb=18.0; double parsat=2.5; if(mchr==-100.0) mchr=-99.8 ; if(mchr==100.0) mchr=99.9; - if(alg==3 || alg==2) ColorTemp::curvecolor(mchr, Mp , sres, parsat); else ColorTemp::curvecolor(0.0, Mp , sres, parsat);//colorfullness + if(alg==3 || alg==2) Ciecam02::curvecolor(mchr, Mp , sres, parsat); else Ciecam02::curvecolor(0.0, Mp , sres, parsat);//colorfullness float dred=100.f;//in C mode float protect_red=80.0f;// in C mode dred *=coe;//in M mode @@ -571,7 +572,7 @@ if(settings->viewinggreySc==1) yb=18.0; parsat=1.5; if(schr==-100.0) schr=-99.; if(schr==100.0) schr=98.; - if(alg==3) ColorTemp::curvecolor(schr, Sp , sres, parsat); else ColorTemp::curvecolor(0.0, Sp , sres, parsat); //saturation + if(alg==3) Ciecam02::curvecolor(schr, Sp , sres, parsat); else Ciecam02::curvecolor(0.0, Sp , sres, parsat); //saturation dred=100.f;// in C mode protect_red=80.0f; // in C mode dred = 100.0 * sqrt((dred*coe)/Q); @@ -584,7 +585,7 @@ if(settings->viewinggreySc==1) yb=18.0; Cp=Cpro/100.0; parsat=1.8;//parsat=1.5 =>saturation ; 1.8 => chroma ; 2.5 => colorfullness (personal evaluation : for not) if(chr==-100.0) chr=-99.8; - if(alg!=2) ColorTemp::curvecolor(chr, Cp , sres, parsat);else ColorTemp::curvecolor(0.0, Cp , sres, parsat); //chroma + if(alg!=2) Ciecam02::curvecolor(chr, Cp , sres, parsat);else Ciecam02::curvecolor(0.0, Cp , sres, parsat); //chroma dred=55.f; protect_red=30.0f; sk=1; @@ -859,7 +860,7 @@ if(settings->viewinggreySc==1) yb=18.0; double xx,yy,zz; //double nj, nbbj, ncbj, flj, czj, dj, awj; //process normal==> viewing - ColorTemp::jch2xyz_ciecam02( xx, yy, zz, + Ciecam02::jch2xyz_ciecam02( xx, yy, zz, J, C, h, xw2, yw2, zw2, yb2, la2, @@ -1066,7 +1067,7 @@ if((params->colorappearance.tonecie || (params->colorappearance.tonecie && param } //end histograms // double nd, nbbd, ncbd, fld, czd, dd, awd; - ColorTemp::jch2xyz_ciecam02( xx, yy, zz, + Ciecam02::jch2xyz_ciecam02( xx, yy, zz, ncie->J_p[i][j], ncie->C_p[i][j], ncie->h_p[i][j], xw2, yw2, zw2, yb2, la2, @@ -1367,22 +1368,22 @@ if(params->colorappearance.enabled) { //evaluate lightness, contrast if (needJ) { if (!CAMBrightCurveJ) { - CAMBrightCurveJ(32768,0); + CAMBrightCurveJ(32768,LUT_CLIP_ABOVE); CAMBrightCurveJ.dirty = false; } float jli=params->colorappearance.jlight; float contra=params->colorappearance.contrast; - ColorTemp::curveJfloat (jli, contra, 1, CAMBrightCurveJ, hist16J);//lightness and contrast J + Ciecam02::curveJfloat (jli, contra, 1, CAMBrightCurveJ, hist16J);//lightness and contrast J } if (needQ) { if (!CAMBrightCurveQ) { - CAMBrightCurveQ(32768,0); + CAMBrightCurveQ(32768,LUT_CLIP_ABOVE); CAMBrightCurveQ.clear(); CAMBrightCurveQ.dirty = false; } float qbri=params->colorappearance.qbright; float qcontra=params->colorappearance.qcontrast; - ColorTemp::curveJfloat (qbri, qcontra, 1, CAMBrightCurveQ, hist16Q);//brightness and contrast Q + Ciecam02::curveJfloat (qbri, qcontra, 1, CAMBrightCurveQ, hist16Q);//brightness and contrast Q } } @@ -1413,10 +1414,10 @@ if(settings->viewinggreySc==1) yb=18.0f;//fixed if(params->colorappearance.wbmodel=="RawT") {xw1=96.46f;yw1=100.0f;zw1=82.445f;xw2=xwd;yw2=ywd;zw2=zwd;} //use RT WB; CAT 02 is used for output device (see prefreneces) else if(params->colorappearance.wbmodel=="RawTCAT02") {xw1=xw;yw1=yw;zw1=zw;xw2=xwd;yw2=ywd;zw2=zwd;} // Settings RT WB are used for CAT02 => mix , CAT02 is use for output device (screen: D50 D65, projector: lamp, LED) see preferences float cz,wh, pfl; - ColorTemp::initcam1float(gamu, yb, pilot, f, la,xw, yw, zw, n, d, nbb, ncb,cz, aw, wh, pfl, fl, c); + Ciecam02::initcam1float(gamu, yb, pilot, f, la,xw, yw, zw, n, d, nbb, ncb,cz, aw, wh, pfl, fl, c); const float pow1 = pow_F( 1.64f - pow_F( 0.29f, n ), 0.73f ); float nj,dj,nbbj,ncbj,czj,awj,flj; - ColorTemp::initcam2float(gamu,yb2, f2, la2, xw2, yw2, zw2, nj, dj, nbbj, ncbj,czj, awj, flj); + Ciecam02::initcam2float(gamu,yb2, f2, la2, xw2, yw2, zw2, nj, dj, nbbj, ncbj,czj, awj, flj); const float pow1n = pow_F( 1.64f - pow_F( 0.29f, nj ), 0.73f ); const float epsil=0.0001f; @@ -1466,7 +1467,7 @@ if(settings->viewinggreySc==1) yb=18.0f;//fixed y=(float)y1/655.35f; z=(float)z1/655.35f; //process source==> normal - ColorTemp::xyz2jchqms_ciecam02float( J, C, h, + Ciecam02::xyz2jchqms_ciecam02float( J, C, h, Q, M, s, aw, fl, wh, x, y, z, xw1, yw1, zw1, @@ -1479,15 +1480,22 @@ if(settings->viewinggreySc==1) yb=18.0f;//fixed Mpro=M; spro=s; // we cannot have all algorithms with all chroma curves + if(alg==0) { + Jpro = CAMBrightCurveJ[Jpro*327.68f]/327.68f; //lightness CIECAM02 + contrast + Qpro = QproFactor * sqrtf(Jpro); + float Cp = (spro*spro*Qpro)/(1000000.f); + Cpro = Cp * 100.f; + float sres; + Ciecam02::curvecolorfloat(chr, Cp , sres, 1.8f); + Color::skinredfloat(Jpro, hpro, sres, Cp, 55.f, 30.f,1,rstprotection, 100.f, Cpro); + } if(alg==1) { // Lightness saturation - if(Jpro > 99.9f) - Jpro = 99.9f; - Jpro=(CAMBrightCurveJ[(float)(Jpro*327.68f)])/327.68f;//lightness CIECAM02 + contrast + Jpro = CAMBrightCurveJ[Jpro*327.68f]/327.68f;//lightness CIECAM02 + contrast float sres; float Sp=spro/100.0f; - float parsat=1.5f; //parsat=1.5 =>saturation ; 1.8 => chroma ; 2.5 => colorfulness (personal evaluation) - ColorTemp::curvecolorfloat(schr, Sp , sres, parsat); + float parsat=1.5f; //parsat=1.5 =>saturation ; 1.8 => chroma ; 2.5 => colorfullness (personal evaluation) + Ciecam02::curvecolorfloat(schr, Sp , sres, parsat); float dred=100.f;// in C mode float protect_red=80.0f; // in C mode dred = 100.0f * sqrtf((dred*coe)/Qpro); @@ -1495,18 +1503,13 @@ if(settings->viewinggreySc==1) yb=18.0f;//fixed Color::skinredfloat(Jpro, hpro, sres, Sp, dred, protect_red,0,rstprotection,100.f, spro); Qpro = QproFactor * sqrtf(Jpro); Cpro=(spro*spro*Qpro)/(10000.0f); - } - else if(alg==3 || alg==0 || alg==2) { - if(alg==3 || alg==2) { - float coef=32760.f/wh; - if(Qpro*coef >= 32767.0f) - Qpro=(CAMBrightCurveQ[32767])/coef;//brightness and contrast - else - Qpro=(CAMBrightCurveQ[(float)(Qpro*coef)])/coef;//brightness and contrast - } + } + else if(alg==2) { + float coef=32767.f/wh; + Qpro=(CAMBrightCurveQ[(float)(Qpro*coef)])/coef;//brightness and contrast float Mp, sres; Mp=Mpro/100.0f; - ColorTemp::curvecolorfloat(mchr, Mp , sres, 2.5f); + Ciecam02::curvecolorfloat(mchr, Mp , sres, 2.5f); float dred=100.f;//in C mode float protect_red=80.0f;// in C mode dred *=coe;//in M mode @@ -1516,13 +1519,30 @@ if(settings->viewinggreySc==1) yb=18.0f;//fixed Cpro= Mpro/coe; Qpro = (Qpro == 0.f ? epsil : Qpro); // avoid division by zero spro = 100.0f * sqrtf( Mpro / Qpro ); - if(alg!=2) { - if(Jpro > 99.9f) - Jpro = 99.9f; - Jpro=(CAMBrightCurveJ[(float)(Jpro*327.68f)])/327.68f;//lightness CIECAM02 + contrast - } + } + else if(alg==3) { + float coef=32760.f/wh; + if(Qpro*coef >= 32767.0f) + Qpro=(CAMBrightCurveQ[32767])/coef;//brightness and contrast + else + Qpro=(CAMBrightCurveQ[(float)(Qpro*coef)])/coef;//brightness and contrast + float Mp, sres; + Mp=Mpro/100.0f; + Ciecam02::curvecolorfloat(mchr, Mp , sres, 2.5f); + float dred=100.f;//in C mode + float protect_red=80.0f;// in C mode + dred *=coe;//in M mode + protect_red *=coe;//M mode + Color::skinredfloat(Jpro, hpro, sres, Mp, dred, protect_red,0,rstprotection,100.f, Mpro); + Jpro = SQR((10.f*Qpro)/wh); + Cpro= Mpro/coe; + Qpro = (Qpro == 0.f ? epsil : Qpro); // avoid division by zero + spro = 100.0f * sqrtf( Mpro / Qpro ); + if(Jpro > 99.9f) + Jpro = 99.9f; + Jpro=(CAMBrightCurveJ[(float)(Jpro*327.68f)])/327.68f;//lightness CIECAM02 + contrast float Sp=spro/100.0f; - ColorTemp::curvecolorfloat(schr, Sp , sres, 1.5f); + Ciecam02::curvecolorfloat(schr, Sp , sres, 1.5f); dred=100.f;// in C mode protect_red=80.0f; // in C mode dred = 100.0f * sqrtf((dred*coe)/Q); @@ -1531,19 +1551,17 @@ if(settings->viewinggreySc==1) yb=18.0f;//fixed Qpro = QproFactor * sqrtf(Jpro); float Cp=(spro*spro*Qpro)/(1000000.f); Cpro = Cp*100.f; - ColorTemp::curvecolorfloat(chr, Cp , sres, 1.8f); + Ciecam02::curvecolorfloat(chr, Cp , sres, 1.8f); Color::skinredfloat(Jpro, hpro, sres, Cp, 55.f, 30.f,1,rstprotection, 100.f, Cpro); // disabled this code, Issue 2690 // if(Jpro < 1.f && Cpro > 12.f) Cpro=12.f;//reduce artifacts by "pseudo gamut control CIECAM" // else if(Jpro < 2.f && Cpro > 15.f) Cpro=15.f; // else if(Jpro < 4.f && Cpro > 30.f) Cpro=30.f; // else if(Jpro < 7.f && Cpro > 50.f) Cpro=50.f; - if(alg==3) { - hpro=hpro+hue; - if( hpro < 0.0f ) - hpro += 360.0f;//hue + hpro=hpro+hue; + if( hpro < 0.0f ) + hpro += 360.0f;//hue } - } if (hasColCurve1) {//curve 1 with Lightness and Brightness if (curveMode==ColorAppearanceParams::TC_MODE_LIGHT){ @@ -1782,7 +1800,7 @@ if(LabPassOne){ float xx,yy,zz; //process normal==> viewing - ColorTemp::jch2xyz_ciecam02float( xx, yy, zz, + Ciecam02::jch2xyz_ciecam02float( xx, yy, zz, J, C, h, xw2, yw2, zw2, yb2, la2, @@ -2048,7 +2066,7 @@ if((params->colorappearance.tonecie && (epdEnabled)) || (params->sharpening.enab } //end histograms - ColorTemp::jch2xyz_ciecam02float( xx, yy, zz, + Ciecam02::jch2xyz_ciecam02float( xx, yy, zz, ncie->J_p[i][j], ncie_C_p, ncie->h_p[i][j], xw2, yw2, zw2, yb2, la2, diff --git a/rtengine/rawimage.cc b/rtengine/rawimage.cc index 6931f32b6..563a0c980 100755 --- a/rtengine/rawimage.cc +++ b/rtengine/rawimage.cc @@ -6,7 +6,6 @@ #include "rawimage.h" #include "settings.h" -#include "colortemp.h" #include "camconst.h" #include "utils.h" #ifdef WIN32 diff --git a/rtengine/rawimagesource.h b/rtengine/rawimagesource.h index cfebfc1cb..f711509bb 100644 --- a/rtengine/rawimagesource.h +++ b/rtengine/rawimagesource.h @@ -20,15 +20,11 @@ #define _RAWIMAGESOURCE_ #include "imagesource.h" -#include #include "dcp.h" #include "array2D.h" #include "curves.h" -#include #include "color.h" #include "iimage.h" -#include "../rtgui/cacheimagedata.h" -#include "../rtgui/threadutils.h" #define HR_SCALE 2 @@ -223,8 +219,6 @@ class RawImageSource : public ImageSource { int findHotDeadPixels( PixelsMap &bpMap, float thresh, bool findHotPixels, bool findDeadPixels ); void cfa_linedn (float linenoiselevel);//Emil's line denoise - void cfa_tile_denoise (fftwf_complex * fcfablox, int vblk, int hblk, int numblox_H, int numblox_W, float noisevar, float * rolloff ); - void cfa_output_tile_row (float *cfabloxrow, float ** cfahipassdn, float ** tilemask_out, int height, int width, int top, int blkrad); void green_equilibrate (float greenthresh);//Emil's green equilibration diff --git a/rtexif/rtexif.cc b/rtexif/rtexif.cc index 46ed2317f..9701b9022 100644 --- a/rtexif/rtexif.cc +++ b/rtexif/rtexif.cc @@ -229,7 +229,7 @@ bool TagDirectory::CPBDump (const Glib::ustring &commFName, const Glib::ustring std::vector tagDirList; std::vector tagDirPaths; - FILE *f; + FILE *f = NULL; if (!keyFile) { // open the file in write mode f = safe_g_fopen (commFName, "wt"); @@ -658,8 +658,9 @@ Tag::Tag (TagDirectory* p, FILE* f, int base) if( tag == 0xc634 ){ // DNGPrivateData int currPos = ftell(f); - char buffer[32],*p=buffer; - while( fread (p, 1, 1, f ) && *p != 0 && p-buffer126) isstring = false; @@ -1709,7 +1710,7 @@ parse_leafdata(TagDirectory* root, ByteOrder order) { int iso_speed = 0; int rotation_angle = 0; int found_count = 0; - while (pos + sizeof(hdr) <= valuesize && found_count < 2) { + while (pos + (int)sizeof(hdr) <= valuesize && found_count < 2) { hdr = (char *)&value[pos]; if (strncmp(hdr, PKTS_tag, 4) != 0) { // in a few cases the header can be offset a few bytes, don't know why diff --git a/rtexif/rtexif.h b/rtexif/rtexif.h index d89d6d5ea..f5f238640 100644 --- a/rtexif/rtexif.h +++ b/rtexif/rtexif.h @@ -27,6 +27,7 @@ #include #include #include +#include #include "../rtengine/procparams.h" #include "../rtengine/safekeyfile.h" diff --git a/rtgui/cacheimagedata.h b/rtgui/cacheimagedata.h index d08ee8cf5..a04bd7300 100644 --- a/rtgui/cacheimagedata.h +++ b/rtgui/cacheimagedata.h @@ -69,7 +69,7 @@ class CacheImageData { enum { FULL_THUMBNAIL = 0, // was the thumbnail generated from whole file - QUICK_THUMBNAIL = 1 // was rhe thumbnail generated from embedded jpeg + QUICK_THUMBNAIL = 1 // was the thumbnail generated from embedded jpeg }; CacheImageData (); diff --git a/rtgui/options.h b/rtgui/options.h index 3f1b3e0c2..cde759862 100644 --- a/rtgui/options.h +++ b/rtgui/options.h @@ -19,7 +19,7 @@ #ifndef _OPTIONS_ #define _OPTIONS_ -#include +#include #include "../rtengine/rtengine.h" #define STARTUPDIR_CURRENT 0