/* * 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 #include "alignedbuffer.h" #include "color.h" #include "iccmatrices.h" #include "iccstore.h" #include "image8.h" #include "imagefloat.h" #include "improcfun.h" #include "labimage.h" #include "procparams.h" #include "rtengine.h" #include "settings.h" #include "utils.h" namespace rtengine { namespace { inline void copyAndClampLine(const float *src, unsigned char *dst, const int W) { for (int j = 0; j < W * 3; ++j) { dst[j] = uint16ToUint8Rounded(CLIP(src[j] * MAXVALF)); } } inline void copyAndClamp(const LabImage *src, unsigned char *dst, const double rgb_xyz[3][3], bool multiThread, int pro)//int pro to switch to Working Profile histogram mode without gamma { const int W = src->W; const int H = src->H; float rgb_xyzf[3][3]; for (int i = 0; i < 3; i++) { for (int j = 0; j < 3; j++) { rgb_xyzf[i][j] = rgb_xyz[i][j]; } } #ifdef __SSE2__ vfloat rgb_xyzv[3][3]; for (int i = 0; i < 3; i++) { for (int j = 0; j < 3; j++) { rgb_xyzv[i][j] = F2V(rgb_xyzf[i][j]); } } #endif #ifdef _OPENMP #pragma omp parallel for schedule(dynamic,16) if (multiThread) #endif for (int i = 0; i < H; ++i) { float* rL = src->L[i]; float* ra = src->a[i]; float* rb = src->b[i]; int ix = i * 3 * W; #ifdef __SSE2__ float rbuffer[W] ALIGNED16; float gbuffer[W] ALIGNED16; float bbuffer[W] ALIGNED16; int j = 0; for (; j < W - 3; j += 4) { vfloat R, G, B; vfloat x_, y_, z_; Color::Lab2XYZ(LVFU(rL[j]), LVFU(ra[j]), LVFU(rb[j]), x_, y_, z_); Color::xyz2rgb(x_, y_, z_, R, G, B, rgb_xyzv); if(pro == 0) { STVF(rbuffer[j], Color::gamma2curve[R]); STVF(gbuffer[j], Color::gamma2curve[G]); STVF(bbuffer[j], Color::gamma2curve[B]); } else {//Working profile and gamma=1 STVF(rbuffer[j], R); STVF(gbuffer[j], G); STVF(bbuffer[j], B); } } for (; j < W; ++j) { float R, G, B; float x_, y_, z_; Color::Lab2XYZ(rL[j], ra[j], rb[j], x_, y_, z_); Color::xyz2rgb(x_, y_, z_, R, G, B, rgb_xyzf); if(pro == 0) { rbuffer[j] = Color::gamma2curve[R]; gbuffer[j] = Color::gamma2curve[G]; bbuffer[j] = Color::gamma2curve[B]; } else {//Working profile and gamma=1 rbuffer[j] = R; gbuffer[j] = G; bbuffer[j] = B; } } for (j = 0; j < W; ++j) { dst[ix++] = uint16ToUint8Rounded(rbuffer[j]); dst[ix++] = uint16ToUint8Rounded(gbuffer[j]); dst[ix++] = uint16ToUint8Rounded(bbuffer[j]); } #else for (int j = 0; j < W; ++j) { float R, G, B; float x_, y_, z_; Color::Lab2XYZ(rL[j], ra[j], rb[j], x_, y_, z_); Color::xyz2rgb(x_, y_, z_, R, G, B, rgb_xyzf); if(pro == 0) { dst[ix++] = uint16ToUint8Rounded(Color::gamma2curve[R]); dst[ix++] = uint16ToUint8Rounded(Color::gamma2curve[G]); dst[ix++] = uint16ToUint8Rounded(Color::gamma2curve[B]); } else {//Working profile and gamma=1 dst[ix++] = uint16ToUint8Rounded(R); dst[ix++] = uint16ToUint8Rounded(G); dst[ix++] = uint16ToUint8Rounded(B); } } #endif } } } // namespace float gammalog(float x, float p, float s, float g3, float g4) { return x <= g3 ? x * s : (1.f + g4) * xexpf(xlogf(x) / p) - g4;//continuous } #ifdef __SSE2__ vfloat gammalog(vfloat x, vfloat p, vfloat s, vfloat g3, vfloat g4) { return vself(vmaskf_le(x, g3), x * s, (F2V(1.f) + g4) * xexpf(xlogf(x) / p) - g4);//continuous } #endif // Used in ImProcCoordinator::updatePreviewImage (rtengine/improccoordinator.cc) // Crop::update (rtengine/dcrop.cc) // Thumbnail::processImage (rtengine/rtthumbnail.cc) // // If monitorTransform, divide by 327.68 then apply monitorTransform (which can integrate soft-proofing) // otherwise divide by 327.68, convert to xyz and apply the sRGB transform, before converting with gamma2curve void ImProcFunctions::lab2monitorRgb(LabImage* lab, Image8* image) { if (monitorTransform) { const int W = lab->W; const int H = lab->H; unsigned char * data = image->data; // cmsDoTransform is relatively expensive #ifdef _OPENMP #pragma omp parallel firstprivate(lab, data, W, H) #endif { AlignedBuffer pBuf(3 * lab->W); AlignedBuffer mBuf; AlignedBuffer gwBuf1; AlignedBuffer gwBuf2; if (gamutWarning) { gwBuf1.resize(3 * lab->W); gwBuf2.resize(3 * lab->W); mBuf.resize(3 * lab->W); } float *buffer = pBuf.data; float *outbuffer = gamutWarning ? mBuf.data : pBuf.data; // make in place transformations when gamutWarning is not needed #ifdef _OPENMP #pragma omp for schedule(dynamic,16) #endif for (int i = 0; i < H; i++) { const int ix = i * 3 * W; int iy = 0; float* rL = lab->L[i]; float* ra = lab->a[i]; float* rb = lab->b[i]; for (int j = 0; j < W; j++) { buffer[iy++] = rL[j] / 327.68f; buffer[iy++] = ra[j] / 327.68f; buffer[iy++] = rb[j] / 327.68f; } cmsDoTransform(monitorTransform, buffer, outbuffer, W); copyAndClampLine(outbuffer, data + ix, W); if (gamutWarning) { gamutWarning->markLine(image, i, buffer, gwBuf1.data, gwBuf2.data); } } } // End of parallelization } else { copyAndClamp(lab, image->data, sRGB_xyz, multiThread, 0);//int pro = 0 always sent 'normal' to monitor } } // Used in ImProcCoordinator::updatePreviewImage (rtengine/improccoordinator.cc) // Crop::update (rtengine/dcrop.cc) // // Generate an Image8 // // If output profile used, divide by 327.68 then apply the "profile" profile (eventually with a standard gamma) // otherwise divide by 327.68, convert to xyz and apply the RGB transform, before converting with gamma2curve Image8* ImProcFunctions::lab2rgb(LabImage* lab, int cx, int cy, int cw, int ch, const procparams::ColorManagementParams &icm, bool consider_histogram_settings) { if (cx < 0) { cx = 0; } if (cy < 0) { cy = 0; } if (cx + cw > lab->W) { cw = lab->W - cx; } if (cy + ch > lab->H) { ch = lab->H - cy; } Image8* image = new Image8(cw, ch); Glib::ustring profile; cmsHPROFILE oprof = nullptr; int pro = 0;//int pro to switch to Working Profile histogram mode without gamma if (settings->HistogramWorking && consider_histogram_settings) { profile = icm.workingProfile; pro = 1;//no gamma for histogram and navigator panel and Lockable color picker } else { profile = icm.outputProfile; //pro = 0 histogram and navigator panel and Lockable color picker with gamma if (icm.outputProfile.empty() || icm.outputProfile == ColorManagementParams::NoICMString) { profile = "sRGB"; } oprof = ICCStore::getInstance()->getProfile(profile); } if (oprof) { const cmsUInt32Number flags = cmsFLAGS_NOOPTIMIZE | cmsFLAGS_NOCACHE | (icm.outputBPC ? cmsFLAGS_BLACKPOINTCOMPENSATION : 0); // NOCACHE is important for thread safety lcmsMutex->lock(); cmsHPROFILE LabIProf = cmsCreateLab4Profile(nullptr); cmsHTRANSFORM hTransform = cmsCreateTransform(LabIProf, TYPE_Lab_DBL, oprof, TYPE_RGB_FLT, icm.outputIntent, flags); cmsCloseProfile(LabIProf); lcmsMutex->unlock(); unsigned char *data = image->data; // cmsDoTransform is relatively expensive #ifdef _OPENMP #pragma omp parallel #endif { AlignedBuffer pBuf(3 * cw); AlignedBuffer oBuf(3 * cw); double *buffer = pBuf.data; float *outbuffer = oBuf.data; int condition = cy + ch; #ifdef _OPENMP #pragma omp for firstprivate(lab) schedule(dynamic,16) #endif for (int i = cy; i < condition; i++) { const int ix = i * 3 * cw; int iy = 0; float* rL = lab->L[i]; float* ra = lab->a[i]; float* rb = lab->b[i]; for (int j = cx; j < cx + cw; j++) { buffer[iy++] = rL[j] / 327.68f; buffer[iy++] = ra[j] / 327.68f; buffer[iy++] = rb[j] / 327.68f; } cmsDoTransform(hTransform, buffer, outbuffer, cw); copyAndClampLine(outbuffer, data + ix, cw); } } // End of parallelization cmsDeleteTransform(hTransform); } else { const auto xyz_rgb = ICCStore::getInstance()->workingSpaceInverseMatrix(profile); copyAndClamp(lab, image->data, xyz_rgb, multiThread, pro);//int pro = 1 to switch to Working Profile for histogram and navigator panel and Lockable color picker whitout gamma } return image; } /** @brief Convert the final Lab image to the output RGB color space * * Used in processImage (rtengine/simpleprocess.cc) * * Provide a pointer to a 7 floats array for "ga" (uninitialized ; this array will be filled with the gamma values) if you want * to use the custom gamma scenario. Those gamma values will correspond to the ones of the chosen standard output profile * (Prophoto if non standard output profile given) * * If "ga" is NULL, then we're considering standard gamma with the chosen output profile. * * Generate an Image16 * * If a custom gamma profile can be created, divide by 327.68, convert to xyz and apply the custom gamma transform * otherwise divide by 327.68, convert to xyz and apply the sRGB transform, before converting with gamma2curve */ Imagefloat* ImProcFunctions::lab2rgbOut(LabImage* lab, int cx, int cy, int cw, int ch, const procparams::ColorManagementParams &icm) { if (cx < 0) { cx = 0; } if (cy < 0) { cy = 0; } if (cx + cw > lab->W) { cw = lab->W - cx; } if (cy + ch > lab->H) { ch = lab->H - cy; } Imagefloat* image = new Imagefloat(cw, ch); cmsHPROFILE oprof = ICCStore::getInstance()->getProfile(icm.outputProfile); if (oprof) { cmsUInt32Number flags = cmsFLAGS_NOOPTIMIZE | cmsFLAGS_NOCACHE; if (icm.outputBPC) { flags |= cmsFLAGS_BLACKPOINTCOMPENSATION; } lcmsMutex->lock(); cmsHPROFILE iprof = cmsCreateLab4Profile(nullptr); cmsHTRANSFORM hTransform = cmsCreateTransform(iprof, TYPE_Lab_FLT, oprof, TYPE_RGB_FLT, icm.outputIntent, flags); lcmsMutex->unlock(); image->ExecCMSTransform(hTransform, *lab, cx, cy); cmsDeleteTransform(hTransform); cmsCloseProfile(iprof); image->normalizeFloatTo65535(); } else { #ifdef _OPENMP #pragma omp parallel for schedule(dynamic,16) if (multiThread) #endif for (int i = cy; i < cy + ch; i++) { float R, G, B; float* rL = lab->L[i]; float* ra = lab->a[i]; float* rb = lab->b[i]; for (int j = cx; j < cx + cw; j++) { float fy = (Color::c1By116 * rL[j]) / 327.68f + Color::c16By116; // (L+16)/116 float fx = (0.002f * ra[j]) / 327.68f + fy; float fz = fy - (0.005f * rb[j]) / 327.68f; float LL = rL[j] / 327.68f; float x_ = 65535.0f * Color::f2xyz(fx) * Color::D50x; //float y_ = 65535.0 * Color::f2xyz(fy); float z_ = 65535.0f * Color::f2xyz(fz) * Color::D50z; float y_ = (LL > (float)Color::epskap) ? 65535.0f * fy * fy * fy : 65535.0f * LL / (float)Color::kappa; Color::xyz2srgb(x_, y_, z_, R, G, B); image->r(i - cy, j - cx) = Color::gamma2curve[CLIP(R)]; image->g(i - cy, j - cx) = Color::gamma2curve[CLIP(G)]; image->b(i - cy, j - cx) = Color::gamma2curve[CLIP(B)]; } } } return image; } void ImProcFunctions::preserv(LabImage *nprevl, LabImage *provis, int cw, int ch) { //avoid too strong in middle values chroma when changing primaries float pres = 0.01f * params->icm.preser; float neutral = 2000000000.f;//if a2 + b2 < 200000000 scale 0..100 a and b about : 140 > a & b > -140 decrease effect float medneutral = 10000000.f;//plein effect 10 > a & b > -10 float aaneu = 1.f / (medneutral - neutral); float bbneu = - aaneu * neutral; #ifdef _OPENMP #pragma omp for schedule(dynamic, 16) nowait #endif for (int i = 0; i < ch; ++i) for (int j = 0; j < cw; ++j) { float neu = SQR(provis->a[i][j]) + SQR(provis->b[i][j]); if (neu < medneutral) {//plein effect nprevl->a[i][j] = intp(pres, provis->a[i][j], nprevl->a[i][j]); nprevl->b[i][j] = intp(pres, provis->b[i][j], nprevl->b[i][j]); } else if (neu < neutral) {//decrease effect float presred = aaneu * neu + bbneu; nprevl->a[i][j] = intp(pres * presred, provis->a[i][j], nprevl->a[i][j]); nprevl->b[i][j] = intp(pres * presred, provis->b[i][j], nprevl->b[i][j]); } } } // ACES-style gamut compression // // tweaked from the original from https://github.com/jedypod/gamut-compress // tweaked from CTL in ART thanks to Alberto Griggio //from ACES https://docs.acescentral.com/specifications/rgc/#appendix-c-illustrations // https://docs.acescentral.com/specifications/rgc/#appendix-d-ctl-reference-implementation // https://docs.acescentral.com/specifications/rgc/ // Distance from achromatic which will be compressed to the gamut boundary // Values calculated to encompass the encoding gamuts of common digital cinema cameras //const float LIM_CYAN = 1.147; //const float LIM_MAGENTA = 1.264; //const float LIM_YELLOW = 1.312; //Percentage of the core gamut to protect // Values calculated to protect all the colors of the ColorChecker Classic 24 as given by // ISO 17321-1 and Ohta (1997) //const float THR_CYAN = 0.815; //const float THR_MAGENTA = 0.803; //const float THR_YELLOW = 0.880; // Aggressiveness of the compression curve //const float PWR = 1.2; void ImProcFunctions::gamutcompr( Imagefloat *src, Imagefloat *dst) const { if (settings->verbose) { printf("Apply compression gamut \n"); } using Triple = std::array; using Matrix = std::array; const TMatrix wprof = ICCStore::getInstance()->workingSpaceMatrix(params->icm.workingProfile); Matrix wpro = {}; //working profile set in Matrix format for (int r = 0; r < 3; ++r) { for (int c = 0; c < 3; ++c) { wpro[r][c] = wprof[r][c]; } } //dcip3 Rec2020, srgb, prophoto, acesp1 - Compression gamut matrix profile Matrix dcip3 = {}; //in fact the white point is "special" - 0.314 - 0.351 Theater dcip3[0][0] = 0.4861607;//0.4451698;(original) //0.4861607 with chromatic adaptation D63 => D50 dcip3[0][1] = 0.3238514;//0.2771344;(original)//0.3238514 with chromatic adaptation D63 => D50 dcip3[0][2] = 0.1541879;//0.1722827;(original//0.1541879 with chromatic adaptation D63 => D50 dcip3[1][0] = 0.2266839;//0.2094917;(original//0.2266839 with chromatic adaptation D63 => D50 dcip3[1][1] = 0.7103336;//0.7215953;(original//0.7103336 with chromatic adaptation D63 => D50 dcip3[1][2] = 0.0629826;//0.0689131;(original// 0.0629826 with chromatic adaptation D63 => D50 dcip3[2][0] = -0.0008016;//0.0;(original//-0.0008016 with chromatic adaptation D63 => D50 dcip3[2][1] = 0.0432353;//0.0470606;(original// 0.0432353 with chromatic adaptation D63 => D50 dcip3[2][2] = 0.7824663;//0.9073554;(original// 0.7824663 with chromatic adaptation D63 => D50 Matrix Rec2020 = {}; Rec2020[0][0] = 0.6734241; Rec2020[0][1] = 0.1656411; Rec2020[0][2] = 0.1251286; Rec2020[1][0] = 0.2790177; Rec2020[1][1] = 0.6753402; Rec2020[1][2] = 0.0456377; Rec2020[2][0] = -0.0019300; Rec2020[2][1] = 0.0299784; Rec2020[2][2] = 0.7973330; Matrix srgb = {}; srgb[0][0] = 0.4360747; srgb[0][1] = 0.3850649; srgb[0][2] = 0.1430804; srgb[1][0] = 0.2225045; srgb[1][1] = 0.7168786; srgb[1][2] = 0.0606169; srgb[2][0] = 0.0139322; srgb[2][1] = 0.0971045; srgb[2][2] = 0.7141733; Matrix adobe = {}; adobe[0][0] = 0.6097559; adobe[0][1] = 0.2052401; adobe[0][2] = 0.1492240; adobe[1][0] = 0.3111242; adobe[1][1] = 0.6256560; adobe[1][2] = 0.0632197; adobe[2][0] = 0.0194811; adobe[2][1] = 0.0608902; adobe[2][2] = 0.7448387; Matrix prophoto = {};//prophoto prophoto[0][0] = 0.7976749; prophoto[0][1] = 0.1351917; prophoto[0][2] = 0.0313534; prophoto[1][0] = 0.2880402; prophoto[1][1] = 0.7118741; prophoto[1][2] = 0.0000857; prophoto[2][0] = 0.0; prophoto[2][1] = 0.0; prophoto[2][2] = 1.2118128; Matrix acesp1 = {};//aces P1 acesp1[0][0] = 0.689697; acesp1[0][1] = 0.149944; acesp1[0][2] = 0.124559; acesp1[1][0] = 0.284448; acesp1[1][1] = 0.671758; acesp1[1][2] = 0.043794; acesp1[2][0] = -0.006043; acesp1[2][1] = 0.009998; acesp1[2][2] = 0.820945; Matrix out = {}; if (params->cg.colorspace == "rec2020") { out = Rec2020; } else if (params->cg.colorspace == "prophoto") { out = prophoto; } else if (params->cg.colorspace == "adobe") { out = adobe; } else if (params->cg.colorspace == "srgb") { out = srgb; } else if (params->cg.colorspace == "dcip3") { out = dcip3; } else if (params->cg.colorspace == "acesp1") { out = acesp1; } else { out = acesp1; // Should never happen, but just in case. } Matrix inv_out = {}; if (!rtengine::invertMatrix(out, inv_out)) {//invert matrix printf("Matrix is not invertible, skipping\n"); } Matrix Rprov = {}; Color::multip(inv_out, wpro, Rprov);//multiply matrix Matrix to_out = {}; Color::transpose(Rprov, to_out);//transpose Matrix for output Matrix from_out = {};//inverse to output if (!rtengine::invertMatrix(to_out, from_out)) { printf("Matrix is not invertible, skipping\n"); } //parameters from GUI const auto thc = static_cast(params->cg.th_c); const auto thm = static_cast(params->cg.th_m); const auto thy = static_cast(params->cg.th_y); const auto dc = static_cast(params->cg.d_c); const auto dm = static_cast(params->cg.d_m); const auto dy = static_cast(params->cg.d_y); const auto pw = static_cast(params->cg.pwr); const bool roll = params->cg.rolloff; const std::array th{thc, thm, thy};//set parameter GUI in th const std::array dl{dc, dm, dy};//set parameter GUI in dl const int height = src->getHeight(); const int width = src->getWidth(); constexpr float range = 65535.f; #ifdef _OPENMP # pragma omp parallel for schedule(dynamic,16) if (multiThread) #endif for (int i = 0; i < height; ++i) for (int j = 0; j < width; ++j) { const float r = src->r(i, j) / range;//in interval 0.. 1 const float g = src->g(i, j) / range; const float b = src->b(i, j) / range; std::array rgb_in{r, g, b}; float rout = 0.f; float gout = 0.f; float bout = 0.f; Color::aces_reference_gamut_compression(rgb_in, th, dl, to_out, from_out, pw, roll, rout, gout, bout); dst->r(i, j) = range * rout;//in interval 0..65535 dst->g(i, j) = range * gout; dst->b(i, j) = range * bout; } } void ImProcFunctions::workingtrc(int sp, Imagefloat* src, Imagefloat* dst, int cw, int ch, int mul, Glib::ustring &profile, double gampos, double slpos, int cat, int &illum, int prim, int locprim, float &rdx, float &rdy, float &grx, float &gry, float &blx, float &bly, float &meanx, float &meany, float &meanxe, float &meanye, cmsHTRANSFORM &transform, bool normalizeIn, bool normalizeOut, bool keepTransForm, bool gamutcontrol) const { const TMatrix wprof = ICCStore::getInstance()->workingSpaceMatrix(params->icm.workingProfile); double wb2[3][3]; float epsilon = 0.000001f; // if(gamutcontrol) { #ifdef _OPENMP #pragma omp parallel for #endif for (int i = 0; i < ch; ++i) for (int j = 0; j < cw; ++j) { src->r(i, j) = (float) rtengine::max(src->r(i, j), epsilon); src->g(i, j) = (float) rtengine::max(src->g(i, j), epsilon); src->b(i, j) = (float) rtengine::max(src->b(i, j), epsilon); } // } if (mul == 5) {//only second pass workingtrc - avoid this code first pass for (int r = 0; r < 3; ++r) { for (int c = 0; c < 3; ++c) { wb2[r][c] = wprof[r][c]; } } //provis - samme approach as in WB itcwb //we can add others function on colors ...others than mean (actually) int precision = 3; const int bfw = cw / precision + ((cw % precision) > 0 ? 1 : 0); const int bfh = ch / precision + ((ch % precision) > 0 ? 1 : 0); Imagefloat *provis = nullptr; provis = new Imagefloat(bfw, bfh);//cw, ch #ifdef _OPENMP #pragma omp parallel for #endif for (int i = 0; i < bfh ; ++i) { const int ii = i * precision; if (ii < ch) { for (int j = 0, jj = 0; j < bfw ; ++j, jj += precision) { provis->r(i, j) = src->r(ii, jj); provis->g(i, j) = src->g(ii, jj); provis->b(i, j) = src->b(ii, jj); } } } // I try to find the dominant color by a simple way (average of x and y) // It is probably intellectually more relevant to place this algorithm at the end, but it is complex at the GUI level (at least for me). // The errors made are relatively minimal and result seems good enough meanx = 0.f; meany = 0.f; #ifdef _OPENMP #pragma omp parallel for reduction(+:meanx, meany) if(multiThread) #endif for (int y = 0; y < bfh ; ++y) { for (int x = 0; x < bfw ; ++x) { const float RR = provis->r(y, x); const float GG = provis->g(y, x); const float BB = provis->b(y, x); float xcb, ycb, zcb; Color::rgbxyz(RR, GG, BB, xcb, ycb, zcb, wb2); float X_r = xcb; float Y_r = ycb; float Z_r = zcb; if(gamutcontrol) { Color::gamutmap(X_r, Y_r, Z_r, wb2);//gamut control } const float som = X_r + Y_r + Z_r; X_r = X_r / som; Y_r = Y_r / som; meanx += X_r; meany += Y_r; } } meanx /= (bfh * bfw); meany /= (bfh * bfw); meanx += 0.005f; meany += 0.005f; //ampirical mean delta with value end in process if (settings->verbose) { printf("Estimation dominant color : x=%f y=%f\n", (double) meanx, (double) meany); } delete provis; } double wprofprim[3][3];//store primaries to XYZ const float toxyz[3][3] = { { static_cast(wprof[0][0] / ((normalizeIn ? 65535.0 : 1.0))), //I have suppressed / Color::D50x static_cast(wprof[0][1] / ((normalizeIn ? 65535.0 : 1.0))), static_cast(wprof[0][2] / ((normalizeIn ? 65535.0 : 1.0))) }, { static_cast(wprof[1][0] / (normalizeIn ? 65535.0 : 1.0)), static_cast(wprof[1][1] / (normalizeIn ? 65535.0 : 1.0)), static_cast(wprof[1][2] / (normalizeIn ? 65535.0 : 1.0)) }, { static_cast(wprof[2][0] / ((normalizeIn ? 65535.0 : 1.0))), //I have suppressed / Color::D50z static_cast(wprof[2][1] / ((normalizeIn ? 65535.0 : 1.0))), static_cast(wprof[2][2] / ((normalizeIn ? 65535.0 : 1.0))) } }; if (profile == "sRGB" || profile == "Adobe RGB" || profile == "ProPhoto" || profile == "WideGamut" || profile == "BruceRGB" || profile == "Beta RGB" || profile == "BestRGB" || profile == "Rec2020" || profile == "ACESp0" || profile == "ACESp1" || profile == "JDCmax" || profile == "JDCmax stdA") { if (settings->verbose) { printf("Profile=%s\n", profile.c_str()); } } else { if (settings->verbose) { printf("profile not accepted\n"); } return; } if (mul == -5 && gampos == 2.4 && slpos == 12.92310) {//must be change if we change settings RT sRGB //only in this case we can shortcut..all process..no gamut control..because we reduce...leads to very small differences, but big speedup #ifdef _OPENMP #pragma omp parallel for schedule(dynamic, 16) if (multiThread) #endif for (int i = 0; i < ch; ++i) for (int j = 0; j < cw; ++j) { double r = (double) src->r(i, j); double g = (double) src->g(i, j); double b = (double) src->b(i, j); r = (Color::igammatab_srgb[r]) / 65535.; g = (Color::igammatab_srgb[g]) / 65535.; b = (Color::igammatab_srgb[b]) / 65535.; dst->r(i, j) = r; dst->g(i, j) = g; dst->b(i, j) = b; } return; } if (mul == 1) { // || (params->icm.wprim == ColorManagementParams::Primaries::DEFAULT && params->icm.will == ColorManagementParams::Illuminant::DEFAULT)) { //shortcut and speedup when no call primaries and illuminant - no gamut control...in this case be careful GammaValues g_a; //gamma parameters double pwr = 1.0 / static_cast(gampos); Color::calcGamma(pwr, slpos, g_a); // call to calcGamma with selected gamma and slope #ifdef _OPENMP # pragma omp parallel for schedule(dynamic,16) if (multiThread) #endif for (int y = 0; y < ch; ++y) { int x = 0; #ifdef __SSE2__ for (; x < cw - 3; x += 4) { STVFU(dst->r(y, x), F2V(65536.f) * gammalog(LVFU(src->r(y, x)), F2V(gampos), F2V(slpos), F2V(g_a[3]), F2V(g_a[4]))); STVFU(dst->g(y, x), F2V(65536.f) * gammalog(LVFU(src->g(y, x)), F2V(gampos), F2V(slpos), F2V(g_a[3]), F2V(g_a[4]))); STVFU(dst->b(y, x), F2V(65536.f) * gammalog(LVFU(src->b(y, x)), F2V(gampos), F2V(slpos), F2V(g_a[3]), F2V(g_a[4]))); } #endif for (; x < cw; ++x) { dst->r(y, x) = 65536.f * gammalog(src->r(y, x), gampos, slpos, g_a[3], g_a[4]); dst->g(y, x) = 65536.f * gammalog(src->g(y, x), gampos, slpos, g_a[3], g_a[4]); dst->b(y, x) = 65536.f * gammalog(src->b(y, x), gampos, slpos, g_a[3], g_a[4]); } } return; } float redxx = params->icm.redx; float redyy = params->icm.redy; float bluxx = params->icm.blux; float bluyy = params->icm.bluy; float grexx = params->icm.grex; float greyy = params->icm.grey; float epsil = 0.0001f; double Wx = 1.0; double Wz = 1.0; cmsCIExyY xyD; if (locprim == 1 && mul == 5) { rdx = params->locallab.spots.at(sp).redxl; rdy = params->locallab.spots.at(sp).redyl; grx = params->locallab.spots.at(sp).grexl; gry = params->locallab.spots.at(sp).greyl; blx = params->locallab.spots.at(sp).bluxl; bly = params->locallab.spots.at(sp).bluyl; if (params->locallab.spots.at(sp).illMethod == "d50") { illum = 2; xyD = {0.3457, 0.3585, 1.0}; // near LCMS values but not perfect... it's a compromise!! Wx = 0.964295676; Wz = 0.825104603; } else if (params->locallab.spots.at(sp).illMethod == "d60") { illum = 4; Wx = 0.952646075; Wz = 1.008825184; xyD = {0.32168, 0.33767, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "d65") { illum = 5; Wx = 0.95045471; Wz = 1.08905029; xyD = {0.312700492, 0.329000939, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "d41") { illum = 1; Wx = 0.991488263; Wz = 0.631604625; xyD = {0.376137, 0.374021, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "d55") { illum = 3; Wx = 0.956565934; Wz = 0.920253249; xyD = {0.332424, 0.347426, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "d80") { illum = 6; Wx = 0.950095542; Wz = 1.284213976; xyD = {0.293756, 0.309185, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "d120") { illum = 7; Wx = 0.979182; Wz = 1.623623; xyD = {0.269669, 0.28078, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "stda") { illum = 8; Wx = 1.098500393; Wz = 0.355848714; xyD = {0.447573, 0.407440, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "T2000") { illum = 9; Wx = 1.274335; Wz = 0.145233; xyD = {0.526591, 0.41331, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "T1500") { illum = 10; Wx = 1.489921; Wz = 0.053826; xyD = {0.585703, 0.393157, 1.0}; } else if (params->locallab.spots.at(sp).illMethod == "iE") { illum = 20; Wx = 1.; Wz = 1.; xyD = {0.333333, 0.333333, 1.0}; } } if (prim == 14 && locprim == 0 && mul == 5) {//convert datas area to xy float redgraphx = params->icm.labgridcieALow; float redgraphy = params->icm.labgridcieBLow; float blugraphx = params->icm.labgridcieAHigh; float blugraphy = params->icm.labgridcieBHigh; float gregraphx = params->icm.labgridcieGx; float gregraphy = params->icm.labgridcieGy; redxx = 0.55f * (redgraphx + 1.f) - 0.1f; redxx = rtengine::LIM(redxx, 0.41f, 1.f);//limit values for xy (arbitrary) redyy = 0.55f * (redgraphy + 1.f) - 0.1f; redyy = rtengine::LIM(redyy, 0.f, 0.7f); bluxx = 0.55f * (blugraphx + 1.f) - 0.1f; bluxx = rtengine::LIM(bluxx, -0.1f, 0.5f); bluyy = 0.55f * (blugraphy + 1.f) - 0.1f; bluyy = rtengine::LIM(bluyy, -0.1f, 0.49f); grexx = 0.55f * (gregraphx + 1.f) - 0.1f; grexx = rtengine::LIM(grexx, -0.1f, 0.4f); greyy = 0.55f * (gregraphy + 1.f) - 0.1f; greyy = rtengine::LIM(greyy, 0.5f, 1.f); } //fixed crash when there is no space or too small..just a line...Possible if bx, by aligned with Gx,Gy Rx,Ry //fix crash if user select 0 for redyy, bluyy, greyy if (redyy == 0.f) { redyy = epsil; } if (bluyy == 0.f) { bluyy = epsil; } if (greyy == 0.f) { greyy = epsil; } //fix crash if grexx - redxx = 0 float grered = 1.f; grered = grexx - redxx; if (grered == 0.f) { grered = epsil; } float ac = (greyy - redyy) / grered; float bc = greyy - ac * grexx; float yc = ac * bluxx + bc; if ((bluyy < yc + 0.0004f) && (bluyy > yc - 0.0004f)) { //under 0.0004 in some case crash because space too small return; } enum class ColorTemp { D50 = 5003, // for Widegamut, ProPhoto Best, Beta -> D50 D65 = 6504, // for sRGB, AdobeRGB, Bruce Rec2020 -> D65 D60 = 6005 // for ACES AP0 and AP1 }; double tempv4 = 5003.; double p[6]; //primaries if (locprim == 0 && mul == 5) { switch (ColorManagementParams::Primaries(prim)) { case ColorManagementParams::Primaries::DEFAULT: { break; } case ColorManagementParams::Primaries::SRGB: { profile = "sRGB"; break; } case ColorManagementParams::Primaries::ADOBE_RGB: { profile = "Adobe RGB"; break; } case ColorManagementParams::Primaries::PRO_PHOTO: { profile = "ProPhoto"; break; } case ColorManagementParams::Primaries::REC2020: { profile = "Rec2020"; break; } case ColorManagementParams::Primaries::ACES_P1: { profile = "ACESp1"; break; } case ColorManagementParams::Primaries::WIDE_GAMUT: { profile = "WideGamut"; break; } case ColorManagementParams::Primaries::ACES_P0: { profile = "ACESp0"; break; } case ColorManagementParams::Primaries::JDC_MAX: { profile = "JDCmax"; break; } case ColorManagementParams::Primaries::JDC_MAXSTDA: { profile = "JDCmax stdA"; break; } case ColorManagementParams::Primaries::BRUCE_RGB: { profile = "BruceRGB"; break; } case ColorManagementParams::Primaries::BETA_RGB: { profile = "Beta RGB"; break; } case ColorManagementParams::Primaries::BEST_RGB: { profile = "BestRGB"; break; } case ColorManagementParams::Primaries::CUSTOM: { profile = "Custom"; break; } case ColorManagementParams::Primaries::CUSTOM_GRID: { profile = "Custom"; break; } } } else if (locprim == 1 && mul == 5) { //local primaries if (prim == 1) { p[0] = 0.6400; // sRGB primaries p[1] = 0.3300; p[2] = 0.3000; p[3] = 0.6000; p[4] = 0.1500; p[5] = 0.0600; tempv4 = 6504.; Wx = 0.95045471; Wz = 1.08905029; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 2) { p[0] = 0.6400; //Adobe primaries p[1] = 0.3300; p[2] = 0.2100; p[3] = 0.7100; p[4] = 0.1500; p[5] = 0.0600; tempv4 = 6504.; Wx = 0.95045471; Wz = 1.08905029; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 3) { p[0] = 0.7347; //ProPhoto and default primaries p[1] = 0.2653; p[2] = 0.1596; p[3] = 0.8404; p[4] = 0.0366; p[5] = 0.0001; Wx = 0.964295676; Wz = 0.825104603; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 4) { p[0] = 0.7080; // Rec2020 primaries p[1] = 0.2920; p[2] = 0.1700; p[3] = 0.7970; p[4] = 0.1310; p[5] = 0.0460; tempv4 = 6504.; Wx = 0.95045471; Wz = 1.08905029; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 5) { p[0] = 0.713; // ACES P1 primaries p[1] = 0.293; p[2] = 0.165; p[3] = 0.830; p[4] = 0.128; p[5] = 0.044; tempv4 = 6004.; Wx = 0.952646075; Wz = 1.008825184; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 6) { p[0] = 0.7350; //Widegamut primaries p[1] = 0.2650; p[2] = 0.1150; p[3] = 0.8260; p[4] = 0.1570; p[5] = 0.0180; Wx = 0.964295676; Wz = 0.825104603; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 7) { p[0] = 0.7347; //ACESp0 primaries p[1] = 0.2653; p[2] = 0.; p[3] = 1.0; p[4] = 0.0001; p[5] = -0.0770; tempv4 = 6004.; Wx = 0.952646075; Wz = 1.008825184; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 8) { p[0] = 0.734702; // JDC max primaries p[1] = 0.265302; p[2] = 0.021908; p[3] = 0.930288; p[4] = 0.120593; p[5] = 0.001583; Wx = 0.964295676; Wz = 0.825104603; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 9) { p[0] = 0.734702; // JDC max primaries p[1] = 0.265302; p[2] = 0.021908; p[3] = 0.930288; p[4] = 0.120593; p[5] = 0.001583; Wx = 1.098500393; Wz = 0.355848714; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 10) { p[0] = 0.64; // Bruce primaries p[1] = 0.33; p[2] = 0.28; p[3] = 0.65; p[4] = 0.15; p[5] = 0.06; Wx = 0.95045471; Wz = 1.08905029; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 11) { p[0] = 0.6888; // Beta primaries p[1] = 0.3112; p[2] = 0.1986; p[3] = 0.7551; p[4] = 0.1265; p[5] = 0.0352; Wx = 0.964295676; Wz = 0.825104603; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 12) { p[0] = 0.7347; // Best primaries p[1] = 0.2653; p[2] = 0.2150; p[3] = 0.7750; p[4] = 0.13; p[5] = 0.0350; Wx = 0.964295676; Wz = 0.825104603; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } else if (prim == 15) { p[0] = rdx; p[1] = rdy; p[2] = grx; p[3] = gry; p[4] = blx; p[5] = bly; } else { p[0] = 0.7347; //ProPhoto and default primaries p[1] = 0.2653; p[2] = 0.1596; p[3] = 0.8404; p[4] = 0.0366; p[5] = 0.0001; Wx = 0.964295676; Wz = 0.825104603; rdx = p[0]; rdy = p[1]; grx = p[2]; gry = p[3]; blx = p[4]; bly = p[5]; } } if (settings->verbose && prim != 0) { printf("prim=%i Profile Destination=%s\n", prim, profile.c_str()); } cmsHTRANSFORM hTransform = nullptr; if (transform) { hTransform = transform; } else { double pwr = 1.0 / gampos; double ts = slpos; int five = mul; if (gampos < 1.0) { pwr = gampos; gampos = 1. / gampos; five = -mul; } // int select_temp = 1; //5003K constexpr double eps = 0.000000001; // not divide by zero //primaries for 10 working profiles ==> output profiles if (locprim == 0 && mul ==5) { if (profile == "WideGamut") { p[0] = 0.7350; //Widegamut primaries p[1] = 0.2650; p[2] = 0.1150; p[3] = 0.8260; p[4] = 0.1570; p[5] = 0.0180; illum = toUnderlying(ColorManagementParams::Illuminant::D50); Wx = 0.964295676; Wz = 0.825104603; } else if (profile == "Adobe RGB") { p[0] = 0.6400; //Adobe primaries p[1] = 0.3300; p[2] = 0.2100; p[3] = 0.7100; p[4] = 0.1500; p[5] = 0.0600; tempv4 = 6504.; illum = toUnderlying(ColorManagementParams::Illuminant::D65); Wx = 0.95045471; Wz = 1.08905029; } else if (profile == "sRGB") { p[0] = 0.6400; // sRGB primaries p[1] = 0.3300; p[2] = 0.3000; p[3] = 0.6000; p[4] = 0.1500; p[5] = 0.0600; tempv4 = 6504.; illum = toUnderlying(ColorManagementParams::Illuminant::D65); Wx = 0.95045471; Wz = 1.08905029; } else if (profile == "BruceRGB") { p[0] = 0.6400; // Bruce primaries p[1] = 0.3300; p[2] = 0.2800; p[3] = 0.6500; p[4] = 0.1500; p[5] = 0.0600; tempv4 = 6504.; illum = toUnderlying(ColorManagementParams::Illuminant::D65); Wx = 0.95045471; Wz = 1.08905029; } else if (profile == "Beta RGB") { p[0] = 0.6888; // Beta primaries p[1] = 0.3112; p[2] = 0.1986; p[3] = 0.7551; p[4] = 0.1265; p[5] = 0.0352; illum = toUnderlying(ColorManagementParams::Illuminant::D50); Wx = 0.964295676; Wz = 0.825104603; } else if (profile == "BestRGB") { p[0] = 0.7347; // Best primaries p[1] = 0.2653; p[2] = 0.2150; p[3] = 0.7750; p[4] = 0.1300; p[5] = 0.0350; illum = toUnderlying(ColorManagementParams::Illuminant::D50); Wx = 0.964295676; Wz = 0.825104603; } else if (profile == "Rec2020") { p[0] = 0.7080; // Rec2020 primaries p[1] = 0.2920; p[2] = 0.1700; p[3] = 0.7970; p[4] = 0.1310; p[5] = 0.0460; tempv4 = 6504.; illum = toUnderlying(ColorManagementParams::Illuminant::D65); Wx = 0.95045471; Wz = 1.08905029; } else if (profile == "ACESp0") { p[0] = 0.7347; // ACES P0 primaries p[1] = 0.2653; p[2] = 0.0000; p[3] = 1.0; p[4] = 0.0001; p[5] = -0.0770; tempv4 = 6004.; illum = toUnderlying(ColorManagementParams::Illuminant::D60); Wx = 0.952646075; Wz = 1.008825184; } else if (profile == "JDCmax") { p[0] = 0.734702; // JDC max primaries p[1] = 0.265302; p[2] = 0.021908; p[3] = 0.930288; p[4] = 0.120593; p[5] = 0.001583; illum = toUnderlying(ColorManagementParams::Illuminant::D50); Wx = 0.964295676; Wz = 0.825104603; } else if (profile == "JDCmax stdA") { p[0] = 0.734702; // JDC max primaries and stdA p[1] = 0.265302; p[2] = 0.021908; p[3] = 0.930288; p[4] = 0.120593; p[5] = 0.001583; illum = toUnderlying(ColorManagementParams::Illuminant::STDA); Wx = 1.098500393; Wz = 0.355848714; } else if (profile == "ACESp1") { p[0] = 0.713; // ACES P1 primaries p[1] = 0.293; p[2] = 0.165; p[3] = 0.830; p[4] = 0.128; p[5] = 0.044; tempv4 = 6004.; illum = toUnderlying(ColorManagementParams::Illuminant::D60); Wx = 0.952646075; Wz = 1.008825184; } else if (profile == "ProPhoto") { p[0] = 0.7347; //ProPhoto and default primaries p[1] = 0.2653; p[2] = 0.1596; p[3] = 0.8404; p[4] = 0.0366; p[5] = 0.0001; illum = toUnderlying(ColorManagementParams::Illuminant::D50); Wx = 0.964295676; Wz = 0.825104603; } else if (profile == "Custom") { p[0] = redxx; p[1] = redyy; p[2] = grexx; p[3] = greyy; p[4] = bluxx; p[5] = bluyy; } else { p[0] = 0.7347; //default primaries always unused p[1] = 0.2653; p[2] = 0.1596; p[3] = 0.8404; p[4] = 0.0366; p[5] = 0.0001; } } if (slpos == 0) { slpos = eps; } GammaValues g_a; //gamma parameters Color::calcGamma(pwr, ts, g_a); // call to calcGamma with selected gamma and slope : return parameters for LCMS2 cmsFloat64Number gammaParams[7]; gammaParams[4] = g_a[3] * ts; gammaParams[0] = gampos; gammaParams[1] = 1. / (1.0 + g_a[4]); gammaParams[2] = g_a[4] / (1.0 + g_a[4]); gammaParams[3] = 1. / slpos; gammaParams[5] = 0.0; gammaParams[6] = 0.0; if(rtengine::settings->verbose) { printf("ga0=%f ga1=%f ga2=%f ga3=%f ga4=%f\n", gammaParams[0], gammaParams[1], gammaParams[2], gammaParams[3], gammaParams[4]); } // 7 parameters for smoother curves // cmsCIExyY xyD; Glib::ustring ills = "D50"; if (locprim == 0 && mul == 5) { switch (ColorManagementParams::Illuminant(illum)) { case ColorManagementParams::Illuminant::DEFAULT: case ColorManagementParams::Illuminant::STDA: case ColorManagementParams::Illuminant::TUNGSTEN_2000K: case ColorManagementParams::Illuminant::TUNGSTEN_1500K: case ColorManagementParams::Illuminant::E:{ break; } case ColorManagementParams::Illuminant::D41: { tempv4 = 4100.; ills = "D41"; break; } case ColorManagementParams::Illuminant::D50: { tempv4 = 5003.; ills = "D50"; break; } case ColorManagementParams::Illuminant::D55: { tempv4 = 5500.; ills = "D55"; break; } case ColorManagementParams::Illuminant::D60: { tempv4 = 6004.; ills = "D60"; break; } case ColorManagementParams::Illuminant::D65: { tempv4 = 6504.; ills = "D65"; break; } case ColorManagementParams::Illuminant::D80: { tempv4 = 8000.; ills = "D80"; break; } case ColorManagementParams::Illuminant::D120: { tempv4 = 12000.; ills = "D120"; break; } } cmsWhitePointFromTemp(&xyD, tempv4); switch (ColorManagementParams::Illuminant(illum)) { case ColorManagementParams::Illuminant::DEFAULT: { break; } case ColorManagementParams::Illuminant::D55: { Wx = 0.956565934; Wz = 0.920253249; break; } case ColorManagementParams::Illuminant::D80: { Wx = 0.950095542; Wz = 1.284213976; break; } case ColorManagementParams::Illuminant::D41: { Wx = 0.991488263; Wz = 0.631604625; break; } case ColorManagementParams::Illuminant::D50: { xyD = {0.3457, 0.3585, 1.0}; // near LCMS values but not perfect... it's a compromise!! Wx = 0.964295676; Wz = 0.825104603; break; } case ColorManagementParams::Illuminant::D60: { Wx = 0.952646075; Wz = 1.008825184; xyD = {0.32168, 0.33767, 1.0}; break; } case ColorManagementParams::Illuminant::D65: { Wx = 0.95045471; Wz = 1.08905029; xyD = {0.312700492, 0.329000939, 1.0}; break; } case ColorManagementParams::Illuminant::D120: { Wx = 0.979182; Wz = 1.623623; xyD = {0.269669, 0.28078, 1.0}; break; } case ColorManagementParams::Illuminant::STDA: { Wx = 1.098500393; Wz = 0.355848714; xyD = {0.447573, 0.407440, 1.0}; ills = "stdA 2875K"; break; } case ColorManagementParams::Illuminant::TUNGSTEN_2000K: { Wx = 1.274335; Wz = 0.145233; xyD = {0.526591, 0.41331, 1.0}; ills = "Tungsten 2000K"; break; } case ColorManagementParams::Illuminant::TUNGSTEN_1500K: { Wx = 1.489921; Wz = 0.053826; xyD = {0.585703, 0.393157, 1.0}; ills = "Tungsten 1500K"; break; } case ColorManagementParams::Illuminant::E: { Wx = 1.; Wz = 1.; xyD = {0.33333, 0.33333, 1.0}; ills = "E"; break; } } } //xyD //meanx, meany // adjust refinement (purity) with a simple algorithm if (mul == 5) { double refin = 0.; if (locprim == 1) { refin = params->locallab.spots.at(sp).refi; meanx += params->locallab.spots.at(sp).shiftxl; meany += params->locallab.spots.at(sp).shiftyl; } else if (locprim == 0) { refin = params->icm.refi; meanx += params->icm.shiftx; meany += params->icm.shifty; } double arefi = (xyD.y - meany) / (xyD.x - meanx); double brefi = xyD.y - arefi * xyD.x; double scalrefi = 0.98 * (meanx - xyD.x); xyD.x = xyD.x + scalrefi * refin; xyD.y = xyD.x * arefi + brefi; // recalculate Wx Wy Wx = xyD.x / xyD.y; Wz = (1. - xyD.x - xyD.y) / xyD.y; } double wprofpri[9]; //xyz in function primaries and illuminant Color::primaries_to_xyz(p, Wx, Wz, wprofpri, cat); for (int i = 0; i < 3; ++i) { for (int j = 0; j < 3; ++j) { wprofprim[i][j] = (double) wprofpri[j * 3 + i]; //xyz in TMatrix format } } //D41 0.377984 0.381229 //D55 0.332424 0.347426 //D80 0.293755 0.309185 //D75 0.299021 0.314852 cmsToneCurve* GammaTRC[3]; GammaTRC[0] = GammaTRC[1] = GammaTRC[2] = cmsBuildParametricToneCurve(NULL, five, gammaParams);//5 = more smoother than 4 cmsHPROFILE oprofdef = nullptr; const cmsCIExyYTRIPLE Primaries = { {p[0], p[1], 1.0}, // red {p[2], p[3], 1.0}, // green {p[4], p[5], 1.0} // blue }; oprofdef = cmsCreateRGBProfile(&xyD, &Primaries, GammaTRC); cmsWriteTag(oprofdef, cmsSigRedTRCTag, GammaTRC[0]); cmsWriteTag(oprofdef, cmsSigGreenTRCTag, GammaTRC[1]); cmsWriteTag(oprofdef, cmsSigBlueTRCTag, GammaTRC[2]); //to read XYZ values and illuminant if (rtengine::settings->verbose && mul == 5) { cmsCIEXYZ *redT = static_cast(cmsReadTag(oprofdef, cmsSigRedMatrixColumnTag)); cmsCIEXYZ *greenT = static_cast(cmsReadTag(oprofdef, cmsSigGreenMatrixColumnTag)); cmsCIEXYZ *blueT = static_cast(cmsReadTag(oprofdef, cmsSigBlueMatrixColumnTag)); if (locprim == 0) { printf("Illuminant=%s\n", ills.c_str()); } else { printf("Illuminant=%s\n", params->locallab.spots.at(sp).illMethod.c_str()); } printf("rX=%f gX=%f bX=%f\n", redT->X, greenT->X, blueT->X); printf("rY=%f gY=%f bY=%f\n", redT->Y, greenT->Y, blueT->Y); printf("rZ=%f gZ=%f bZ=%f\n", redT->Z, greenT->Z, blueT->Z); } cmsFreeToneCurve(GammaTRC[0]); if (oprofdef) { constexpr cmsUInt32Number flags = cmsFLAGS_NOOPTIMIZE | cmsFLAGS_NOCACHE | cmsFLAGS_BLACKPOINTCOMPENSATION | cmsFLAGS_GAMUTCHECK; const cmsHPROFILE iprof = ICCStore::getInstance()->getXYZProfile(); lcmsMutex->lock(); hTransform = cmsCreateTransform(iprof, TYPE_RGB_FLT, oprofdef, TYPE_RGB_FLT, params->icm.aRendIntent, flags); lcmsMutex->unlock(); } } if (hTransform) { #ifdef _OPENMP #pragma omp parallel if (multiThread) #endif { AlignedBuffer pBuf(cw * 3); const float normalize = normalizeOut ? 65535.f : 1.f; #ifdef _OPENMP #pragma omp for schedule(dynamic, 16) nowait #endif for (int i = 0; i < ch; ++i) { float *p = pBuf.data; for (int j = 0; j < cw; ++j) { const float r = src->r(i, j); const float g = src->g(i, j); const float b = src->b(i, j); float X = toxyz[0][0] * r + toxyz[0][1] * g + toxyz[0][2] * b; float Y = toxyz[1][0] * r + toxyz[1][1] * g + toxyz[1][2] * b; float Z = toxyz[2][0] * r + toxyz[2][1] * g + toxyz[2][2] * b; if (gamutcontrol) { Color::gamutmap(X, Y, Z, wprofprim);//gamut control } *(p++) = X; *(p++) = Y; *(p++) = Z; } p = pBuf.data; cmsDoTransform(hTransform, p, p, cw); for (int j = 0; j < cw; ++j) { dst->r(i, j) = *(p++) * normalize; dst->g(i, j) = *(p++) * normalize; dst->b(i, j) = *(p++) * normalize; } } } // alternative to find dominant color xy // Not use : // 1) GUI complex at least for mean // 2) small difference for meanxe, meanye with meanx , meany above in most cases /* if (locprim == 1) { meanxe = 0.f; meanye = 0.f; #ifdef _OPENMP #pragma omp parallel for reduction(+:meanxe, meanye) if(multiThread) #endif for (int y = 0; y < ch ; ++y) { for (int x = 0; x < cw ; ++x) { const float RR = dst->r(y,x); const float GG = dst->g(y,x); const float BB = dst->b(y,x); float xcb, ycb, zcb; Color::rgbxyz(RR, GG, BB, xcb, ycb, zcb, wb2);//use sRGB Adobe Rec2020 ACESp0 float X_r = xcb; float Y_r = ycb; float Z_r = zcb; Color::gamutmap(X_r, Y_r, Z_r, wb2);//gamut control const float som = X_r + Y_r + Z_r; X_r = X_r / som; Y_r = Y_r / som; meanxe += X_r; meanye += Y_r; } } meanxe /= (ch*cw); meanye /= (ch*cw); printf("DiffmeanxE=%f DiffmeanyE=%f \n", (double) (meanxe - meanx), (double) (meanye - meany)); } */ if (!keepTransForm) { cmsDeleteTransform(hTransform); hTransform = nullptr; } transform = hTransform; } } }