diff --git a/rtengine/improccoordinator.cc b/rtengine/improccoordinator.cc index b6862ab66..7434776d5 100644 --- a/rtengine/improccoordinator.cc +++ b/rtengine/improccoordinator.cc @@ -55,7 +55,7 @@ ImProcCoordinator::ImProcCoordinator () lhist16(65536), lhist16Cropped(65536), lhist16CAM(65536), lhist16CroppedCAM(65536), lhist16CCAM(65536), - lhist16RETI(65536), + lhist16RETI(), histCropped(65536), lhist16Clad(65536), lhist16CLlad(65536), lhist16LClad(65536), lhist16LLClad(65536), @@ -236,6 +236,7 @@ void ImProcCoordinator::updatePreviewImage (int todo, Crop* cropCall) } if (params.retinex.enabled) { + lhist16RETI(32768); lhist16RETI.clear(); imgsrc->retinexPrepareBuffers(params.icm, params.retinex, conversionBuffer, lhist16RETI); diff --git a/rtengine/ipretinex.cc b/rtengine/ipretinex.cc index f66fb39af..6773e7a12 100644 --- a/rtengine/ipretinex.cc +++ b/rtengine/ipretinex.cc @@ -47,7 +47,6 @@ #include "opthelper.h" #include "StopWatch.h" -#define MAX_RETINEX_SCALES 8 #define clipretinex( val, minv, maxv ) (( val = (val < minv ? minv : val ) ) > maxv ? maxv : val ) #define med3(a0,a1,a2,a3,a4,a5,a6,a7,a8,median) { \ @@ -60,13 +59,9 @@ PIX_SORT(pp[3],pp[6]); PIX_SORT(pp[1],pp[4]); PIX_SORT(pp[2],pp[5]); \ PIX_SORT(pp[4],pp[7]); PIX_SORT(pp[4],pp[2]); PIX_SORT(pp[6],pp[4]); \ PIX_SORT(pp[4],pp[2]); median=pp[4];} //pp4 = median -namespace rtengine + +namespace { - -extern const Settings* settings; - -static float RetinexScales[MAX_RETINEX_SCALES]; - void retinex_scales( float* scales, int nscales, int mode, int s, float high) { if ( nscales == 1 ) { @@ -102,6 +97,7 @@ void retinex_scales( float* scales, int nscales, int mode, int s, float high) } } } + void mean_stddv2( float **dst, float &mean, float &stddv, int W_L, int H_L, float &maxtr, float &mintr) { // summation using double precision to avoid too large summation error for large pictures @@ -123,14 +119,8 @@ void mean_stddv2( float **dst, float &mean, float &stddv, int W_L, int H_L, floa sum += dst[i][j]; vsquared += (dst[i][j] * dst[i][j]); - if ( dst[i][j] > lmax) { - lmax = dst[i][j] ; - } - - if ( dst[i][j] < lmin) { - lmin = dst[i][j] ; - } - + lmax = dst[i][j] > lmax ? dst[i][j] : lmax; + lmin = dst[i][j] < lmin ? dst[i][j] : lmin; } #ifdef _OPENMP @@ -148,81 +138,30 @@ void mean_stddv2( float **dst, float &mean, float &stddv, int W_L, int H_L, floa stddv = (float)sqrt(stddv); } - - - - - -void mean_stddv( float **dst, float &mean, float &stddv, int W_L, int H_L, const float factor, float &maxtr, float &mintr) - -{ - // summation using double precision to avoid too large summation error for large pictures - double vsquared = 0.f; - double sum = 0.f; - maxtr = 0.f; - mintr = 0.f; -#ifdef _OPENMP - #pragma omp parallel -#endif - { - float lmax = 0.f, lmin = 0.f; - -#ifdef _OPENMP - #pragma omp for reduction(+:sum,vsquared) // this can lead to differences, but parallel summation is more accurate -#endif - - for (int i = 0; i < H_L; i++ ) - for (int j = 0; j < W_L; j++) { - sum += dst[i][j]; - vsquared += (dst[i][j] * dst[i][j]); - - if ( dst[i][j] > lmax) { - lmax = dst[i][j] ; - } - - if ( dst[i][j] < lmin) { - lmin = dst[i][j] ; - } - - } - -#ifdef _OPENMP - #pragma omp critical -#endif - { - maxtr = maxtr > lmax ? maxtr : lmax; - mintr = mintr < lmin ? mintr : lmin; - } - - } - - sum *= factor; - maxtr *= factor; - mintr *= factor; - vsquared *= (factor * factor); - mean = sum / (float) (W_L * H_L); - vsquared /= (float) W_L * H_L; - stddv = ( vsquared - (mean * mean) ); - stddv = (float)sqrt(stddv); } + +namespace rtengine +{ + +extern const Settings* settings; + void RawImageSource::MSR(float** luminance, float** originalLuminance, float **exLuminance, LUTf & mapcurve, bool &mapcontlutili, int width, int height, RetinexParams deh, const RetinextransmissionCurve & dehatransmissionCurve, const RetinexgaintransmissionCurve & dehagaintransmissionCurve, float &minCD, float &maxCD, float &mini, float &maxi, float &Tmean, float &Tsigma, float &Tmin, float &Tmax) { if (deh.enabled) {//enabled float mean, stddv, maxtr, mintr; - //float mini, delta, maxi; float delta; - float eps = 2.f; + constexpr float eps = 2.f; bool useHsl = deh.retinexcolorspace == "HSLLOG"; bool useHslLin = deh.retinexcolorspace == "HSLLIN"; float gain2 = (float) deh.gain / 100.f; //def =1 not use gain2 = useHslLin ? gain2 * 0.5f : gain2; float offse = (float) deh.offs; //def = 0 not use int iter = deh.iter; - int gradient = deh.scal; + int gradient = deh.scal; int scal = 3;//disabled scal - int nei = (int) 2.8f * deh.neigh; //def = 220 + int nei = (int) (2.8f * deh.neigh); //def = 220 float vart = (float)deh.vart / 100.f;//variance float gradvart = (float)deh.grad; float gradstr = (float)deh.grads; @@ -231,9 +170,7 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e limD = pow(limD, 1.7f);//about 2500 enough limD *= useHslLin ? 10.f : 1.f; float ilimD = 1.f / limD; - int moderetinex = 2; // default to 2 ( deh.retinexMethod == "high" ) float hig = ((float) deh.highl) / 100.f; - bool higplus = false ; float elogt; float hl = deh.baselog; scal = deh.skal; @@ -248,50 +185,43 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e int W_L = width; float *tran[H_L] ALIGNED16; - float *tranBuffer; - int viewmet = 0; + float *tranBuffer = nullptr; elogt = 2.71828f;//disabled baselog - FlatCurve* shcurve = NULL;//curve L=f(H) bool lhutili = false; - if (deh.enabled) { - shcurve = new FlatCurve(deh.lhcurve); + FlatCurve* shcurve = new FlatCurve(deh.lhcurve); //curve L=f(H) - if (!shcurve || shcurve->isIdentity()) { - if (shcurve) { - delete shcurve; - shcurve = NULL; - } - } else { - lhutili = true; + if (!shcurve || shcurve->isIdentity()) { + if (shcurve) { + delete shcurve; + shcurve = nullptr; } + } else { + lhutili = true; } + bool higplus = false ; + int moderetinex = 2; // default to 2 ( deh.retinexMethod == "high" ) + if(deh.retinexMethod == "highliplus") { higplus = true; - } - - if (deh.retinexMethod == "uni") { + moderetinex = 3; + } else if (deh.retinexMethod == "uni") { moderetinex = 0; - } - - if (deh.retinexMethod == "low") { + } else if (deh.retinexMethod == "low") { moderetinex = 1; - } - - if (deh.retinexMethod == "highli" || deh.retinexMethod == "highliplus") { + } else { /*if (deh.retinexMethod == "highli") */ moderetinex = 3; } - for(int it = 1; it < iter + 1; it++) { //iter nb max of iterations - float aahi = 49.f / 99.f; ////reduce sensibility 50% - float bbhi = 1.f - aahi; - float high; - high = bbhi + aahi * (float) deh.highl; + constexpr float aahi = 49.f / 99.f; ////reduce sensibility 50% + constexpr float bbhi = 1.f - aahi; + + for(int it = 1; it < iter + 1; it++) { //iter nb max of iterations + float high = bbhi + aahi * (float) deh.highl; - float grads; float grad = 1.f; float sc = scal; @@ -382,7 +312,6 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e } scal = round(sc); - float strengthx; float ks = 1.f; if(gradstr != 0) { @@ -413,8 +342,11 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e } } - strengthx = ks * strength; - //printf("scale=%d\n", scal); + float strengthx = ks * strength; + + constexpr auto maxRetinexScales = 8; + float RetinexScales[maxRetinexScales]; + retinex_scales( RetinexScales, scal, moderetinex, nei / grad, high ); float *src[H_L] ALIGNED16; @@ -428,39 +360,28 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e int shHighlights = deh.highlights; int shShadows = deh.shadows; + int mapmet = 0; if(deh.mapMethod == "map") { mapmet = 2; - } - - if(deh.mapMethod == "mapT") { + } else if(deh.mapMethod == "mapT") { mapmet = 3; - } - - /*if(deh.mapMethod == "curv") { - mapmet = 1; - }*/ - - if(deh.mapMethod == "gaus") { + } else if(deh.mapMethod == "gaus") { mapmet = 4; } const double shradius = mapmet == 4 ? (double) deh.radius : 40.; + int viewmet = 0; + if(deh.viewMethod == "mask") { viewmet = 1; - } - - if(deh.viewMethod == "tran") { + } else if(deh.viewMethod == "tran") { viewmet = 2; - } - - if(deh.viewMethod == "tran2") { + } else if(deh.viewMethod == "tran2") { viewmet = 3; - } - - if(deh.viewMethod == "unsharp") { + } else if(deh.viewMethod == "unsharp") { viewmet = 4; } @@ -554,6 +475,12 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e #endif if(mapmet > 0 && mapcontlutili && it == 1) { + // TODO: When rgbcurvespeedup branch is merged into master we can simplify the code by + // 1) in rawimagesource.retinexPrepareCurves() insert + // mapcurve *= 0.5f; + // after + // CurveFactory::mapcurve (mapcontlutili, retinexParams.mapcurve, mapcurve, 1, lhist16RETI, histLRETI); + // 2) remove the division by 2.f from the code 7 lines below this line #ifdef _OPENMP #pragma omp parallel for #endif @@ -567,21 +494,21 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e } if(((mapmet == 2 && scale > 2) || mapmet == 3 || mapmet == 4) && it == 1) { - - + float hWeight = (100.f - shHighlights) / 100.f; + float sWeight = (100.f - shShadows) / 100.f; #ifdef _OPENMP - #pragma omp parallel for + #pragma omp parallel for schedule(dynamic,16) #endif for (int i = 0; i < H_L; i++) { for (int j = 0; j < W_L; j++) { - double mapval = 1.0 + shmap->map[i][j]; - double factor = 1.0; + float mapval = 1.f + shmap->map[i][j]; + float factor = 1.f; if (mapval > h_th) { - factor = (h_th + (100.0 - shHighlights) * (mapval - h_th) / 100.0) / mapval; + factor = (h_th + hWeight * (mapval - h_th)) / mapval; } else if (mapval < s_th) { - factor = (s_th - (100.0 - shShadows) * (s_th - mapval) / 100.0) / mapval; + factor = (s_th - sWeight * (s_th - mapval)) / mapval; } out[i][j] *= factor; @@ -629,10 +556,9 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e } } - shmap = NULL; + shmap = nullptr; delete [] buffer; - //delete [] outBuffer; delete [] srcBuffer; mean = 0.f; @@ -657,6 +583,7 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e bmax *= 500.f; amin *= 500.f; bmin *= 500.f; + #ifdef _OPENMP #pragma omp parallel #endif @@ -676,6 +603,7 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e absciss = amin * luminance[i][j] + bmin; } + //TODO : move multiplication by 4.f and subtraction of 1.f inside the curve luminance[i][j] *= (-1.f + 4.f * dehatransmissionCurve[absciss]); //new transmission if(viewmet == 3 || viewmet == 2) { @@ -797,10 +725,9 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e #pragma omp parallel #endif { - float absciss; float cdmax = -999999.f, cdmin = 999999.f; #ifdef _OPENMP - #pragma omp for + #pragma omp for schedule(dynamic,16) nowait #endif for ( int i = 0; i < H_L; i ++ ) @@ -808,6 +735,7 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e float gan; if (dehagaintransmissionCurve && mean != 0.f && stddv != 0.f) { + float absciss; if (LIKELY(fabsf(luminance[i][j] - mean) < stddv)) { absciss = asig * luminance[i][j] + bsig; @@ -819,6 +747,7 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e // float cd = cdfactor * ( luminance[i][j] - mini ) + offse; + // TODO : move multiplication by 2.f inside the curve gan = 2.f * (dehagaintransmissionCurve[absciss]); //new gain function transmission } else { gan = 0.5f; @@ -826,17 +755,12 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e float cd = gan * cdfactor * ( luminance[i][j] ) + offse; - if(cd > cdmax) { - cdmax = cd; - } - - if(cd < cdmin) { - cdmin = cd; - } + cdmax = cd > cdmax ? cd : cdmax; + cdmin = cd < cdmin ? cd : cdmin; float str = strengthx; - if(lhutili && it == 1) { // S=f(H) + if(lhutili && it == 1) { // S=f(H) { float HH = exLuminance[i][j]; float valparam; @@ -851,31 +775,23 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e } } - if(exLuminance[i][j] > 65535.f * hig && higplus) { + if(higplus && exLuminance[i][j] > 65535.f * hig) { str *= hig; } if(viewmet == 0) { - luminance[i][j] = clipretinex( cd, 0.f, 32768.f ) * str + (1.f - str) * originalLuminance[i][j]; - } - - if(viewmet == 1) { + luminance[i][j] = intp(str, clipretinex( cd, 0.f, 32768.f ), originalLuminance[i][j]); + } else if(viewmet == 1) { luminance[i][j] = out[i][j]; - } - - if(viewmet == 4) { - luminance[i][j] = (1.f + str) * originalLuminance[i][j] - str * out[i][j]; //unsharp - } - - if(viewmet == 2) { + } else if(viewmet == 4) { + luminance[i][j] = originalLuminance[i][j] + str * (originalLuminance[i][j] - out[i][j]);//unsharp + } else if(viewmet == 2) { if(tran[i][j] <= mean) { luminance[i][j] = azb + aza * tran[i][j]; //auto values } else { luminance[i][j] = bzb + bza * tran[i][j]; } - } - - if(viewmet == 3) { + } else { /*if(viewmet == 3) */ luminance[i][j] = 1000.f + tran[i][j] * 700.f; //arbitrary values to help display log values which are between -20 to + 30 - usage values -4 + 5 } @@ -890,23 +806,23 @@ void RawImageSource::MSR(float** luminance, float** originalLuminance, float **e } } + delete [] outBuffer; - outBuffer = NULL; + outBuffer = nullptr; //printf("cdmin=%f cdmax=%f\n",minCD, maxCD); Tmean = mean; Tsigma = stddv; Tmin = mintr; Tmax = maxtr; - - if (shcurve && it == 1) { + if (shcurve) { delete shcurve; + shcurve = nullptr; } } - if(viewmet == 3 || viewmet == 2) { + if(tranBuffer) { delete [] tranBuffer; - tranBuffer = NULL; } } diff --git a/rtengine/rawimagesource.cc b/rtengine/rawimagesource.cc index 27ddb6b3c..2780c033f 100644 --- a/rtengine/rawimagesource.cc +++ b/rtengine/rawimagesource.cc @@ -1764,7 +1764,7 @@ void RawImageSource::preprocess (const RAWParams &raw, const LensProfParams &le LCPProfile *pLCPProf = lcpStore->getProfile(lensProf.lcpFile); if (pLCPProf) { // don't check focal length to allow distortion correction for lenses without chip, also pass dummy focal length 1 in case of 0 - LCPMapper map(pLCPProf, max(idata->getFocalLen(),1.0), idata->getFocalLen35mm(), idata->getFocusDist(), idata->getFNumber(), true, false, W, H, coarse, -1); + LCPMapper map(pLCPProf, max(idata->getFocalLen(), 1.0), idata->getFocalLen35mm(), idata->getFocusDist(), idata->getFNumber(), true, false, W, H, coarse, -1); #ifdef _OPENMP #pragma omp parallel for @@ -1970,7 +1970,6 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar LUTf *retinexgamtab;//gamma before and after Retinex to restore tones LUTf lutTonereti; - lutTonereti(65536); if(retinexParams.gammaretinex == "low") { retinexgamtab = &(Color::gammatab_115_2); @@ -1986,24 +1985,33 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar double gamm2 = retinexParams.gam; if(gamm2 < 1.) { - pwr = 1. / pwr; - gamm = 1. / gamm; + std::swap(pwr, gamm); } int mode = 0, imax = 0; Color::calcGamma(pwr, ts, mode, imax, g_a0, g_a1, g_a2, g_a3, g_a4, g_a5); // call to calcGamma with selected gamma and slope // printf("g_a0=%f g_a1=%f g_a2=%f g_a3=%f g_a4=%f\n", g_a0,g_a1,g_a2,g_a3,g_a4); + double start; + double add; + + if(gamm2 < 1.) { + start = g_a2; + add = g_a4; + } else { + start = g_a3; + add = g_a4; + } + + double mul = 1. + g_a4; + + lutTonereti(65536); + for (int i = 0; i < 65536; i++) { double val = (i) / 65535.; - double start = g_a3; - double add = g_a4; - double mul = 1. + g_a4; double x; if(gamm2 < 1.) { - start = g_a2; - add = g_a4; x = Color::igammareti (val, gamm, start, ts, mul , add); } else { x = Color::gammareti (val, gamm, start, ts, mul , add); @@ -2055,6 +2063,7 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar } */ if(retinexParams.gammaretinex != "none" && retinexParams.str != 0) {//gamma + #ifdef _OPENMP #pragma omp parallel for #endif @@ -2083,7 +2092,7 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar if(lhist16RETI) { - lhist16RETIThr(32769, 0); + lhist16RETIThr(lhist16RETI.getSize()); lhist16RETIThr.clear(); } @@ -2101,7 +2110,6 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar for (; j < W - border - 3; j += 4) { vfloat H, S, L; - int pos; Color::rgb2hsl(LVFU(red[i][j]), LVFU(green[i][j]), LVFU(blue[i][j]), H, S, L); STVFU(conversionBuffer[0][i - border][j - border], H); STVFU(conversionBuffer[1][i - border][j - border], S); @@ -2111,7 +2119,7 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar if(lhist16RETI) { for(int p = 0; p < 4; p++) { - pos = (int) clipretinex( conversionBuffer[2][i - border][j - border + p], 0.f, 32768.f);//histogram in curve HSL + int pos = ( conversionBuffer[2][i - border][j - border + p]);//histogram in curve HSL lhist16RETIThr[pos]++; } } @@ -2121,14 +2129,13 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar for (; j < W - border; j++) { float H, S, L; - int pos; //rgb=>lab Color::rgb2hslfloat(red[i][j], green[i][j], blue[i][j], conversionBuffer[0][i - border][j - border], conversionBuffer[1][i - border][j - border], L); L *= 32768.f; conversionBuffer[2][i - border][j - border] = L; if(lhist16RETI) { - pos = (int) clipretinex(L, 0, 32768); + int pos = L; lhist16RETIThr[pos]++; } } @@ -2139,10 +2146,7 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar { if(lhist16RETI) { - // Add per Thread LUT to global LUT - for(int i = 0; i < 32769; i++) { - lhist16RETI[i] += lhist16RETIThr[i]; - } + lhist16RETI += lhist16RETIThr; // Add per Thread LUT to global LUT } } #endif @@ -2150,10 +2154,10 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar } } else { TMatrix wprof = iccStore->workingSpaceMatrix (cmp.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]} + float wp[3][3] = { + {static_cast(wprof[0][0]), static_cast(wprof[0][1]), static_cast(wprof[0][2])}, + {static_cast(wprof[1][0]), static_cast(wprof[1][1]), static_cast(wprof[1][2])}, + {static_cast(wprof[2][0]), static_cast(wprof[2][1]), static_cast(wprof[2][2])} }; // Conversion rgb -> lab is hard to vectorize because it uses a lut (that's not the main problem) @@ -2166,25 +2170,19 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar LUTu lhist16RETIThr; if(lhist16RETI) { - lhist16RETIThr(32769, 0); + lhist16RETIThr(lhist16RETI.getSize()); lhist16RETIThr.clear(); } #ifdef _OPENMP - #pragma omp for + #pragma omp for schedule(dynamic,16) #endif for (int i = border; i < H - border; i++ ) for (int j = border; j < W - border; j++) { float X, Y, Z, L, aa, bb; - int pos; - float R_, G_, B_; - R_ = red[i][j]; - G_ = green[i][j]; - B_ = blue[i][j]; - float k; //rgb=>lab - Color::rgbxyz(R_, G_, B_, X, Y, Z, wp); + Color::rgbxyz(red[i][j], green[i][j], blue[i][j], X, Y, Z, wp); //convert Lab Color::XYZ2Lab(X, Y, Z, L, aa, bb); conversionBuffer[0][i - border][j - border] = aa; @@ -2195,7 +2193,7 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar // if(R_>40000.f && G_ > 30000.f && B_ > 30000.f) conversionBuffer[3][i - border][j - border] = R_; // else conversionBuffer[3][i - border][j - border] = 0.f; if(lhist16RETI) { - pos = (int) clipretinex(L, 0, 32768); + int pos = L; lhist16RETIThr[pos]++;//histogram in Curve Lab } } @@ -2204,10 +2202,7 @@ void RawImageSource::retinexPrepareBuffers(ColorManagementParams cmp, RetinexPar #pragma omp critical { if(lhist16RETI) { - // Add per Thread LUT to global LUT - for(int i = 0; i < 32769; i++) { - lhist16RETI[i] += lhist16RETIThr[i]; - } + lhist16RETI += lhist16RETIThr; // Add per Thread LUT to global LUT } } #endif @@ -2263,23 +2258,29 @@ void RawImageSource::retinex(ColorManagementParams cmp, RetinexParams deh, ToneC int mode = 0, imax = 0; if(gamm2 < 1.) { - pwr = 1. / pwr; - gamm = 1. / gamm; + std::swap(pwr, gamm); } Color::calcGamma(pwr, ts, mode, imax, g_a0, g_a1, g_a2, g_a3, g_a4, g_a5); // call to calcGamma with selected gamma and slope + double mul = 1. + g_a4; + double add; + double start; + + if(gamm2 < 1.) { + start = g_a3; + add = g_a3; + } else { + add = g_a4; + start = g_a2; + } + // printf("g_a0=%f g_a1=%f g_a2=%f g_a3=%f g_a4=%f\n", g_a0,g_a1,g_a2,g_a3,g_a4); for (int i = 0; i < 65536; i++) { double val = (i) / 65535.; double x; - double mul = 1. + g_a4; - double add = g_a4; - double start = g_a2; if(gamm2 < 1.) { - start = g_a3; - add = g_a3; x = Color::gammareti (val, gamm, start, ts, mul , add); } else { x = Color::igammareti (val, gamm, start, ts, mul , add); @@ -2303,11 +2304,10 @@ void RawImageSource::retinex(ColorManagementParams cmp, RetinexParams deh, ToneC float val; if(dehacontlutili && histLRETI) { - hist16RET(32769, 0); + hist16RET(32768); hist16RET.clear(); histLRETI.clear(); - dLcurve(32769, 0); - dLcurve.clear(); + dLcurve(32768); } FlatCurve* chcurve = NULL;//curve c=f(H) @@ -2335,8 +2335,8 @@ void RawImageSource::retinex(ColorManagementParams cmp, RetinexParams deh, ToneC // one LUT per thread LUTu hist16RETThr; - if(dehacontlutili && histLRETI) { - hist16RETThr(32769, 0); + if(hist16RET) { + hist16RETThr(hist16RET.getSize()); hist16RETThr.clear(); } @@ -2350,7 +2350,7 @@ void RawImageSource::retinex(ColorManagementParams cmp, RetinexParams deh, ToneC LBuffer[i][j] = cdcurve[2.f * temp[i][j]] / 2.f; if(histLRETI) { - int pos = (int) clipretinex(LBuffer[i][j], 0.f, 32768.f); + int pos = LBuffer[i][j]; hist16RETThr[pos]++; //histogram in Curve } } @@ -2363,24 +2363,24 @@ void RawImageSource::retinex(ColorManagementParams cmp, RetinexParams deh, ToneC #pragma omp critical #endif { - if(dehacontlutili && histLRETI) { - // Add per Thread LUT to global LUT - for(int i = 0; i < 32769; i++) { - hist16RET[i] += hist16RETThr[i]; - } + if(hist16RET) { + hist16RET += hist16RETThr; // Add per Thread LUT to global LUT } } } - if(dehacontlutili && histLRETI) {//update histogram + if(hist16RET) {//update histogram + // TODO : When rgbcurvesspeedup branch is merged into master, replace this by the following 1-liner + // hist16RET.compressTo(histLRETI); + // also remove declaration and init of dLcurve some lines above then and finally remove this comment :) for (int i = 0; i < 32768; i++) { val = (double)i / 32767.0; - dLcurve[i] = CLIPD(val); + dLcurve[i] = val; } for (int i = 0; i < 32768; i++) { float hval = dLcurve[i]; - int hi = (int)(255.0f * CLIPD(hval)); + int hi = (int)(255.0f * hval); histLRETI[hi] += hist16RET[i]; } } @@ -2399,7 +2399,6 @@ void RawImageSource::retinex(ColorManagementParams cmp, RetinexParams deh, ToneC for (; j < W - border; j++) { float valp; - float chr; // if(chutili) { // c=f(H) { valp = float((chcurve->getVal(conversionBuffer[3][i - border][j - border]) - 0.5f)); @@ -4365,12 +4364,12 @@ void RawImageSource::hlRecovery (std::string method, float* red, float* green, f void RawImageSource::getAutoExpHistogram (LUTu & histogram, int& histcompr) { -BENCHFUN + BENCHFUN histcompr = 3; histogram(65536 >> histcompr); histogram.clear(); - const float refwb[3] = {static_cast(refwb_red),static_cast(refwb_green),static_cast(refwb_blue)}; + const float refwb[3] = {static_cast(refwb_red), static_cast(refwb_green), static_cast(refwb_blue)}; #ifdef _OPENMP #pragma omp parallel @@ -4388,11 +4387,11 @@ BENCHFUN if (ri->getSensorType() == ST_BAYER) { for (int j = start; j < end; j++) { - tmphistogram[(int)(refwb[ri->FC(i,j)] * rawData[i][j]) >> histcompr] += 4; + tmphistogram[(int)(refwb[ri->FC(i, j)] * rawData[i][j]) >> histcompr] += 4; } } else if (ri->getSensorType() == ST_FUJI_XTRANS) { for (int j = start; j < end; j++) { - tmphistogram[(int)(refwb[ri->XTRANSFC(i,j)] * rawData[i][j]) >> histcompr] += 4; + tmphistogram[(int)(refwb[ri->XTRANSFC(i, j)] * rawData[i][j]) >> histcompr] += 4; } } else if (ri->get_colors() == 1) { for (int j = start; j < end; j++) { @@ -4419,7 +4418,7 @@ BENCHFUN // Histogram MUST be 256 in size; gamma is applied, blackpoint and gain also void RawImageSource::getRAWHistogram (LUTu & histRedRaw, LUTu & histGreenRaw, LUTu & histBlueRaw) { -BENCHFUN + BENCHFUN histRedRaw.clear(); histGreenRaw.clear(); histBlueRaw.clear(); @@ -4429,17 +4428,19 @@ BENCHFUN 65535.0f / ri->get_white(3) }; - const bool fourColours = ri->getSensorType() == ST_BAYER && ((mult[1] != mult[3] || cblacksom[1] != cblacksom[3]) || FC(0,0) == 3 || FC(0,1) == 3 || FC(1,0) == 3 || FC(1,1) == 3); + const bool fourColours = ri->getSensorType() == ST_BAYER && ((mult[1] != mult[3] || cblacksom[1] != cblacksom[3]) || FC(0, 0) == 3 || FC(0, 1) == 3 || FC(1, 0) == 3 || FC(1, 1) == 3); LUTu hist[4]; hist[0](65536); hist[0].clear(); + if (ri->get_colors() > 1) { hist[1](65536); hist[1].clear(); hist[2](65536); hist[2].clear(); } + if (fourColours) { hist[3](65536); hist[3].clear(); @@ -4458,16 +4459,19 @@ BENCHFUN LUTu tmphist[4]; tmphist[0](65536); tmphist[0].clear(); + if (ri->get_colors() > 1) { tmphist[1](65536); tmphist[1].clear(); tmphist[2](65536); tmphist[2].clear(); + if (fourColours) { tmphist[3](65536); tmphist[3].clear(); } } + #ifdef _OPENMP #pragma omp for nowait #endif @@ -4514,9 +4518,11 @@ BENCHFUN #endif { hist[0] += tmphist[0]; + if (ri->get_colors() > 1) { hist[1] += tmphist[1]; hist[2] += tmphist[2]; + if (fourColours) { hist[3] += tmphist[3]; } @@ -4528,13 +4534,16 @@ BENCHFUN int idx; idx = CLIP((int)Color::gamma(mult[0] * (i - (cblacksom[0]/*+black_lev[0]*/)))); histRedRaw[idx >> 8] += hist[0][i]; + if (ri->get_colors() > 1) { idx = CLIP((int)Color::gamma(mult[1] * (i - (cblacksom[1]/*+black_lev[1]*/)))); histGreenRaw[idx >> 8] += hist[1][i]; + if (fourColours) { idx = CLIP((int)Color::gamma(mult[3] * (i - (cblacksom[3]/*+black_lev[3]*/)))); histGreenRaw[idx >> 8] += hist[3][i]; } + idx = CLIP((int)Color::gamma(mult[2] * (i - (cblacksom[2]/*+black_lev[2]*/)))); histBlueRaw[idx >> 8] += hist[2][i]; } @@ -4574,6 +4583,7 @@ void RawImageSource::getAutoWBMultipliers (double &rm, double &gm, double &bm) { BENCHFUN constexpr double clipHigh = 64000.0; + if (ri->get_colors() == 1) { rm = gm = bm = 1; return; @@ -4653,6 +4663,7 @@ void RawImageSource::getAutoWBMultipliers (double &rm, double &gm, double &bm) #ifdef _OPENMP #pragma omp for schedule(dynamic,16) nowait #endif + for (int i = 32; i < H - 32; i++) { for (int j = 32; j < W - 32; j++) { // each loop read 1 rgb triplet value