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

Fixed bad behavior sigma FFTW

Browse Source
This commit is contained in:
Desmis
2019-06-29 15:53:26 +02:00
parent dc1d1a561b
commit fa92226ba3
4 changed files with 32 additions and 20 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
rtengine/improcfun.h
View File

@@ -300,8 +300,8 @@ public:
//locallab
void normalize_mean_dt(float *data, const float *ref, size_t size);
void retinex_pde(float *datain, float * dataout, int bfw, int bfh, float thresh, float multy, float *dE);
void fftw_convol_blur(float *input, float *output, int bfw, int bfh, float radius, int fftkern);
void fftw_convol_blur2(float **input2, float **output2, int bfw, int bfh, float radius, int fftkern);
void fftw_convol_blur(float *input, float *output, int bfw, int bfh, float radius, int fftkern, int algo);
void fftw_convol_blur2(float **input2, float **output2, int bfw, int bfh, float radius, int fftkern, int algo);
void fftw_tile_blur(int GW, int GH, int tilssize , int max_numblox_W, int min_numblox_W, float **tmp1, int numThreads, double radius);
void MSRLocal(int sp, bool fftw, int lum, LabImage * bufreti, LabImage * bufmask, LabImage * buforig, LabImage * buforigmas, float** luminance, float** templ, const float* const *originalLuminance, const int width, const int height, const procparams::LocallabParams &loc, const int skip, const LocretigainCurve &locRETgainCcurve, const int chrome, const int scall, const float krad, float &minCD, float &maxCD, float &mini, float &maxi, float &Tmean, float &Tsigma, float &Tmin, float &Tmax,

4
rtengine/iplocalcontrast.cc
View File

@@ -55,8 +55,8 @@ void ImProcFunctions::localContrast(LabImage *lab, float **destination, const Lo
#endif
gaussianBlur(lab->L, buf, width, height, sigma);
} else {
sigma *= 50.f;//approximation to convert sigma "gaussianBlur" to FFTW blur
ImProcFunctions::fftw_convol_blur2(lab->L, buf, width, height, sigma, 0);
sigma *= 1.f;
ImProcFunctions::fftw_convol_blur2(lab->L, buf, width, height, sigma, 0, 0);
}
#ifdef _OPENMP
#pragma omp parallel for if(multiThread)

38
rtengine/iplocallab.cc
View File

@@ -3994,7 +3994,7 @@ void ImProcFunctions::retinex_pde(float *datain, float * dataout, int bfw, int b
void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int bfh, float radius, int fftkern)
void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int bfh, float radius, int fftkern, int algo)
{
/*
** Jacques Desmis june 2019 - inspired by Copyright 2013 IPOL Image Processing On Line http://www.ipol.im/
@@ -4026,7 +4026,9 @@ void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int
fftwf_plan p;
fftwf_plan pkern;//plan for FFT
int image_size, image_sizechange;
float n_x, n_y;//relative coordonates for kernel Gauss
float n_x = 1.f;
float n_y = 1.f;//relative coordonates for kernel Gauss
float radsig = 1.f;
out = (float*) fftwf_malloc(sizeof(float) * (bfw * bfh));//allocate real datas for FFT
@@ -4043,14 +4045,17 @@ void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int
fftwf_execute(p);
fftwf_destroy_plan(p);
/*define the gaussian constants for the convolution kernel*/
// n_x = rtengine::RT_PI/(double) bfw;//ipol
// n_y = rtengine::RT_PI/(double) bfh;
n_x = 1.f/(float) bfw;//gauss
n_y = 1.f/(float) bfh;
if(algo == 0) {
n_x = rtengine::RT_PI/(double) bfw;//ipol
n_y = rtengine::RT_PI/(double) bfh;
} else if(algo == 1) {
n_x = 1.f/(float) bfw;//gauss
n_y = 1.f/(float) bfh;
radsig = 1.f / (2.f * rtengine::RT_PI * radius * radius);//gauss
}
n_x = n_x * n_x;
n_y = n_y * n_y;
float radsig = 1.f / (2.f * rtengine::RT_PI * radius * radius);//gauss
image_size = bfw * bfh;
image_sizechange = 4 * image_size;
@@ -4061,8 +4066,11 @@ void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int
for(int j = 0; j < bfh; j++){
int index = j * bfw;
for(int i = 0; i < bfw; i++)
// kern[ i+ index] = exp((float)(-radius)*(n_x * i * i + n_y * j * j));//calculate Gauss kernel Ipol formula
kern[ i+ index] = radsig * exp((float)(-(n_x * i * i + n_y * j * j)/ (2.f * radius * radius)));//calculate Gauss kernel with Gauss formula
if(algo == 0) {
kern[ i+ index] = exp((float)(-radius)*(n_x * i * i + n_y * j * j));//calculate Gauss kernel Ipol formula
} else if(algo == 1) {
kern[ i+ index] = radsig * exp((float)(-(n_x * i * i + n_y * j * j)/ (2.f * radius * radius)));//calculate Gauss kernel with Gauss formula
}
}
/*compute the Fourier transform of the kernel data*/
pkern = fftwf_plan_r2r_2d(bfh, bfw, kern, outkern, FFTW_REDFT10, FFTW_REDFT10,FFTW_ESTIMATE);//FFT 2 dimensions forward
@@ -4088,7 +4096,11 @@ void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int
for(int j = 0; j < bfh; j++){
int index = j * bfw;
for(int i = 0; i < bfw; i++)
out[i + index] *= exp((float)(-radius)*(n_x * i * i + n_y * j * j));//apply Gauss kernel whithout FFT
if(algo == 0) {
out[i + index] *= exp((float)(-radius)*(n_x * i * i + n_y * j * j));//apply Gauss kernel whithout FFT
} else if(algo == 1) {
out[i + index] *= radsig * exp((float)(-(n_x * i * i + n_y * j * j)/ (2.f * radius * radius)));//calculate Gauss kernel with Gauss formula
}
}
}
@@ -4107,7 +4119,7 @@ void ImProcFunctions::fftw_convol_blur(float *input, float *output, int bfw, int
}
}
void ImProcFunctions::fftw_convol_blur2(float **input2, float **output2, int bfw, int bfh, float radius, int fftkern)
void ImProcFunctions::fftw_convol_blur2(float **input2, float **output2, int bfw, int bfh, float radius, int fftkern, int algo)
{
MyMutex::MyLock lock(*fftwMutex);
@@ -4131,7 +4143,7 @@ void ImProcFunctions::fftw_convol_blur2(float **input2, float **output2, int bfw
input[y * bfw + x] = input2[y][x];
}
}
ImProcFunctions::fftw_convol_blur(input, output, bfw, bfh, radius, fftkern);
ImProcFunctions::fftw_convol_blur(input, output, bfw, bfh, radius, fftkern, algo);
#ifdef _OPENMP
#pragma omp parallel for schedule(dynamic,16)

6
rtengine/ipretinex.cc
View File

@@ -961,7 +961,7 @@ void ImProcFunctions::MSRLocal(int sp, bool fftw, int lum, LabImage * bufreti, L
}
float *buffer = new float[W_L * H_L];
float mulradiusfftw = 40.f;
float mulradiusfftw = 1.f;
for (int scale = scal - 1; scale >= 0; scale--) {
// printf("retscale=%f scale=%i \n", mulradiusfftw * RetinexScales[scale], scale);
if(!fftw) {
@@ -982,11 +982,11 @@ void ImProcFunctions::MSRLocal(int sp, bool fftw, int lum, LabImage * bufreti, L
} else {
if (scale == scal - 1)
{
ImProcFunctions::fftw_convol_blur2(src, out, W_L, H_L, mulradiusfftw * RetinexScales[scale], 0);
ImProcFunctions::fftw_convol_blur2(src, out, W_L, H_L, mulradiusfftw * RetinexScales[scale], 0, 0);
} else // reuse result of last iteration
{
// out was modified in last iteration => restore it
ImProcFunctions::fftw_convol_blur2(out, out, W_L, H_L,sqrtf(SQR(mulradiusfftw * RetinexScales[scale]) - SQR(mulradiusfftw * RetinexScales[scale + 1])), 0);
ImProcFunctions::fftw_convol_blur2(out, out, W_L, H_L,sqrtf(SQR(mulradiusfftw * RetinexScales[scale]) - SQR(mulradiusfftw * RetinexScales[scale + 1])), 0, 0);
}
}