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

Optimization for Flat Field, Issue 2262

Browse Source
This commit is contained in:
Ingo
2014-02-27 15:24:37 +01:00
parent 6c089a8654
commit 6014c2fde9
1 changed files with 113 additions and 21 deletions
Download Patch File Download Diff File Expand all files Collapse all files

134
rtengine/rawimagesource.cc
View File

@@ -39,6 +39,7 @@
#ifdef _OPENMP #ifdef _OPENMP
#include <omp.h> #include <omp.h>
#endif #endif
#include "opthelper.h"
namespace rtengine { namespace rtengine {
@@ -995,7 +996,7 @@ void RawImageSource::preprocess (const RAWParams &raw, const LensProfParams &le
printf( "Subtracting Darkframe:%s\n",rid->get_filename().c_str()); printf( "Subtracting Darkframe:%s\n",rid->get_filename().c_str());
} }
//copyOriginalPixels(ri, rid); //copyOriginalPixels(ri, rid);
//FLATFIELD start //FLATFIELD start
Glib::ustring newFF = raw.ff_file; Glib::ustring newFF = raw.ff_file;
RawImage *rif=NULL; RawImage *rif=NULL;
@@ -1273,7 +1274,7 @@ void RawImageSource::copyOriginalPixels(const RAWParams &raw, RawImage *src, Raw
else //(raw.ff_BlurType == RAWParams::ff_BlurTypestring[RAWParams::area_ff]) else //(raw.ff_BlurType == RAWParams::ff_BlurTypestring[RAWParams::area_ff])
cfaboxblur(riFlatFile, cfablur, BS, BS); cfaboxblur(riFlatFile, cfablur, BS, BS);
float refcolor[2][2],vignettecorr; float refcolor[2][2];
//find center ave values by channel //find center ave values by channel
for (int m=0; m<2; m++) for (int m=0; m<2; m++)
for (int n=0; n<2; n++) { for (int n=0; n<2; n++) {
@@ -1286,13 +1287,17 @@ void RawImageSource::copyOriginalPixels(const RAWParams &raw, RawImage *src, Raw
for (int m=0; m<2; m++) for (int m=0; m<2; m++)
for (int n=0; n<2; n++) { for (int n=0; n<2; n++) {
for (int row = 0; row+m < H; row+=2) #pragma omp parallel for
for (int col = 0; col+n < W; col+=2) { for (int row = 0; row< H-m; row+=2) {
int c = FC(row, col); int c = FC(row, 0);
int c4 = ( c == 1 && !(row&1) ) ? 3 : c; int c4 = ( c == 1 && !(row&1) ) ? 3 : c;
vignettecorr = ( refcolor[m][n]/max(1e-5f,cfablur[(row+m)*W+col+n]-black[c4]) ); for (int col = 0; col < W-n; col+=2) {
// int c = FC(row, col);
// int c4 = ( c == 1 && !(row&1) ) ? 3 : c;
float vignettecorr = ( refcolor[m][n]/max(1e-5f,cfablur[(row+m)*W+col+n]-black[c4]) );
rawData[row+m][col+n] = (rawData[row+m][col+n]-black[c4]) * vignettecorr + black[c4]; rawData[row+m][col+n] = (rawData[row+m][col+n]-black[c4]) * vignettecorr + black[c4];
} }
}
} }
if (raw.ff_BlurType == RAWParams::ff_BlurTypestring[RAWParams::vh_ff]) { if (raw.ff_BlurType == RAWParams::ff_BlurTypestring[RAWParams::vh_ff]) {
@@ -1304,18 +1309,18 @@ void RawImageSource::copyOriginalPixels(const RAWParams &raw, RawImage *src, Raw
cfaboxblur(riFlatFile, cfablur1, 0, 2*BS);//now do horizontal blur cfaboxblur(riFlatFile, cfablur1, 0, 2*BS);//now do horizontal blur
cfaboxblur(riFlatFile, cfablur2, 2*BS, 0);//now do vertical blur cfaboxblur(riFlatFile, cfablur2, 2*BS, 0);//now do vertical blur
float vlinecorr, hlinecorr;
for (int m=0; m<2; m++) for (int m=0; m<2; m++)
for (int n=0; n<2; n++) { for (int n=0; n<2; n++) {
for (int row = 0; row+m < H; row+=2) #pragma omp parallel for
for (int col = 0; col+n < W; col+=2) { for (int row = 0; row< H-m; row+=2) {
int c = FC(row, col); int c = FC(row, 0);
int c4 = ( c == 1 && !(row&1) ) ? 3 : c; int c4 = ( c == 1 && !(row&1) ) ? 3 : c;
hlinecorr = (max(1e-5f,cfablur[(row+m)*W+col+n]-black[c4])/max(1e-5f,cfablur1[(row+m)*W+col+n]-black[c4]) ); for (int col = 0; col < W-n; col+=2) {
vlinecorr = (max(1e-5f,cfablur[(row+m)*W+col+n]-black[c4])/max(1e-5f,cfablur2[(row+m)*W+col+n]-black[c4]) ); float hlinecorr = (max(1e-5f,cfablur[(row+m)*W+col+n]-black[c4])/max(1e-5f,cfablur1[(row+m)*W+col+n]-black[c4]) );
float vlinecorr = (max(1e-5f,cfablur[(row+m)*W+col+n]-black[c4])/max(1e-5f,cfablur2[(row+m)*W+col+n]-black[c4]) );
rawData[row+m][col+n] = ((rawData[row+m][col+n]-black[c4]) * hlinecorr * vlinecorr + black[c4]); rawData[row+m][col+n] = ((rawData[row+m][col+n]-black[c4]) * hlinecorr * vlinecorr + black[c4]);
} }
}
} }
free (cfablur1); free (cfablur1);
free (cfablur2); free (cfablur2);
@@ -1354,15 +1359,18 @@ void RawImageSource::copyOriginalPixels(const RAWParams &raw, RawImage *src, Raw
} }
} }
void RawImageSource::cfaboxblur(RawImage *riFlatFile, float* cfablur, int boxH, int boxW ) { SSEFUNCTION void RawImageSource::cfaboxblur(RawImage *riFlatFile, float* cfablur, const int boxH, const int boxW ) {
float (*cfatmp); float (*cfatmp);
cfatmp = (float (*)) calloc (H*W, sizeof *cfatmp); cfatmp = (float (*)) calloc (H*W, sizeof *cfatmp);
float hotdeadthresh = 0.5; // const float hotdeadthresh = 0.5;
#pragma omp parallel
{
#pragma omp for
for (int i=0; i<H; i++) { for (int i=0; i<H; i++) {
int iprev,inext,jprev,jnext; int iprev,inext,jprev,jnext;
float p[9],temp, median; int p[5],temp, median;
if (i<2) {iprev=i+2;} else {iprev=i-2;} if (i<2) {iprev=i+2;} else {iprev=i-2;}
if (i>H-3) {inext=i-2;} else {inext=i+2;} if (i>H-3) {inext=i-2;} else {inext=i+2;}
for (int j=0; j<W; j++) { for (int j=0; j<W; j++) {
@@ -1373,7 +1381,9 @@ void RawImageSource::cfaboxblur(RawImage *riFlatFile, float* cfablur, int boxH,
// riFlatFile->data[inext][jprev], riFlatFile->data[inext][j], riFlatFile->data[inext][jnext], cfatmp[i*W+j]); // riFlatFile->data[inext][jprev], riFlatFile->data[inext][j], riFlatFile->data[inext][jnext], cfatmp[i*W+j]);
med5(riFlatFile->data[iprev][j], riFlatFile->data[i][jprev],riFlatFile->data[i][j], med5(riFlatFile->data[iprev][j], riFlatFile->data[i][jprev],riFlatFile->data[i][j],
riFlatFile->data[i][jnext], riFlatFile->data[inext][j],median); riFlatFile->data[i][jnext], riFlatFile->data[inext][j],median);
if (riFlatFile->data[i][j]>hotdeadthresh*median || median>hotdeadthresh*riFlatFile->data[i][j]) {
// if (riFlatFile->data[i][j]>hotdeadthresh*median || median>hotdeadthresh*riFlatFile->data[i][j]) {
if (((int)riFlatFile->data[i][j]<<1)>median || (median<<1) > riFlatFile->data[i][j]) {
cfatmp[i*W+j] = median; cfatmp[i*W+j] = median;
} else { } else {
cfatmp[i*W+j] = riFlatFile->data[i][j]; cfatmp[i*W+j] = riFlatFile->data[i][j];
@@ -1384,6 +1394,7 @@ void RawImageSource::cfaboxblur(RawImage *riFlatFile, float* cfablur, int boxH,
//box blur cfa image; box size = BS //box blur cfa image; box size = BS
//horizontal blur //horizontal blur
#pragma omp for
for (int row = 0; row < H; row++) { for (int row = 0; row < H; row++) {
int len = boxW/2 + 1; int len = boxW/2 + 1;
cfatmp[row*W+0] = cfatmp[row*W+0]/len; cfatmp[row*W+0] = cfatmp[row*W+0]/len;
@@ -1410,6 +1421,86 @@ void RawImageSource::cfaboxblur(RawImage *riFlatFile, float* cfablur, int boxH,
} }
//vertical blur //vertical blur
#ifdef __SSE2__
__m128 leninitv = _mm_set1_ps( (float)((int)(boxH/2 + 1)));
__m128 onev = _mm_set1_ps( 1.0f );
__m128 temp1v,temp2v,lenv,lenp1v,lenm1v;
int row;
#pragma omp for
for (int col = 0; col < W-3; col+=4) {
lenv = leninitv;
temp1v = LVFU(cfatmp[0*W+col]) / lenv;
temp2v = LVFU(cfatmp[1*W+col]) / lenv;
for (int i=2; i<boxH+2; i+=2) {
temp1v += LVFU(cfatmp[i*W+col]) / lenv;
temp2v += LVFU(cfatmp[(i+1)*W+col]) / lenv;
}
_mm_storeu_ps(&cfablur[0*W+col], temp1v);
_mm_storeu_ps(&cfablur[1*W+col], temp2v);
for (row=2; row<boxH+2; row+=2) {
lenp1v = lenv + onev;
temp1v = (temp1v * lenv + LVFU(cfatmp[(row+boxH)*W+col])) / lenp1v;
temp2v = (temp2v * lenv + LVFU(cfatmp[(row+boxH+1)*W+col])) / lenp1v;
_mm_storeu_ps( &cfablur[row*W+col], temp1v);
_mm_storeu_ps( &cfablur[(row+1)*W+col], temp2v);
lenv = lenp1v;
}
for (; row < H-boxH-1; row+=2) {
temp1v = temp1v + (LVFU(cfatmp[(row+boxH)*W+col]) - LVFU(cfatmp[(row-boxH-2)*W+col]))/lenv;
temp2v = temp2v + (LVFU(cfatmp[(row+1+boxH)*W+col]) - LVFU(cfatmp[(row+1-boxH-2)*W+col]))/lenv;
_mm_storeu_ps(&cfablur[row*W+col], temp1v);
_mm_storeu_ps(&cfablur[(row+1)*W+col], temp2v);
}
for(; row < H-boxH; row++) {
temp1v = temp1v + (LVFU(cfatmp[(row+boxH)*W+col]) - LVFU(cfatmp[(row-boxH-2)*W+col]))/lenv;
_mm_storeu_ps(&cfablur[row*W+col], temp1v);
__m128 swapv = temp1v;
temp1v = temp2v;
temp2v = swapv;
}
for (; row<H-1; row+=2) {
lenm1v = lenv - onev;
temp1v = (temp1v * lenv - LVFU(cfatmp[(row-boxH-2)*W+col])) / lenm1v;
temp2v = (temp2v * lenv - LVFU(cfatmp[(row-boxH-1)*W+col])) / lenm1v;
_mm_storeu_ps(&cfablur[row*W+col], temp1v);
_mm_storeu_ps(&cfablur[(row+1)*W+col], temp2v);
lenv = lenm1v;
}
for(; row < H; row++) {
lenm1v = lenv - onev;
temp1v = (temp1v * lenv - LVFU(cfatmp[(row-boxH-2)*W+col])) / lenm1v;
_mm_storeu_ps(&cfablur[(row)*W+col], temp1v);
}
}
for (int col = W-(W%4); col < W; col++) {
int len = boxH/2 + 1;
cfablur[0*W+col] = cfatmp[0*W+col]/len;
cfablur[1*W+col] = cfatmp[1*W+col]/len;
for (int i=2; i<boxH+2; i+=2) {
cfablur[0*W+col] += cfatmp[i*W+col]/len;
cfablur[1*W+col] += cfatmp[(i+1)*W+col]/len;
}
for (int row=2; row<boxH+2; row+=2) {
cfablur[row*W+col] = (cfablur[(row-2)*W+col]*len + cfatmp[(row+boxH)*W+col])/(len+1);
cfablur[(row+1)*W+col] = (cfablur[(row-1)*W+col]*len + cfatmp[(row+boxH+1)*W+col])/(len+1);
len ++;
}
for (int row = boxH+2; row < H-boxH; row++) {
cfablur[row*W+col] = cfablur[(row-2)*W+col] + (cfatmp[(row+boxH)*W+col] - cfatmp[(row-boxH-2)*W+col])/len;
}
for (int row=H-boxH; row<H; row+=2) {
cfablur[row*W+col] = (cfablur[(row-2)*W+col]*len - cfatmp[(row-boxH-2)*W+col])/(len-1);
if (row+1<H)
cfablur[(row+1)*W+col] = (cfablur[(row-1)*W+col]*len - cfatmp[(row-boxH-1)*W+col])/(len-1);
len --;
}
}
#else
#pragma omp for
for (int col = 0; col < W; col++) { for (int col = 0; col < W; col++) {
int len = boxH/2 + 1; int len = boxH/2 + 1;
cfablur[0*W+col] = cfatmp[0*W+col]/len; cfablur[0*W+col] = cfatmp[0*W+col]/len;
@@ -1433,8 +1524,9 @@ void RawImageSource::cfaboxblur(RawImage *riFlatFile, float* cfablur, int boxH,
len --; len --;
} }
} }
#endif
}
free (cfatmp); free (cfatmp);
} }