liz/rawTherapee
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Retinex: some small improvements and some cleanup

Browse Source
This commit is contained in:
heckflosse
2016-05-11 21:04:32 +02:00
parent a742e3a7ba
commit 1265679651
3 changed files with 158 additions and 230 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
rtengine/improccoordinator.cc
View File

@@ -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);

246
rtengine/ipretinex.cc
View File

@@ -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;
}
}

139
rtengine/rawimagesource.cc
View File

@@ -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<float>(wprof[0][0]), static_cast<float>(wprof[0][1]), static_cast<float>(wprof[0][2])},
{static_cast<float>(wprof[1][0]), static_cast<float>(wprof[1][1]), static_cast<float>(wprof[1][2])},
{static_cast<float>(wprof[2][0]), static_cast<float>(wprof[2][1]), static_cast<float>(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<float>(refwb_red),static_cast<float>(refwb_green),static_cast<float>(refwb_blue)};
const float refwb[3] = {static_cast<float>(refwb_red), static_cast<float>(refwb_green), static_cast<float>(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