From cd91e7890e526d170c3396962a762cdf1041e410 Mon Sep 17 00:00:00 2001 From: gabor Date: Fri, 21 May 2010 16:14:18 +0200 Subject: [PATCH] ImProcFunctions cleanup and transition to OpenMP -- phase 1 --- CMakeLists.txt | 5 + rtengine/CMakeLists.txt | 5 +- rtengine/bilateral2.h | 458 ++---- rtengine/dcrop.cc | 18 +- rtengine/gauss.cc | 49 - rtengine/gauss.h | 480 +----- rtengine/improccoordinator.cc | 21 +- rtengine/improcfun.cc | 1929 +----------------------- rtengine/improcfun.h | 86 +- rtengine/iplab2rgb.cc | 259 ++++ rtengine/ipresize.cc | 293 ++++ rtengine/ipsharpen.cc | 201 +++ rtengine/iptransform.cc | 821 ++++++++++ rtengine/labimage.cc | 42 + rtengine/{bilateral2.cc => labimage.h} | 36 +- rtengine/rtthumbnail.cc | 16 +- rtengine/shmap.cc | 61 +- rtengine/shmap.h | 3 +- rtengine/simpleprocess.cc | 17 +- 19 files changed, 1944 insertions(+), 2856 deletions(-) delete mode 100644 rtengine/gauss.cc create mode 100644 rtengine/iplab2rgb.cc create mode 100644 rtengine/ipresize.cc create mode 100644 rtengine/ipsharpen.cc create mode 100644 rtengine/iptransform.cc create mode 100644 rtengine/labimage.cc rename rtengine/{bilateral2.cc => labimage.h} (54%) diff --git a/CMakeLists.txt b/CMakeLists.txt index a67151f21..5e6abd6c3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -116,6 +116,11 @@ if (WITH_RAWZOR) endif (WIN32) endif (WITH_RAWZOR) +find_package(OpenMP) +if (OPENMP_FOUND) + set(CMAKE_CXX_FLAGS ${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}) +endif (OPENMP_FOUND) + add_subdirectory (rtexif) add_subdirectory (rtengine) add_subdirectory (rtgui) diff --git a/rtengine/CMakeLists.txt b/rtengine/CMakeLists.txt index f27265f7a..d8f51a125 100644 --- a/rtengine/CMakeLists.txt +++ b/rtengine/CMakeLists.txt @@ -6,11 +6,12 @@ link_directories (${CMAKE_CURRENT_SOURCE_DIR}/../rtexif ${EXTRA_LIBDIR} ${GTHREA ${GOBJECT_LIBRARY_DIRS} ${GLIB2_LIBRARY_DIRS} ${GLIBMM_LIBRARY_DIRS} ${IPTCDATA_LIBRARY_DIRS} ${LCMS_LIBRARY_DIRS}) -set (RTENGINESOURCEFILES colortemp.cc curves.cc dcraw.cc gauss.cc iccstore.cc +set (RTENGINESOURCEFILES colortemp.cc curves.cc dcraw.cc iccstore.cc image8.cc image16.cc imagedata.cc imageio.cc improcfun.cc init.cc dcrop.cc loadinitial.cc procparams.cc rawimagesource.cc shmap.cc simpleprocess.cc refreshmap.cc stdimagesource.cc myfile.cc iccjpeg.c hlmultipliers.cc improccoordinator.cc - processingjob.cc rtthumbnail.cc utils.cc bilateral2.cc) + processingjob.cc rtthumbnail.cc utils.cc labimage.cc + iplab2rgb.cc ipsharpen.cc iptransform.cc ipresize.cc) add_library (rtengine ${RTENGINESOURCEFILES}) #It may be nice to store library version too diff --git a/rtengine/bilateral2.h b/rtengine/bilateral2.h index 531c248cc..0c59f9db5 100644 --- a/rtengine/bilateral2.h +++ b/rtengine/bilateral2.h @@ -26,6 +26,9 @@ #include #include #include +#include + +// This seems ugly, but way faster than any other solutions I tried #define ELEM(a,b) (src[i - a][j - b] * ec[src[i - a][j - b]-src[i][j]+0x10000]) #define SULY(a,b) (ec[src[i - a][j - b]-src[i][j]+0x10000]) @@ -34,22 +37,22 @@ int* ec = new int [0x20000]; \ for (int i=0; i<0x20000; i++) \ ec[i] = (int)(exp(-(double)(i-0x10000)*(double)(i-0x10000) / (2.0*sens*sens))*scale); \ - int start = row_from; \ - if (start<(b)) start = (b); \ - int end = row_to; \ - if (end>H-(b)) end = H-(b); \ - for (int i=start; i=end || j>=W-(b)) \ + if (i=rend || j>=cend) \ dst[i][j] = src[i][j]; \ else \ dst[i][j] = buffer[i][j]; -#define BL_OPER3(a11,a12,a21,a22) A v = a11*ELEM(-1,-1) + a12*ELEM(-1,0) + a11*ELEM(-1,1) + \ +#define BL_OPER3(a11,a12,a21,a22) for (int i=rstart; i void bilateral05 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral05 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(318,1) + #pragma omp parallel for if (multiThread) BL_OPER3(1,7,7,55) BL_END(1) } // sigma = 0.6 -template void bilateral06 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral06 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(768,1) + #pragma omp parallel for if (multiThread) BL_OPER3(1,4,4,16) BL_END(1) } // sigma = 0.7 -template void bilateral07 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral07 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(366,2) + #pragma omp parallel for if (multiThread) BL_OPER5(0,0,1,0,8,21,1,21,59) BL_END(2) } // sigma = 0.8 -template void bilateral08 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral08 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(753,2) + #pragma omp parallel for if (multiThread) BL_OPER5(0,0,1,0,5,10,1,10,23) BL_END(2) } // sigma = 0.9 -template void bilateral09 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral09 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(595,2) + #pragma omp parallel for if (multiThread) BL_OPER5(0,1,2,1,6,12,2,12,22) BL_END(2) } // sigma = 1.0 -template void bilateral10 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral10 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(910,2) + #pragma omp parallel for if (multiThread) BL_OPER5(0,1,2,1,4,7,2,7,12) BL_END(2) } // sigma = 1.1 -template void bilateral11 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral11 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(209,3) + #pragma omp parallel for if (multiThread) BL_OPER7(0,0,1,1,0,2,5,8,1,5,18,27,1,8,27,41) BL_END(3) } // sigma = 1.2 -template void bilateral12 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral12 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(322,3) + #pragma omp parallel for if (multiThread) BL_OPER7(0,0,1,1,0,1,4,6,1,4,11,16,1,6,16,23) BL_END(3) } // sigma = 1.3 -template void bilateral13 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral13 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(336,3) + #pragma omp parallel for if (multiThread) BL_OPER7(0,0,1,1,0,2,4,6,1,4,11,14,1,6,14,19) BL_END(3) } // sigma = 1.4 -template void bilateral14 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral14 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(195,3) + #pragma omp parallel for if (multiThread) BL_OPER7(0,1,2,3,1,4,8,10,2,8,17,21,3,10,21,28) BL_END(3) } // sigma = 1.5 -template void bilateral15 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral15 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(132,4) + #pragma omp parallel for if (multiThread) BL_OPER9(0,0,0,1,1,0,1,2,4,5,0,2,6,12,14,1,4,12,22,28,1,5,14,28,35) BL_END(4) } // sigma = 1.6 -template void bilateral16 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral16 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(180,4) + #pragma omp parallel for if (multiThread) BL_OPER9(0,0,0,1,1,0,1,2,3,4,0,2,5,9,10,1,3,9,15,19,1,4,10,19,23) BL_END(4) } // sigma = 1.7 -template void bilateral17 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral17 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(195,4) + #pragma omp parallel for if (multiThread) BL_OPER9(0,0,1,1,1,0,1,2,3,4,1,2,5,8,9,1,3,8,13,16,1,4,9,16,19) BL_END(4) } // sigma = 1.8 -template void bilateral18 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral18 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(151,4) + #pragma omp parallel for if (multiThread) BL_OPER9(0,0,1,2,2,0,1,3,5,5,1,3,6,10,12,2,5,10,16,19,2,5,12,19,22) BL_END(4) } // sigma = 1.9 -template void bilateral19 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral19 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(151,4) + #pragma omp parallel for if (multiThread) BL_OPER9(0,0,1,2,2,0,1,3,4,5,1,3,5,8,9,2,4,8,12,14,2,5,9,14,16) BL_END(4) } // sigma = 2 -template void bilateral20 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral20 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(116,5) + #pragma omp parallel for if (multiThread) BL_OPER11(0,0,0,1,1,1,0,0,1,2,3,3,0,1,2,4,7,7,1,2,4,8,12,14,1,3,7,12,18,20,1,3,7,14,20,23) BL_END(5) } // sigma = 2.1 -template void bilateral21 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral21 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(127,5) + #pragma omp parallel for if (multiThread) BL_OPER11(0,0,0,1,1,1,0,0,1,2,3,3,0,1,2,4,6,7,1,2,4,8,11,12,1,3,6,11,15,17,1,3,7,12,17,19) BL_END(5) } // sigma = 2.2 -template void bilateral22 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral22 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(109,5) + #pragma omp parallel for if (multiThread) BL_OPER11(0,0,0,1,1,2,0,1,2,3,3,4,1,2,3,5,7,8,1,3,5,9,12,13,1,3,7,12,16,18,2,4,8,13,18,20) BL_END(5) } // sigma = 2.3 -template void bilateral23 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral23 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(132,5) + #pragma omp parallel for if (multiThread) BL_OPER11(0,0,1,1,1,1,0,1,1,2,3,3,1,1,3,5,6,7,1,2,5,7,10,11,1,3,6,10,13,14,1,3,7,11,14,16) BL_END(5) } // sigma = 2.4 -template void bilateral24 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral24 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { BL_BEGIN(156,5) + #pragma omp parallel for if (multiThread) BL_OPER11(0,0,1,1,1,1,0,1,1,2,3,3,1,1,3,4,5,6,1,2,4,6,8,9,1,3,5,8,10,11,1,3,6,9,11,12) BL_END(5) } // sigma = 2.5 -template void bilateral25 (T** src, T** dst, T** buffer, int W, int H, int row_from, int row_to, double sens) { +template void bilateral25 (T** src, T** dst, T** buffer, int W, int H, double sens, bool multiThread) { - BL_BEGIN(173,5) + BL_BEGIN(173,5) + #pragma omp parallel for if (multiThread) BL_OPER11(0,0,1,1,1,1,0,1,1,2,3,3,1,1,2,4,5,5,1,2,4,5,7,7,1,3,5,7,9,9,1,3,5,7,9,10) BL_END(5) } -class Dim { - - public: - int W, H, row_from, row_to; - - Dim (int w, int h, int rf, int rt) : W(w), H(h), row_from(rf), row_to(rt) {} - -}; - // main bilateral filter -template void bilateral (T** src, T** dst, T** buffer, Dim dim, double sigma, double sens) { - - int W = dim.W; - int H = dim.H; - int row_from = dim.row_from; - int row_to = dim.row_to; +template void bilateral (T** src, T** dst, T** buffer, int W, int H, double sigma, double sens, bool multiThread) { if (sigma<0.45) - for (int i=row_from; i (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral05 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<0.65) - bilateral06 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral06 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<0.75) - bilateral07 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral07 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<0.85) - bilateral08 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral08 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<0.95) - bilateral09 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral09 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.05) - bilateral10 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral10 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.15) - bilateral11 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral11 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.25) - bilateral12 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral12 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.35) - bilateral13 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral13 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.45) - bilateral14 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral14 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.55) - bilateral15 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral15 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.65) - bilateral16 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral16 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.75) - bilateral17 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral17 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.85) - bilateral18 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral18 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<1.95) - bilateral19 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral19 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<2.05) - bilateral20 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral20 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<2.15) - bilateral21 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral21 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<2.25) - bilateral22 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral22 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<2.35) - bilateral23 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral23 (src, dst, buffer, W, H, sens, multiThread); else if (sigma<2.45) - bilateral24 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral24 (src, dst, buffer, W, H, sens, multiThread); else - bilateral25 (src, dst, buffer, W, H, row_from, row_to, sens); + bilateral25 (src, dst, buffer, W, H, sens, multiThread); } -void bilateral_unsigned (unsigned short** src, unsigned short** dst, unsigned short** buffer, Dim dim, double sigma, double sens); -void bilateral_signed (short** src, short** dst, short** buffer, Dim dim, double sigma, double sens); - -/* -template void bilateral (T** src, int** dst, int W, int H, int sigmar, double sigmas) { - - time_t t1 = clock (); - - int r = 2.0*sigmas + 0.5; - - double scaleg = 1024.0; - - double MAXINT = 65536.0*65536.0; - double sgmax = 275000/(MAXINT/2.0/sigmar/sigmar+(r+1)*(r+1)/sigmas/sigmas); - printf ("sgmax = %lf\n", sgmax); - if (scaleg>sgmax) - scaleg = sgmax; - - // kernel vector for the spatial gaussian filter * 1024 - int* kernel = new int[1+2*r]; - for (int i=0; i<2*r+1; i++) { - kernel[i] = scaleg / (2.0*sigmas*sigmas) * (i-r)*(i-r) + 0.25; - printf ("kerneli = %d\n", kernel[i]); - } - - // exponential lookup table - int scale = (2.0*sigmar*sigmar) / scaleg; - int scalem = 65535/((1+2*r)*(1+2*r)); - int ec[256000]; - for (int i=0; i<256000; i++) - ec[i] = exp (-i/scaleg) * scalem; - - for (int i=r; i0?((a) void bilateral (T** src, T** dst, AlignedBuffer* buffer, int W, int H, double sigmar, double sigmas) { - - time_t t1 = clock (); - - float alpha = 0.5 / sigmar / sigmar; - - // buffer for the final image - float** buff_final = new float*[H]; - for (int i=0; i void bilateral (T** src, T** dst, int W, int H, int sigmar, double sigmas, int row_from, int row_to) { // range weights @@ -654,95 +469,4 @@ template void bilateral (T** src, T** dst, int W, int H, int sigmar, do delete [] buff_final[i]; } -struct bilateralparams { - - int row_from, row_to; -}; - -void bilateral_box_unsigned (unsigned short** src, unsigned short** dst, int W, int H, int sigmar, double sigmas, bilateralparams row); - -/* -#define BINBIT 12 -#define TRANSBIT 4 - -#define ADDHIST(a,b) { for (int k=0; k<(1< void bilateral (T** src, T** dst, AlignedBuffer* buffer, int W, int H, int sigmar, double sigmas) { - - time_t t1 = clock (); - - unsigned short** rowhist = new unsigned short* [H]; - for (int i=0; i>TRANSBIT]++; - - // let the game begin... - for (int j=r+1; j>TRANSBIT]--; - rowhist[i][src[i][j+r]>>TRANSBIT]++; - } - // sum up upper histograms - memset (hist, 0, (1<>TRANSBIT); k++) { - float w = 1.0 - (double)k/(sigmar>>TRANSBIT); - int v = src[i][j]>>TRANSBIT; - if (v-k >= 0) { - weight += hist [v-k] * w; - buff_final[i][j] += hist [v-k] * w * (src[i][j]-(k<H-r-1 || jW-r-1) - dst[i][j] = src[i][j]; - else - dst[i][j] = (T)CLIP(buff_final[i][j]); - - delete [] hist; - for (int i=0; iipf.setScale (skip); + t1.set (); // give possibility to the listener to modify crop window (as the full image dimensions are already known at this point) int wx, wy, ww, wh, ws; @@ -130,7 +132,7 @@ void Crop::update (int todo, bool internal) { } if (!resizeCrop) resizeCrop = new Image16 (rcw, rch); - parent->ipf.resize (origCrop, resizeCrop, params.resize); + parent->ipf.resize (origCrop, resizeCrop); baseCrop = resizeCrop; } parent->minit.unlock (); @@ -173,7 +175,7 @@ void Crop::update (int todo, bool internal) { if (settings->verbose) printf ("C-BLURMAP: %d\n", t4.etime(t3)); if (todo & M_RGBCURVE) { - parent->ipf.rgbProc (baseCrop, laboCrop, ¶ms, parent->tonecurve, cshmap); + parent->ipf.rgbProc (baseCrop, laboCrop, parent->tonecurve, cshmap); } t5.set (); @@ -182,8 +184,8 @@ void Crop::update (int todo, bool internal) { if (todo & M_LUMINANCE) { parent->ipf.luminanceCurve (laboCrop, labnCrop, parent->lumacurve, 0, croph); if (skip==1) { - parent->ipf.lumadenoise (labnCrop, ¶ms, 1, cbuffer); - parent->ipf.sharpening (labnCrop, ¶ms, 1, (unsigned short**)cbuffer); + parent->ipf.lumadenoise (labnCrop, cbuffer); + parent->ipf.sharpening (labnCrop, (unsigned short**)cbuffer); } } @@ -191,9 +193,9 @@ void Crop::update (int todo, bool internal) { if (settings->verbose) printf ("C-LUMINANCE: %d\n", t6.etime(t5)); if (todo & M_COLOR) { - parent->ipf.colorCurve (laboCrop, labnCrop, ¶ms); + parent->ipf.colorCurve (laboCrop, labnCrop); if (skip==1) - parent->ipf.colordenoise (labnCrop, ¶ms, 1, cbuffer); + parent->ipf.colordenoise (labnCrop, cbuffer); } t7.set (); @@ -290,7 +292,7 @@ if (settings->verbose) printf ("setcropsizes before lock\n"); // determine which part of the source image is required to compute the crop rectangle int orx, ory, orw, orh; ProcParams& params = parent->params; - parent->ipf.transCoord (¶ms, parent->fw, parent->fh, bx1, by1, bw, bh, orx, ory, orw, orh); + ImProcFunctions::transCoord (¶ms, parent->fw, parent->fh, bx1, by1, bw, bh, orx, ory, orw, orh); int tr = TR_NONE; if (params.coarse.rotate==90) tr |= TR_R90; @@ -324,7 +326,7 @@ if (settings->verbose) printf ("setcropsizes before lock\n"); laboCrop = new LabImage (cropw, croph); labnCrop = new LabImage (cropw, croph); cropImg = new Image8 (cropw, croph); - cshmap = new SHMap (cropw, croph); + cshmap = new SHMap (cropw, croph, true); cbuffer = new int*[croph]; for (int i=0; i - * - * 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 - -void gaussHorizontal_unsigned (unsigned short** src, unsigned short** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma) { - - gaussHorizontal (src, dst, buffer, W, row_from, row_to, sigma); -} - -void gaussVertical_unsigned (unsigned short** src, unsigned short** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma) { - - gaussVertical (src, dst, buffer, H, col_from, col_to, sigma); -} - -void gaussHorizontal_signed (short** src, short** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma) { - - gaussHorizontal (src, dst, buffer, W, row_from, row_to, sigma); -} - -void gaussVertical_signed (short** src, short** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma) { - - gaussVertical (src, dst, buffer, H, col_from, col_to, sigma); -} - -void gaussHorizontal_float (float** src, float** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma) { - - gaussHorizontal (src, dst, buffer, W, row_from, row_to, sigma); -} - -void gaussVertical_float (float** src, float** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma) { - - gaussVertical (src, dst, buffer, H, col_from, col_to, sigma); -} diff --git a/rtengine/gauss.h b/rtengine/gauss.h index 3e7fe2252..c9e25a0b0 100644 --- a/rtengine/gauss.h +++ b/rtengine/gauss.h @@ -23,453 +23,45 @@ #include #include #include +#include -#define NOSSE 1 +// classical filtering if the support window is small: -#ifndef NOSSE -#include -template void gaussHorizontal3 (T** src, T** dst, T* buffer, int W, int row_from, int row_to, const float c0, const float c1) { - - typedef float pfloat[4]; - pfloat* temp = (pfloat*)buffer + 1; - pfloat* tmp = (pfloat*)buffer; +template void gaussHorizontal3 (T** src, T** dst, T* buffer, int W, int H, const float c0, const float c1, bool multiThread) { - __m128 xmm1; xmm1 = _mm_load1_ps (&c0); - __m128 xmm2; xmm2 = _mm_load1_ps (&c1); - __m128 xmm3; __m128 xmm4; __m128 xmm5; __m128 xmm6; - __m128 xmm0; - - int i; - for (i = row_from; i void gaussHorizontal3 (T** src, T** dst, T* buffer, int W, int row_from, int row_to, const int c0, const int c1) { +template void gaussVertical3 (T** src, T** dst, T* buffer, int W, int H, const float c0, const float c1, bool multiThread) { - time_t t1 = clock (); - - const int csum = c0 + 2 * c1; - for (int i = row_from; i void gaussVertical3 (T** src, T** dst, T* buffer, int H, int col_from, int col_to, const float c0, const float c1) { - - typedef float pfloat[4]; - pfloat* temp = (pfloat*)buffer + 1; - pfloat* tmp = (pfloat*)buffer; - - __m128 xmm1; xmm1 = _mm_load1_ps (&c0); - __m128 xmm2; xmm2 = _mm_load1_ps (&c1); - __m128 xmm3; __m128 xmm4; __m128 xmm5; __m128 xmm6; - __m128 xmm0; - - int i; - for (i = col_from; i void gaussVertical3 (T** src, T** dst, T* buffer, int H, int col_from, int col_to, const float c0, const float c1) { - - time_t t1 = clock (); - - int stride = dst[1] - dst[0]; - for (int j=col_from; j void gaussHorizontal (T** src, T** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma) { - - if (sigma<0.25) { - // dont perform filtering - if (src!=dst) - for (int i = row_from; i (src, dst, (T*)(buffer->data), W, row_from, row_to, (int) round(c0), 2); - double c1 = exp (-1.0 / (2.0 * sigma * sigma)); - double csum = 2.0 * c1 + 1.0; - c1 /= csum; - double c0 = 1.0 / csum; - gaussHorizontal3 (src, dst, (T*)(buffer->data), W, row_from, row_to, c0, c1); - return; - } - - // horizontal - float q = 0.98711 * sigma - 0.96330; - if (sigma<2.5) - q = 3.97156 - 4.14554 * sqrt (1.0 - 0.26891 * sigma); - float b0 = 1.57825 + 2.44413*q + 1.4281*q*q + 0.422205*q*q*q; - float b1 = 2.44413*q + 2.85619*q*q + 1.26661*q*q*q; - float b2 = -1.4281*q*q - 1.26661*q*q*q; - float b3 = 0.422205*q*q*q; - float B = 1.0 - (b1+b2+b3) / b0; - - b1 /= b0; - b2 /= b0; - b3 /= b0; - - // SSE optimized version: - __m128 xmm1; xmm1 = _mm_load1_ps (&B); - __m128 xmm2; xmm2 = _mm_load1_ps (&b1); - __m128 xmm3; xmm3 = _mm_load1_ps (&b2); - __m128 xmm4; xmm4 = _mm_load1_ps (&b3); - __m128 xmm5; - __m128 xmm6; - - typedef float pfloat[4]; - pfloat* temp = (pfloat*)buffer->data + 1; - pfloat* tmp = (pfloat*)buffer->data; - - memset (temp, 0, W*sizeof(pfloat)); - - int i; - for (i=row_from; i=0; j--) { - xmm5 = _mm_load_ps ((float*)&temp[j]); - xmm5 =_mm_mul_ps (xmm1, xmm5); - xmm6 = _mm_load_ps ((float*)&temp[j+1]); - xmm6 = _mm_mul_ps (xmm2, xmm6); - xmm5 = _mm_add_ps (xmm6, xmm5); - xmm6 = _mm_load_ps ((float*)&temp[j+2]); - xmm6 = _mm_mul_ps (xmm6, xmm3); - xmm5 = _mm_add_ps (xmm5, xmm6); - xmm6 = _mm_load_ps ((float*)&temp[j+3]); - xmm6 = _mm_mul_ps (xmm6, xmm4); - xmm5 = _mm_add_ps (xmm5, xmm6); - _mm_store_ps ((float*)&temp[j], xmm5); - } - for (int j=0; jdata; - for (; i=0; j--) - temp2[j] = B * temp2[j] + b1*temp2[j+1] + b2*temp2[j+2] + b3*temp2[j+3]; - for (int j=0; j void gaussVertical (T** src, T** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma) { - - if (sigma<0.25) { - // dont perform filtering - if (src!=dst) - for (int i = 0; i (src, dst, (T*)(buffer->data), H, col_from, col_to, c0, c1); -// double c0 = 2.0 / exp (-1.0 / (2.0 * sigma * sigma)); -// gaussVertical3 (src, dst, (T*)(buffer->data), H, col_from, col_to, (int) round(c0), 2); - return; - } - - // vertical - float q = 0.98711 * sigma - 0.96330; - if (sigma<2.5) - q = 3.97156 - 4.14554 * sqrt (1.0 - 0.26891 * sigma); - float b0 = 1.57825 + 2.44413*q + 1.4281*q*q + 0.422205*q*q*q; - float b1 = 2.44413*q + 2.85619*q*q + 1.26661*q*q*q; - float b2 = -1.4281*q*q - 1.26661*q*q*q; - float b3 = 0.422205*q*q*q; - float B = 1.0 - (b1+b2+b3) / b0; - - b1 /= b0; - b2 /= b0; - b3 /= b0; - - // SSE optimized version: - __m128 xmm1; xmm1 = _mm_load1_ps (&B); - __m128 xmm2; xmm2 = _mm_load1_ps (&b1); - __m128 xmm3; xmm3 = _mm_load1_ps (&b2); - __m128 xmm4; xmm4 = _mm_load1_ps (&b3); - __m128 xmm5; - __m128 xmm6; - - typedef float pfloat[4]; - pfloat* temp = (pfloat*)buffer->data + 1; - pfloat* tmp = (pfloat*)buffer->data; - - memset (temp, 0, H*sizeof(pfloat)); - - int i; - for (i=col_from; i=0; j--) { - xmm5 = _mm_load_ps ((float*)&temp[j]); - xmm5 =_mm_mul_ps (xmm1, xmm5); - xmm6 = _mm_load_ps ((float*)&temp[j+1]); - xmm6 = _mm_mul_ps (xmm2, xmm6); - xmm5 = _mm_add_ps (xmm6, xmm5); - xmm6 = _mm_load_ps ((float*)&temp[j+2]); - xmm6 = _mm_mul_ps (xmm6, xmm3); - xmm5 = _mm_add_ps (xmm5, xmm6); - xmm6 = _mm_load_ps ((float*)&temp[j+3]); - xmm6 = _mm_mul_ps (xmm6, xmm4); - xmm5 = _mm_add_ps (xmm5, xmm6); - _mm_store_ps ((float*)&temp[j], xmm5); - } - for (int j=0; jdata; - for (; i=0; j--) - temp2[j] = B * temp2[j] + b1*temp2[j+1] + b2*temp2[j+2] + b3*temp2[j+3]; - for (int j=0; j void gaussHorizontal3 (T** src, T** dst, T* buffer, int W, int row_from, int row_to, const float c0, const float c1) { - - for (int i=row_from; i void gaussVertical3 (T** src, T** dst, T* buffer, int H, int col_from, int col_to, const float c0, const float c1) { - - for (int i=col_from; i void gaussHorizontal (T** src, T** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma) { +// fast gaussian approximation if the support window is large + +template void gaussHorizontal (T** src, T** dst, AlignedBuffer* buffer, int W, int H, double sigma, bool multiThread) { if (sigma<0.25) { // dont perform filtering if (src!=dst) - for (int i = row_from; i void gaussHorizontal (T** src, T** dst, AlignedBuffer* double csum = 2.0 * c1 + 1.0; c1 /= csum; double c0 = 1.0 / csum; - gaussHorizontal3 (src, dst, (T*)(buffer->data), W, row_from, row_to, c0, c1); + gaussHorizontal3 (src, dst, (T*)(buffer->data), W, H, c0, c1, multiThread); return; } - // horizontal + // coefficient calculation double q = 0.98711 * sigma - 0.96330; if (sigma<2.5) q = 3.97156 - 4.14554 * sqrt (1.0 - 0.26891 * sigma); @@ -513,8 +105,11 @@ template void gaussHorizontal (T** src, T** dst, AlignedBuffer* for (int j=0; j<3; j++) M[i][j] /= (1.0+b1-b2+b3)*(1.0+b2+(b1-b3)*b3); - double* temp2 = (double*)buffer->data; - for (int i=row_from; idata + omp_get_thread_num() * W; + temp2[0] = B * src[i][0] + b1*src[i][0] + b2*src[i][0] + b3*src[i][0]; temp2[1] = B * src[i][1] + b1*temp2[0] + b2*src[i][0] + b3*src[i][0]; temp2[2] = B * src[i][2] + b1*temp2[1] + b2*temp2[0] + b3*src[i][0]; @@ -537,14 +132,13 @@ template void gaussHorizontal (T** src, T** dst, AlignedBuffer* } } -template void gaussVertical (T** src, T** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma) { +template void gaussVertical (T** src, T** dst, AlignedBuffer* buffer, int W, int H, double sigma, bool multiThread) { if (sigma<0.25) { // dont perform filtering if (src!=dst) - for (int i = 0; i void gaussVertical (T** src, T** dst, AlignedBuffer* b double csum = 2.0 * c1 + 1.0; c1 /= csum; double c0 = 1.0 / csum; - gaussVertical3 (src, dst, (T*)(buffer->data), H, col_from, col_to, c0, c1); + gaussVertical3 (src, dst, (T*)(buffer->data), W, H, c0, c1, multiThread); return; } - // vertical + // coefficient calculation double q = 0.98711 * sigma - 0.96330; if (sigma<2.5) q = 3.97156 - 4.14554 * sqrt (1.0 - 0.26891 * sigma); @@ -587,9 +181,12 @@ template void gaussVertical (T** src, T** dst, AlignedBuffer* b for (int j=0; j<3; j++) M[i][j] /= (1.0+b1-b2+b3)*(1.0+b2+(b1-b3)*b3); - double* temp2 = (double*)buffer->data; - for (int i=col_from; idata + omp_get_thread_num() * H; + + temp2[0] = B * src[0][i] + b1*src[0][i] + b2*src[0][i] + b3*src[0][i]; temp2[1] = B * src[1][i] + b1*temp2[0] + b2*src[0][i] + b3*src[0][i]; temp2[2] = B * src[2][i] + b1*temp2[1] + b2*temp2[0] + b3*src[0][i]; @@ -611,14 +208,13 @@ template void gaussVertical (T** src, T** dst, AlignedBuffer* b dst[j][i] = (T)temp2[j]; } } -#endif - +/* void gaussHorizontal_unsigned (unsigned short** src, unsigned short** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma); void gaussVertical_unsigned (unsigned short** src, unsigned short** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma); void gaussHorizontal_signed (short** src, short** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma); void gaussVertical_signed (short** src, short** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma); void gaussHorizontal_float (float** src, float** dst, AlignedBuffer* buffer, int W, int row_from, int row_to, double sigma); void gaussVertical_float (float** src, float** dst, AlignedBuffer* buffer, int H, int col_from, int col_to, double sigma); - +*/ #endif diff --git a/rtengine/improccoordinator.cc b/rtengine/improccoordinator.cc index 0e8cf81f1..a015da737 100644 --- a/rtengine/improccoordinator.cc +++ b/rtengine/improccoordinator.cc @@ -28,7 +28,7 @@ namespace rtengine { extern Settings* settings; ImProcCoordinator::ImProcCoordinator () - : allocated(false), scale(-1), pW(-1), pH(-1), + : ipf(¶ms, true), allocated(false), scale(-1), pW(-1), pH(-1), plistener(NULL), imageListener(NULL), aeListener(NULL), hListener(NULL), resultValid(false), awbComputed(false), changeSinceLast(0), updaterRunning(false), @@ -54,7 +54,6 @@ ImProcCoordinator::~ImProcCoordinator () { delete toDel[i]; imgsrc->decreaseRef (); - ipf.release (); updaterThreadStart.unlock (); } @@ -79,6 +78,8 @@ void ImProcCoordinator::updatePreviewImage (int todo) { else if (params.resize.dataspec==2) params.resize.scale = (double)params.resize.height / (params.coarse.rotate==90 || params.coarse.rotate==270 ? fw : fh); + ipf.setScale (scale); + progress ("Applying white balance, color correction & sRBG conversion...",100*readyphase/numofphases); if (todo & M_INIT) { minit.lock (); @@ -164,7 +165,7 @@ void ImProcCoordinator::updatePreviewImage (int todo) { progress ("Exposure curve & CIELAB conversion...",100*readyphase/numofphases); if (todo & M_RGBCURVE) { CurveFactory::complexCurve (params.toneCurve.expcomp, params.toneCurve.black/65535.0, params.toneCurve.hlcompr, params.toneCurve.shcompr, params.toneCurve.brightness, params.toneCurve.contrast, imgsrc->getDefGain(), imgsrc->getGamma(), true, params.toneCurve.curve, vhist16, tonecurve, bcrgbhist, scale==1 ? 1 : 1); - ipf.rgbProc (oprevi, oprevl, ¶ms, tonecurve, shmap); + ipf.rgbProc (oprevi, oprevl, tonecurve, shmap); // recompute luminance histogram memset (lhist16, 0, 65536*sizeof(int)); @@ -186,12 +187,12 @@ void ImProcCoordinator::updatePreviewImage (int todo) { readyphase++; if (scale==1) { progress ("Denoising luminance...",100*readyphase/numofphases); - ipf.lumadenoise (nprevl, ¶ms, scale*params.resize.scale, buffer); + ipf.lumadenoise (nprevl, buffer); } readyphase++; if (scale==1) { progress ("Sharpening...",100*readyphase/numofphases); - ipf.sharpening (nprevl, ¶ms, scale*params.resize.scale, (unsigned short**)buffer); + ipf.sharpening (nprevl, (unsigned short**)buffer); } readyphase++; } @@ -201,11 +202,11 @@ void ImProcCoordinator::updatePreviewImage (int todo) { if (todo & M_COLOR) { progress ("Applying Color Boost...",100*readyphase/numofphases); - ipf.colorCurve (oprevl, nprevl, ¶ms); + ipf.colorCurve (oprevl, nprevl); readyphase++; if (scale==1) { progress ("Denoising color...",100*readyphase/numofphases); - ipf.colordenoise (nprevl, ¶ms, scale*params.resize.scale, buffer); + ipf.colordenoise (nprevl, buffer); } readyphase++; } @@ -314,7 +315,7 @@ if (settings->verbose) printf ("setscale before lock\n"); oprevl = new LabImage (pW, pH); nprevl = new LabImage (pW, pH); previmg = new Image8 (pW, pH); - shmap = new SHMap (pW, pH); + shmap = new SHMap (pW, pH, true); buffer = new int*[pH]; for (int i=0; i0) diff --git a/rtengine/improcfun.cc b/rtengine/improcfun.cc index beae7850d..afa99d020 100644 --- a/rtengine/improcfun.cc +++ b/rtengine/improcfun.cc @@ -27,97 +27,36 @@ #include #include #include +#include namespace rtengine { using namespace procparams; +#undef MAXVAL +#undef CMAXVAL +#undef MAXL #undef MAX #undef MIN -#undef MAXVAL +#undef ABS #undef CLIP #undef CLIPS #undef CLIPC #undef CLIPTO -#undef CLIPTOC -#undef THREAD_PRIORITY_NORMAL #define MAXVAL 0xffff -#define CLIP(a) ((a)>0?((a)(b)?(b):(a)) +#define ABS(a) ((a)<0?-(a):(a)) +#define CLIP(a) ((a)>0?((a)-32768?((a)<32767?(a):32767):-32768) #define CLIPC(a) ((a)>-32000?((a)<32000?(a):32000):-32000) -#define MAX(a,b) ((a)<(b)?(b):(a)) -#define MIN(a,b) ((a)>(b)?(b):(a)) #define CLIPTO(a,b,c) ((a)>(b)?((a)<(c)?(a):(c)):(b)) -#define CLIPTOC(a,b,c,d) ((a)>=(b)?((a)<=(c)?(a):((c),d=true)):((b),d=true)) extern const Settings* settings; -// -// STRUCTURES FOR THE SHADOW MAP AND LAB SPACE IMAGE -//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -LabImage::LabImage (int w, int h) : W(w), H(h), fromImage(false) { - - L = new unsigned short*[H]; - for (int i=0; iwidth; - H = im->height; - L = im->r; - a = (short**) im->g; - b = (short**) im->b; - fromImage = true; -} - -LabImage::~LabImage () { - - if (!fromImage) { - for (int i=0; i0?((a)(c)?(c):(a))) - -#define MAX(a,b) ((a)<(b)?(b):(a)) -#define MIN(a,b) ((a)>(b)?(b):(a)) - -#define MAXL 65535 -#define ABS(a) ((a)<0?-(a):(a)) - - int* ImProcFunctions::cacheL; int* ImProcFunctions::cachea; int* ImProcFunctions::cacheb; @@ -126,17 +65,6 @@ int* ImProcFunctions::ycache; int* ImProcFunctions::zcache; unsigned short ImProcFunctions::gamma2curve[65536]; -/*const int c00 = (int) (32768.0 * 0.412453 / 0.950456); -const int c01 = (int) (32768.0 * 0.357580 / 0.950456); -const int c02 = (int) (32768.0 * 0.180423 / 0.950456); -const int c10 = (int) (32768.0 * 0.212671); -const int c11 = (int) (32768.0 * 0.715160); -const int c12 = (int) (32768.0 * 0.072169); -const int c20 = (int) (32768.0 * 0.019334 / 1.088754); -const int c21 = (int) (32768.0 * 0.119193 / 1.088754); -const int c22 = (int) (32768.0 * 0.950227 / 1.088754); -*/ - void ImProcFunctions::initCache () { int maxindex = 2*65536; @@ -180,11 +108,14 @@ void ImProcFunctions::initCache () { } } -void ImProcFunctions::release () { +ImProcFunctions::~ImProcFunctions () { if (monitorTransform!=NULL) cmsDeleteTransform (monitorTransform); - monitorTransform = NULL; +} + +void ImProcFunctions::setScale (double iscale) { + scale = iscale; } void ImProcFunctions::firstAnalysis_ (Image16* original, Glib::ustring wprofile, unsigned int* histogram, int* chroma_radius, int row_from, int row_to) { @@ -244,13 +175,8 @@ void ImProcFunctions::firstAnalysis_ (Image16* original, Glib::ustring wprofile, void ImProcFunctions::firstAnalysis (Image16* original, const ProcParams* params, unsigned int* histogram, double gamma) { - int cr1, cr2; - unsigned int* hist1 = new unsigned int[65536]; memset (hist1, 0, 65536*sizeof(int)); - unsigned int* hist2 = new unsigned int[65536]; memset (hist2, 0, 65536*sizeof(int)); - - int H = original->height; - - Glib::ustring wprofile = params->icm.working; + // set up monitor transform + Glib::ustring wprofile = params->icm.working; if (monitorTransform) cmsDeleteTransform (monitorTransform); monitorTransform = NULL; @@ -264,47 +190,49 @@ void ImProcFunctions::firstAnalysis (Image16* original, const ProcParams* params monitorTransform = cmsCreateTransform (iprof, TYPE_RGB_16, monitor, TYPE_RGB_8, settings->colorimetricIntent, 0); lcmsMutex->unlock (); } - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::firstAnalysis_), original, wprofile, hist1, &cr1, 0, H/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::firstAnalysis_), original, wprofile, hist2, &cr2, H/2, H), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); + + // calculate chroma radius and histogram of the y channel needed for exposure curve calculation + int T = omp_get_max_threads(); + int* cr = new int [T]; + unsigned int** hist = new unsigned int* [T]; + for (int i=0; iheight; + #pragma omp parallel if (multiThread) + { + int tid = omp_get_thread_num(); + int nthreads = omp_get_num_threads(); + int blk = H/nthreads; + + if (tidchroma_radius) + chroma_radius = cr[i]; + memset (histogram, 0, 65536*sizeof(int)); for (int i=0; i<65536; i++) - histogram[i] = hist1[i]+hist2[i]; + for (int j=0; jdualThreadEnabled) { - - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::rgbProc_), working, lab, params, tonecurve, shmap, 0, working->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::rgbProc_), working, lab, params, tonecurve, shmap, working->height/2, working->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - rgbProc_ (working, lab, params, tonecurve, shmap, 0, working->height); -} - -void ImProcFunctions::rgbProc_ (Image16* working, LabImage* lab, const ProcParams* params, int* tonecurve, SHMap* shmap, int row_from, int row_to) { - - int r, g, b; +void ImProcFunctions::rgbProc (Image16* working, LabImage* lab, int* tonecurve, SHMap* shmap) { int h_th, s_th; if (shmap) { @@ -333,7 +261,11 @@ void ImProcFunctions::rgbProc_ (Image16* working, LabImage* lab, const ProcParam int mapval; double factor; int tW = working->width; - for (int i=row_from; iheight; + + #pragma omp parallel for if (multiThread) + for (int i=0; ir[i][j]; @@ -375,9 +307,6 @@ void ImProcFunctions::rgbProc_ (Image16* working, LabImage* lab, const ProcParam g = tonecurve[g]; b = tonecurve[b]; -// int x = (14219 * r + 12328 * g + 6220 * b) >> 15; -// int y = ( 6968 * r + 23434 * g + 2365 * b) >> 15; -// int z = ( 582 * r + 3587 * g + 28598 * b) >> 15; int x = (toxyz[0][0] * r + toxyz[1][0] * g + toxyz[2][0] * b) >> 15; int y = (toxyz[0][1] * r + toxyz[1][1] * g + toxyz[2][1] * b) >> 15; int z = (toxyz[0][2] * r + toxyz[1][2] * g + toxyz[2][2] * b) >> 15; @@ -405,7 +334,7 @@ void ImProcFunctions::luminanceCurve (LabImage* lold, LabImage* lnew, int* curve #include "cubic.cc" -void ImProcFunctions::colorCurve (LabImage* lold, LabImage* lnew, const ProcParams* params) { +void ImProcFunctions::colorCurve (LabImage* lold, LabImage* lnew) { double* cmultiplier = new double [181021]; @@ -443,43 +372,13 @@ void ImProcFunctions::colorCurve (LabImage* lold, LabImage* lnew, const ProcPara } } - if (settings->dualThreadEnabled) { - - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::colorCurve_), lold, lnew, params, 0, lnew->H/2, cmultiplier), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::colorCurve_), lold, lnew, params, lnew->H/2, lnew->H, cmultiplier), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - colorCurve_ (lold, lnew, params, 0, lnew->H, cmultiplier); - - delete [] cmultiplier; -} - -void ImProcFunctions::colorCurve_ (LabImage* lold, LabImage* lnew, const ProcParams* params, int row_from, int row_to, double* cmultiplier) { - - double boost_a = (params->colorBoost.amount + 100.0) / 100.0; - double boost_b = (params->colorBoost.amount + 100.0) / 100.0; - - double c, amul = 1.0, bmul = 1.0; - if (boost_a > boost_b) { - c = boost_a; - if (boost_a > 0) - bmul = boost_b / boost_a; - } - else { - c = boost_b; - if (boost_b > 0) - amul = boost_a / boost_b; - } - - int nna, nnb; double shift_a = params->colorShift.a * chroma_scale, shift_b = params->colorShift.b * chroma_scale; short** oa = lold->a; short** ob = lold->b; - for (int i=row_from; iH; i++) for (int j=0; jW; j++) { double wanted_c = c; @@ -504,1696 +403,36 @@ void ImProcFunctions::colorCurve_ (LabImage* lold, LabImage* lnew, const ProcPar } } - nna = (int)((oa[i][j]+shift_a) * real_c * amul); - nnb = (int)((ob[i][j]+shift_b) * real_c * bmul); - lnew->a[i][j] = CLIPV(nna,-32000,32000); - lnew->b[i][j] = CLIPV(nnb,-32000,32000); + int nna = (int)((oa[i][j]+shift_a) * real_c * amul); + int nnb = (int)((ob[i][j]+shift_b) * real_c * bmul); + lnew->a[i][j] = CLIPTO(nna,-32000,32000); + lnew->b[i][j] = CLIPTO(nnb,-32000,32000); } + + delete [] cmultiplier; } -void blur (float** src, float** dst, int W, int H, int r) { +void ImProcFunctions::lumadenoise (LabImage* lab, int** b2) { - float** tmpI = new float*[H]; - for (int i=0; i=H-r) - tmpI[i][j] = src[i][j]; - else { - int num = 0; - float sum = 0.0; - for (int x=-r; x<=r; x++) - for (int y=-r; y<=r; y++) - if (x*x+y*y<=r*r) { - sum += src[i+x][j+y]; - num++; - } - tmpI[i][j] = sum / num; - } - } - for (int i=0; ilumaDenoise.enabled && lab->W>=8 && lab->H>=8) + bilateral (lab->L, lab->L, (unsigned short**)b2, lab->W, lab->H, params->lumaDenoise.radius / scale, params->lumaDenoise.edgetolerance, multiThread); } -void ImProcFunctions::damping_ (float** aI, unsigned short** aO, float damping, int W, int rowfrom, int rowto) { - - for (int i=rowfrom; isharpening.enabled==false || params->sharpening.deconvamount<1) - return; - - int W = lab->W, H = lab->H; - - float** tmpI = new float*[H]; - for (int i=0; iL[i][j]; - } - - float** tmp = (float**)b2; - - AlignedBuffer* buffer1 = new AlignedBuffer (MAX(W,H)*5); - AlignedBuffer* buffer2 = new AlignedBuffer (MAX(W,H)*5); - - float damping = params->sharpening.deconvdamping / 5.0; - bool needdamp = params->sharpening.deconvdamping > 0; - for (int k=0; ksharpening.deconviter; k++) { - -// gaussHorizontal (tmpI, tmp, buffer1, W, 0, H, params->sharpening.deconvradius / scale); -// gaussVertical (tmp, tmp, buffer1, H, 0, W, params->sharpening.deconvradius / scale); - // apply blur function (gaussian blur) - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_float), tmpI, tmp, buffer1, W, 0, H/2, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_float), tmpI, tmp, buffer2, W, H/2, H, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - gaussHorizontal_float (tmpI, tmp, buffer1, W, 0, H, params->sharpening.deconvradius / scale); - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_float), tmp, tmp, buffer1, H, 0, W/2, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_float), tmp, tmp, buffer2, H, W/2, W, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - gaussVertical_float (tmp, tmp, buffer1, H, 0, W, params->sharpening.deconvradius / scale); - -// blur (tmpI, tmp, W, H, params->sharpening.radius / scale); - if (!needdamp) { - for (int i=0; i0) - tmp[i][j] = (float)lab->L[i][j] / tmp[i][j]; - } - else { - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::damping_), tmp, lab->L, damping, W, 0, H/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::damping_), tmp, lab->L, damping, W, H/2, H), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - damping_ (tmp, lab->L, damping, W, 0, H); - } -// gaussHorizontal (tmp, tmp, buffer1, W, 0, H, params->sharpening.deconvradius / scale); -// gaussVertical (tmp, tmp, buffer1, H, 0, W, params->sharpening.deconvradius / scale); - // apply blur function (gaussian blur) - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_float), tmp, tmp, buffer1, W, 0, H/2, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_float), tmp, tmp, buffer2, W, H/2, H, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - gaussHorizontal_float (tmp, tmp, buffer1, W, 0, H, params->sharpening.deconvradius / scale); - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_float), tmp, tmp, buffer1, H, 0, W/2, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_float), tmp, tmp, buffer2, H, W/2, W, params->sharpening.deconvradius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - gaussVertical_float (tmp, tmp, buffer1, H, 0, W, params->sharpening.deconvradius / scale); - -// blur (tmp, tmp, W, H, params->sharpening.radius / scale); - - for (int i=0; iL[i][j] = lab->L[i][j]*(100-params->sharpening.deconvamount) / 100 + (int)CLIP(tmpI[i][j])*params->sharpening.deconvamount / 100; - - for (int i=0; isharpening.method=="rld") { - deconvsharpening (lab, params, scale, b2); - return; - } - - if (params->sharpening.enabled==false || params->sharpening.amount<1 || lab->W<8 || lab->H<8) - return; - - int W = lab->W, H = lab->H; - unsigned short** b3; - if (params->sharpening.edgesonly==false) { - - AlignedBuffer* buffer1 = new AlignedBuffer (MAX(W,H)*5); - AlignedBuffer* buffer2 = new AlignedBuffer (MAX(W,H)*5); - - MyTime t1, t2, t3; - t1.set (); - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_unsigned), lab->L, b2, buffer1, W, 0, H/2, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_unsigned), lab->L, b2, buffer2, W, H/2, H, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - gaussHorizontal_unsigned (lab->L, b2, buffer1, W, 0, H, params->sharpening.radius / scale); - - t2.set (); -// printf ("Horizontal: %d\n", t2.etime (t1)); - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_unsigned), b2, b2, buffer1, H, 0, W/2, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_unsigned), b2, b2, buffer2, H, W/2, W, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - gaussVertical_unsigned (b2, b2, buffer1, H, 0, W, params->sharpening.radius / scale); - - t3.set (); -// printf ("Vertical: %d\n", t3.etime (t2)); - - delete buffer1; - delete buffer2; - } - else { - b3 = new unsigned short*[H]; - for (int i=0; i* buffer1 = new AlignedBuffer (MAX(W,H)*5); - AlignedBuffer* buffer2 = new AlignedBuffer (MAX(W,H)*5); - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_unsigned), lab->L, b3, b2, dim1, params->sharpening.edges_radius / scale, params->sharpening.edges_tolerance), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_unsigned), lab->L, b3, b2, dim2, params->sharpening.edges_radius / scale, params->sharpening.edges_tolerance), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_unsigned), b3, b2, buffer1, W, 0, H/2, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_unsigned), b3, b2, buffer2, W, H/2, H, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_unsigned), b2, b2, buffer1, H, 0, W/2, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_unsigned), b2, b2, buffer2, H, W/2, W, params->sharpening.radius / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else { - bilateral_unsigned (lab->L, (unsigned short**)b2, b3, dim1, params->sharpening.edges_radius / scale, params->sharpening.edges_tolerance); - bilateral_unsigned (lab->L, (unsigned short**)b2, b3, dim2, params->sharpening.edges_radius / scale, params->sharpening.edges_tolerance); - gaussHorizontal_unsigned (b2, b2, buffer1, W, 0, H, params->sharpening.radius / scale); - gaussVertical_unsigned (b2, b2, buffer1, H, 0, W, params->sharpening.radius / scale); - } - - delete buffer1; - delete buffer2; - } - unsigned short** base = lab->L; - if (params->sharpening.edgesonly) - base = b3; - - if (params->sharpening.halocontrol==false) { - for (int i=0; iparams->sharpening.threshold) { - int val = lab->L[i][j] + params->sharpening.amount * diff / 100; - lab->L[i][j] = CLIP(val); - } - } - } - else { - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::sharpenHaloCtrl), lab, params, b2, base, W, 2, H/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::sharpenHaloCtrl), lab, params, b2, base, W, H/2, H-2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - sharpenHaloCtrl (lab, params, b2, base, W, 2, H-2); - } - - if (params->sharpening.edgesonly) { - for (int i=0; isharpening.halocontrol_amount); - unsigned short** nL = base; - for (int i=row_from; i params->sharpening.threshold) { - // compute maximum/minimum in a delta environment - np1 = 2*(nL[i-2][j] + nL[i-2][j+1] + nL[i-2][j+2] + nL[i-1][j] + nL[i-1][j+1] + nL[i-1][j+2] + nL[i][j] + nL[i][j+1] + nL[i][j+2]) / 27 + nL[i-1][j+1] / 3; - np2 = 2*(nL[i-1][j] + nL[i-1][j+1] + nL[i-1][j+2] + nL[i][j] + nL[i][j+1] + nL[i][j+2] + nL[i+1][j] + nL[i+1][j+1] + nL[i+1][j+2]) / 27 + nL[i][j+1] / 3; - np3 = 2*(nL[i][j] + nL[i][j+1] + nL[i][j+2] + nL[i+1][j] + nL[i+1][j+1] + nL[i+1][j+2] + nL[i+2][j] + nL[i+2][j+1] + nL[i+2][j+2]) / 27 + nL[i+1][j+1] / 3; - MINMAX3(np1,np2,np3,maxn,minn); - MAX3(max1,max2,maxn,max); - MIN3(min1,min2,minn,min); - max1 = max2; max2 = maxn; - min1 = min2; min2 = minn; - if (max < lab->L[i][j]) - max = lab->L[i][j]; - if (min > lab->L[i][j]) - min = lab->L[i][j]; - int val = lab->L[i][j] + params->sharpening.amount * diff / 100; - int newL = CLIP(val); - // applying halo control - if (newL > max) - newL = max + (newL-max) * scale / 10000; - else if (newLL[i][j] = newL; - } - } - } -} - -void ImProcFunctions::lumadenoise (LabImage* lab, const ProcParams* params, double scale, int** b2) { - -// MyTime t1, t2; -// t1.set (); - - if (params->lumaDenoise.enabled && lab->W>=8 && lab->H>=8) { - - Dim dim1 (lab->W, lab->H, 0, lab->H/2); - Dim dim2 (lab->W, lab->H, lab->H/2, lab->H); - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_unsigned), lab->L, lab->L, (unsigned short**)b2, dim1, params->lumaDenoise.radius / scale, params->lumaDenoise.edgetolerance), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_unsigned), lab->L, lab->L, (unsigned short**)b2, dim2, params->lumaDenoise.radius / scale, params->lumaDenoise.edgetolerance), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else { - bilateral_unsigned (lab->L, lab->L, (unsigned short**)b2, dim1, params->lumaDenoise.radius / scale, params->lumaDenoise.edgetolerance); - bilateral_unsigned (lab->L, lab->L, (unsigned short**)b2, dim2, params->lumaDenoise.radius / scale, params->lumaDenoise.edgetolerance); - } - } -// t2.set (); -// printf ("Luminance denoising time = %d\n", t2.etime (t1)); -} - -void ImProcFunctions::colordenoise (LabImage* lab, const ProcParams* params, double scale, int** b2) { +void ImProcFunctions::colordenoise (LabImage* lab, int** b2) { if (params->colorDenoise.enabled && lab->W>=8 && lab->H>=8) { -/* if (params->colorDenoise.edgesensitive) { + AlignedBuffer* buffer = new AlignedBuffer (MAX(lab->W,lab->H)*omp_get_max_threads()); - short** buffer1 = (short**)b2; - short** buffer2 = new short*[lab->H]; - for (int i=0; iH; i++) - buffer2[i] = buffer1[i]+lab->W; - Dim dim (lab->W, lab->H, 0, lab->H); - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_signed), lab->a, lab->a, buffer1, dim, params->colorDenoise.radius / scale, params->colorDenoise.edgetolerance), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_signed), lab->b, lab->b, buffer2, dim, params->colorDenoise.radius / scale, params->colorDenoise.edgetolerance), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else { - bilateral_signed (lab->a, lab->a, buffer1, dim, params->colorDenoise.radius / scale, params->colorDenoise.edgetolerance); - bilateral_signed (lab->b, lab->b, buffer1, dim, params->colorDenoise.radius / scale, params->colorDenoise.edgetolerance); - } - delete [] buffer2; - } - else { -*/ - AlignedBuffer* buffer1 = new AlignedBuffer (MAX(lab->W,lab->H)*5); - AlignedBuffer* buffer2 = new AlignedBuffer (MAX(lab->W,lab->H)*5); - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_signed), lab->a, lab->a, buffer1, lab->W, 0, lab->H, params->colorDenoise.amount / 10.0 / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_signed), lab->b, lab->b, buffer2, lab->W, 0, lab->H, params->colorDenoise.amount / 10.0 / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_signed), lab->a, lab->a, buffer1, lab->H, 0, lab->W, params->colorDenoise.amount / 10.0 / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_signed), lab->b, lab->b, buffer2, lab->H, 0, lab->W, params->colorDenoise.amount / 10.0 / scale), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else { - gaussHorizontal_signed (lab->a, lab->a, buffer1, lab->W, 0, lab->H, params->colorDenoise.amount / 10.0 / scale); - gaussHorizontal_signed (lab->b, lab->b, buffer1, lab->W, 0, lab->H, params->colorDenoise.amount / 10.0 / scale); - gaussVertical_signed (lab->a, lab->a, buffer1, lab->H, 0, lab->W, params->colorDenoise.amount / 10.0 / scale); - gaussVertical_signed (lab->b, lab->b, buffer1, lab->H, 0, lab->W, params->colorDenoise.amount / 10.0 / scale); - } - delete buffer1; - delete buffer2; -// } + gaussHorizontal (lab->a, lab->a, buffer, lab->W, lab->H, params->colorDenoise.amount / 10.0 / scale, multiThread); + gaussHorizontal (lab->b, lab->b, buffer, lab->W, lab->H, params->colorDenoise.amount / 10.0 / scale, multiThread); + gaussVertical (lab->a, lab->a, buffer, lab->W, lab->H, params->colorDenoise.amount / 10.0 / scale, multiThread); + gaussVertical (lab->b, lab->b, buffer, lab->W, lab->H, params->colorDenoise.amount / 10.0 / scale, multiThread); + + delete buffer; } } -void ImProcFunctions::vignetting_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { - - int oW = sizes.oW; - int oH = sizes.oH; - int cx = sizes.cx; - int cy = sizes.cy; - - double w2 = (double) oW / 2.0 - 0.5; - double h2 = (double) oH / 2.0 - 0.5; - - double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2; - - double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; - double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; - - double mul = (1.0-v) / tanh(b); - - int val; - for (int y=row_from; ywidth; x++) { - double x_d = (double) (x + cx) - w2 ; - double r = sqrt(x_d*x_d + y_d*y_d); - double vign = v + mul * tanh (b*(maxRadius-r) / maxRadius); - val = original->r[y][x] / vign; - transformed->r[y][x] = CLIP(val); - val = original->g[y][x] / vign; - transformed->g[y][x] = CLIP(val); - val = original->b[y][x] / vign; - transformed->b[y][x] = CLIP(val); - } - } -} - -void ImProcFunctions::vignetting (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int oW, int oH) { - - STemp sizes; - sizes.cx = cx; - sizes.cy = cy; - sizes.oW = oW; - sizes.oH = oH; - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::vignetting_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::vignetting_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - vignetting_ (original, transformed, params, sizes, 0, transformed->height); -} - -#include "cubint.cc" -void ImProcFunctions::transform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { - - int oW = sizes.oW; - int oH = sizes.oH; - int cx = sizes.cx; - int cy = sizes.cy; - int sx = sizes.sx; - int sy = sizes.sy; - - double w2 = (double) oW / 2.0 - 0.5; - double h2 = (double) oH / 2.0 - 0.5; - - double cost = cos(params->rotate.degree * 3.14/180.0); - double sint = sin(params->rotate.degree * 3.14/180.0); - - double max_x = (double) (sx + original->width - 1); - double max_y = (double) (sy + original->height - 1); - double min_x = (double) sx; - double min_y = (double) sy; - - const int n2 = 2; - const int n = 4; - - int mix = original->width - 1; // maximum x-index src - int miy = original->height - 1;// maximum y-index src - int mix2 = mix +1 - n; - int miy2 = miy +1 - n; - - double scale = (oW>oH) ? (double)oW / 2.0 : (double)oH / 2.0 ; - double radius = sqrt( (double)( oW*oW + oH*oH ) ); - radius /= (oWdistortion.amount; - - double d = 1.0 - a; - - // magnify image to keep size - double rotmagn = 1.0; - if (params->rotate.fill) { - double beta = atan((double)MIN(oH,oW)/MAX(oW,oH)); - rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); - } - // 1. check upper and lower border - double d1 = rotmagn - a*h2/scale; - double d2 = rotmagn - a*w2/scale; - double d3 = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; - d = MIN(d,MIN(d1,MIN(d2,d3))); - - // auxilary variables for vignetting - double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2 / scale; - - double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; - double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; - - double mul = (1.0-v) / tanh(b); - - // main cycle - double eps = 1e-10; - bool calc_r=( (fabs(a)>eps) || (fabs(1.0-v)>eps) ); - bool do_vign = (fabs(1.0-v)>eps); - - for (int y=row_from; ywidth; x++) { - double x_d = (double) (x + cx) - w2 ; - - double r=0.0; - double s = d;//10000.0; - if (calc_r) - { - r=(sqrt(x_d*x_d + y_d*y_d)) / scale; - if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); - - // Convert only valid pixels - if (valid) { - // Extract integer and fractions of source screen coordinates - int xc = (int) (Dx); Dx -= (double)xc; - int yc = (int) (Dy); Dy -= (double)yc; - int ys = yc +1 - n2 - sy; // smallest y-index used for interpolation - int xs = xc +1 - n2 - sx; // smallest x-index used for interpolation - - double vignmul = 1.0; - if (do_vign) vignmul /= (v + mul * tanh (b*(maxRadius-s*r) / maxRadius)); - - if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2) // all interpolation pixels inside image - cubint (original, xs, ys, Dx, Dy, &(transformed->r[y][x]), &(transformed->g[y][x]), &(transformed->b[y][x]), vignmul); - else { // edge pixels - int y1 = (yc>0) ? yc : 0; - if (y1>miy) y1 = miy; - int y2 = (yc0) ? xc : 0; - if (x1>mix) x1 = mix; - int x2 = (xcr[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->r[y1][x2]*Dx*(1.0-Dy) + original->r[y2][x1]*(1.0-Dx)*Dy + original->r[y2][x2]*Dx*Dy); - int g = vignmul*(original->g[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->g[y1][x2]*Dx*(1.0-Dy) + original->g[y2][x1]*(1.0-Dx)*Dy + original->g[y2][x2]*Dx*Dy); - int b = vignmul*(original->b[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->b[y1][x2]*Dx*(1.0-Dy) + original->b[y2][x1]*(1.0-Dx)*Dy + original->b[y2][x2]*Dx*Dy); - transformed->r[y][x] = CLIP(r); - transformed->g[y][x] = CLIP(g); - transformed->b[y][x] = CLIP(b); - } - } - else { - // not valid (source pixel x,y not inside source image, etc.) - transformed->r[y][x] = 0; - transformed->g[y][x] = 0; - transformed->b[y][x] = 0; - } - } - } -} - -void ImProcFunctions::simpltransform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { - - int oW = sizes.oW; - int oH = sizes.oH; - int cx = sizes.cx; - int cy = sizes.cy; - int sx = sizes.sx; - int sy = sizes.sy; - - double w2 = (double) oW / 2.0 - 0.5; - double h2 = (double) oH / 2.0 - 0.5; - - double cost = cos(params->rotate.degree * 3.14/180.0); - double sint = sin(params->rotate.degree * 3.14/180.0); - - double max_x = (double) (sx + original->width - 1); - double max_y = (double) (sy + original->height - 1); - double min_x = (double) sx; - double min_y = (double) sy; - - const int n2 = 2; - const int n = 2; - - int mix = original->width - 1; // maximum x-index src - int miy = original->height - 1;// maximum y-index src - int mix2 = mix +1 - n; - int miy2 = miy +1 - n; - - double scale = (oW>oH) ? (double)oW / 2.0 : (double)oH / 2.0 ; - double radius = sqrt( (double)( oW*oW + oH*oH ) ); - radius /= (oWdistortion.amount; - - double d = 1.0 - a; - - - // magnify image to keep size - double rotmagn = 1.0; - if (params->rotate.fill) { - double beta = atan((double)MIN(oH,oW)/MAX(oW,oH)); - rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); - } - // 1. check upper and lower border - double d1r = rotmagn - a*h2/scale - params->cacorrection.red; - double d2r = rotmagn - a*w2/scale - params->cacorrection.red; - double d3r = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.red; - double dr = MIN(d,MIN(d1r,MIN(d2r,d3r))); - double d1b = rotmagn - a*h2/scale - params->cacorrection.blue; - double d2b = rotmagn - a*w2/scale - params->cacorrection.blue; - double d3b = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.blue; - double db = MIN(d,MIN(d1b,MIN(d2b,d3b))); - double d1g = rotmagn - a*h2/scale; - double d2g = rotmagn - a*w2/scale; - double d3g = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; - double dg = MIN(d,MIN(d1g,MIN(d2g,d3g))); - - d = MIN(dg,MIN(dr,db)); - - // auxilary variables for vignetting - double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2 / scale; - - double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; - double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; - - double mul = (1.0-v) / tanh(b); - - // main cycle - double eps = 1e-10; - bool calc_r=( (fabs(a)>eps) || (fabs(1.0-v)>eps) ); - bool do_vign = (fabs(1.0-v)>eps); - - for (int y=row_from; ywidth; x++) { - double x_d = (double) (x + cx) - w2 ; - - double r=0.0; - double s = d;//10000.0; - if (calc_r) - { - r=(sqrt(x_d*x_d + y_d*y_d)) / scale; - if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); - - // Convert only valid pixels - if (valid) { - // Extract integer and fractions of source screen coordinates - int xc = (int) (Dx); Dx -= (double)xc; - int yc = (int) (Dy); Dy -= (double)yc; - int ys = yc +1 - n2 - sy; // smallest y-index used for interpolation - int xs = xc +1 - n2 - sx; // smallest x-index used for interpolation - - double vignmul = 1.0; - if (do_vign) vignmul /= (v + mul * tanh (b*(maxRadius-s*r) / maxRadius)); - - if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2 && yc < miy-1) { // all interpolation pixels inside image - - int r = vignmul*(original->r[yc][xc]*(1.0-Dx)*(1.0-Dy) + original->r[yc][xc+1]*Dx*(1.0-Dy) + original->r[yc+1][xc]*(1.0-Dx)*Dy + original->r[yc+1][xc+1]*Dx*Dy); - int g = vignmul*(original->g[yc][xc]*(1.0-Dx)*(1.0-Dy) + original->g[yc][xc+1]*Dx*(1.0-Dy) + original->g[yc+1][xc]*(1.0-Dx)*Dy + original->g[yc+1][xc+1]*Dx*Dy); - int b = vignmul*(original->b[yc][xc]*(1.0-Dx)*(1.0-Dy) + original->b[yc][xc+1]*Dx*(1.0-Dy) + original->b[yc+1][xc]*(1.0-Dx)*Dy + original->b[yc+1][xc+1]*Dx*Dy); - transformed->r[y][x] = CLIP(r); - transformed->g[y][x] = CLIP(g); - transformed->b[y][x] = CLIP(b); - } - else { // edge pixels - int y1 = (yc>0) ? yc : 0; - if (y1>miy) y1 = miy; - int y2 = (yc0) ? xc : 0; - if (x1>mix) x1 = mix; - int x2 = (xcr[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->r[y1][x2]*Dx*(1.0-Dy) + original->r[y2][x1]*(1.0-Dx)*Dy + original->r[y2][x2]*Dx*Dy); - int g = vignmul*(original->g[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->g[y1][x2]*Dx*(1.0-Dy) + original->g[y2][x1]*(1.0-Dx)*Dy + original->g[y2][x2]*Dx*Dy); - int b = vignmul*(original->b[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->b[y1][x2]*Dx*(1.0-Dy) + original->b[y2][x1]*(1.0-Dx)*Dy + original->b[y2][x2]*Dx*Dy); - transformed->r[y][x] = CLIP(r); - transformed->g[y][x] = CLIP(g); - transformed->b[y][x] = CLIP(b); - } - } - else { - // not valid (source pixel x,y not inside source image, etc.) - transformed->r[y][x] = 0; - transformed->g[y][x] = 0; - transformed->b[y][x] = 0; - } - } - } -} - - -#include "cubintch.cc" -void ImProcFunctions::transform_sep_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { - - int oW = sizes.oW; - int oH = sizes.oH; - int cx = sizes.cx; - int cy = sizes.cy; - int sx = sizes.sx; - int sy = sizes.sy; - - double w2 = (double) oW / 2.0 - 0.5; - double h2 = (double) oH / 2.0 - 0.5; - - double cost = cos(params->rotate.degree * 3.14/180.0); - double sint = sin(params->rotate.degree * 3.14/180.0); - - double max_x = (double) (sx + original->width - 1); - double max_y = (double) (sy + original->height - 1); - double min_x = (double) sx; - double min_y = (double) sy; - - const int n2 = 2; - const int n = 4; - - int mix = original->width - 1; // maximum x-index src - int miy = original->height - 1;// maximum y-index src - int mix2 = mix +1 - n; - int miy2 = miy +1 - n; - - double scale = (oW>oH) ? (double)oW / 2.0 : (double)oH / 2.0 ; - double radius = sqrt( (double)( oW*oW + oH*oH ) ); - radius /= (oWdistortion.amount; - double d = 1.0 - a; - - double cdist[3]; - cdist[0] = params->cacorrection.red; - cdist[1] = 0.0; - cdist[2] = params->cacorrection.blue; - - // magnify image to keep size - double rotmagn = 1.0; - if (params->rotate.fill) { - double beta = atan((double)MIN(oH,oW)/MAX(oW,oH)); - rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); - } - // 1. check upper and lower border - double d1r = rotmagn - a*h2/scale - params->cacorrection.red; - double d2r = rotmagn - a*w2/scale - params->cacorrection.red; - double d3r = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.red; - double dr = MIN(d,MIN(d1r,MIN(d2r,d3r))); - double d1b = rotmagn - a*h2/scale - params->cacorrection.blue; - double d2b = rotmagn - a*w2/scale - params->cacorrection.blue; - double d3b = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.blue; - double db = MIN(d,MIN(d1b,MIN(d2b,d3b))); - double d1g = rotmagn - a*h2/scale; - double d2g = rotmagn - a*w2/scale; - double d3g = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; - double dg = MIN(d,MIN(d1g,MIN(d2g,d3g))); - - d = MIN(dg,MIN(dr,db)); - - unsigned short** chorig[3]; - chorig[0] = original->r; - chorig[1] = original->g; - chorig[2] = original->b; - - unsigned short** chtrans[3]; - chtrans[0] = transformed->r; - chtrans[1] = transformed->g; - chtrans[2] = transformed->b; - - - // auxilary variables for vignetting - double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2 / scale; - - double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; - double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; - - double mul = (1.0-v) / tanh(b); - - // main cycle - double eps = 1e-10; - for (int y=row_from; ywidth; x++) { - double x_d = (double) (x + cx) - w2 ; - - double r = (sqrt(x_d*x_d + y_d*y_d)) / scale; - double s = 10000.0; - if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); - - // Convert only valid pixels - if (valid) { - // Extract integer and fractions of source screen coordinates - int xc = (int) (Dx); Dx -= (double)xc; - int yc = (int) (Dy); Dy -= (double)yc; - int ys = yc +1 - n2 - sy; // smallest y-index used for interpolation - int xs = xc +1 - n2 - sx; // smallest x-index used for interpolation - - if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2) // all interpolation pixels inside image - cubintch (chorig[c], xs, ys, Dx, Dy, &(chtrans[c][y][x]), vignmul); - else {// edge pixels, linear interpolation - int y1 = (yc>0) ? yc : 0; - if (y1>miy) y1 = miy; - int y2 = (yc0) ? xc : 0; - if (x1>mix) x1 = mix; - int x2 = (xc &src, std::vector &red, std::vector &green, std::vector &blue) { - - bool clipresize = true; - bool clipped = false; - - red.clear (); - green.clear (); - blue.clear (); - bool needstransform = 0;// fabs(params->rotate.degree)>1e-15 || fabs(params->distortion.amount)>1e-15 || fabs(params->cacorrection.red)>1e-15 || fabs(params->cacorrection.blue)>1e-15; - if (!needstransform) { - if (clipresize) { - // Apply resizing - if (fabs(params->resize.scale-1.0)>=1e-7) { - for (int i=0; iresize.scale, src[i].y / params->resize.scale)); - green.push_back (Coord2D (src[i].x / params->resize.scale, src[i].y / params->resize.scale)); - blue.push_back (Coord2D (src[i].x / params->resize.scale, src[i].y / params->resize.scale)); - } - for (int i=0; iresize.scale; - double rH = H*params->resize.scale; - double w2 = (double) rW / 2.0 - 0.5; - double h2 = (double) rH / 2.0 - 0.5; - double cost = cos(params->rotate.degree * 3.14/180.0); - double sint = sin(params->rotate.degree * 3.14/180.0); - - double scale = (rW>rH) ? rW / 2.0 : rH / 2.0 ; - double radius = sqrt ((double)(rW*rW + rH*rH )); - radius /= (rWdistortion.amount; - double d = 1.0 - a; - - // magnify image to keep size - double rotmagn = 1.0; - if (params->rotate.fill) { - double beta = atan(MIN(rH,rW)/MAX(rW,rH)); - rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); - } - if (params->cacorrection.red==0 && params->cacorrection.blue==0) { - // 1. check upper and lower border - double d1 = rotmagn - a*h2/scale; - double d2 = rotmagn - a*w2/scale; - double d3 = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; - d = MIN(d,MIN(d1,MIN(d2,d3))); - - for (int i=0; icacorrection.red; - cdist[1] = 0.0; - cdist[2] = params->cacorrection.blue; - - // 1. check upper and lower border - double d1r = rotmagn - a*h2/scale - params->cacorrection.red; - double d2r = rotmagn - a*w2/scale - params->cacorrection.red; - double d3r = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.red; - double dr = MIN(d,MIN(d1r,MIN(d2r,d3r))); - double d1b = rotmagn - a*h2/scale - params->cacorrection.blue; - double d2b = rotmagn - a*w2/scale - params->cacorrection.blue; - double d3b = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.blue; - double db = MIN(d,MIN(d1b,MIN(d2b,d3b))); - double d1g = rotmagn - a*h2/scale; - double d2g = rotmagn - a*w2/scale; - double d3g = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; - double dg = MIN(d,MIN(d1g,MIN(d2g,d3g))); - - d = MIN(dg,MIN(dr,db)); - - for (int i=0; iresize.scale-1.0)>=1e-7) { - for (int i=0; iresize.scale; - red[i].y /= params->resize.scale; - green[i].x /= params->resize.scale; - green[i].y /= params->resize.scale; - blue[i].x /= params->resize.scale; - blue[i].y /= params->resize.scale; - } - } - for (int i=0; i corners (8); - corners[0].set (x1, y1); - corners[1].set (x1, y2); - corners[2].set (x2, y2); - corners[3].set (x2, y1); - corners[4].set ((x1+x2)/2, y1); - corners[5].set ((x1+x2)/2, y2); - corners[6].set (x1, (y1+y2)/2); - corners[7].set (x2, (y1+y2)/2); - - std::vector r, g, b; - - bool result = transCoord (params, W, H, corners, r, g, b); - - std::vector transCorners; - transCorners.insert (transCorners.end(), r.begin(), r.end()); - transCorners.insert (transCorners.end(), g.begin(), g.end()); - transCorners.insert (transCorners.end(), b.begin(), b.end()); - - double x1d = transCorners[0].x; - for (int i=1; ix2d) - x2d = transCorners[i].x; - int x2v = (int)ceil(x2d); - - double y2d = transCorners[0].y; - for (int i=1; iy2d) - y2d = transCorners[i].y; - int y2v = (int)ceil(y2d); - - xv = x1v; - yv = y1v; - wv = x2v - x1v + 1; - hv = y2v - y1v + 1; - - return result; -} - -void ImProcFunctions::transform (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int sx, int sy, int oW, int oH) { - - STemp sizes; - sizes.cx = 0;//cx; - sizes.cy = 0;//cy; - sizes.oW = oW; - sizes.oH = oH; - sizes.sx = 0;//sx; - sizes.sy = 0;//sy; - - if (params->cacorrection.red==0 && params->cacorrection.blue==0) { - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - transform_ (original, transformed, params, sizes, 0, transformed->height); - } - else { - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_sep_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_sep_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - transform_sep_ (original, transformed, params, sizes, 0, transformed->height); - } -} - -void ImProcFunctions::simpltransform (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int sx, int sy, int oW, int oH) { - - STemp sizes; - sizes.cx = 0;//cx; - sizes.cy = 0;//cy; - sizes.oW = oW; - sizes.oH = oH; - sizes.sx = 0;//sx; - sizes.sy = 0;//sy; - - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::simpltransform_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::simpltransform_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - simpltransform_ (original, transformed, params, sizes, 0, transformed->height); -} -/*void ImProcFunctions::transform (Image16* original, Image16* transformed, const ProcParams* params, int ox, int oy) { - - if (!transformed) - return; - - int oW = W, oH = H, tW = W, tH = H; - - double w2 = (double) tW / 2.0 - 0.5; - double h2 = (double) tH / 2.0 - 0.5; - double sw2 = (double) oW / 2.0 - 0.5; - double sh2 = (double) oH / 2.0 - 0.5; - - double cost = cos(params->rotate_fine * 3.14/180.0); - double sint = sin(params->rotate_fine * 3.14/180.0); - - double max_x = (double) oW; - double max_y = (double) oH; - double min_x = 0.0; - double min_y = 0.0; - - const int n2 = 2; - const int n = 4; - - int mix = oW - 1; // maximum x-index src - int miy = oH - 1;// maximum y-index src - int mix2 = mix +1 - n; - int miy2 = miy +1 - n; - - double scale = (tW>tH) ? (double)tW / 2.0 : (double)tH / 2.0 ; - double radius = sqrt( (double)( tW*tW + tH*tH ) ); - radius /= (tWlens_distortion; - - for (int y=0; yheight; y++) { - double y_d = (double) y + oy - h2 ; - - for (int x=0; xwidth; x++) { - double x_d = (double) x + ox - w2 ; - - double r = (sqrt(x_d*x_d + y_d*y_d)) / scale; - double s = 10000.0; - if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); - - // Convert only valid pixels - if (valid) { - // Extract integer and fractions of source screen coordinates - int xc = (int) floor (Dx) ; Dx -= (double)xc; - int yc = (int) floor (Dy) ; Dy -= (double)yc; - int ys = yc +1 - n2 ; // smallest y-index used for interpolation - int xs = xc +1 - n2 ; // smallest x-index used for interpolation - - unsigned short sr[2][2], sg[2][2], sb[2][2]; - - if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2) // all interpolation pixels inside image - cubint (original, xs, ys, Dx, Dy, &(transformed->r[y][x]), &(transformed->g[y][x]), &(transformed->b[y][x])); - else { // edge pixels - transformed->r[y][x] = 0; - transformed->g[y][x] = 0; - transformed->b[y][x] = 0; - } - } - else { - // not valid (source pixel x,y not inside source image, etc.) - transformed->r[y][x] = 0; - transformed->g[y][x] = 0; - transformed->b[y][x] = 0; - } - } - } -}*/ - -void ImProcFunctions::lab2rgb (LabImage* lab, Image8* image) { - - if (settings->dualThreadEnabled) { - - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::lab2rgb_), lab, image, 0, lab->H/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::lab2rgb_), lab, image, lab->H/2, lab->H), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - lab2rgb_ (lab, image, 0, lab->H); -} - -void ImProcFunctions::lab2rgb_ (LabImage* lab, Image8* image, int row_from, int row_to) { - - int X, Y, Z; - unsigned short** nL = lab->L; - short** na = lab->a; - short** nb = lab->b; - int tW = lab->W; - int ix = row_from*tW*3; - - if (monitorTransform) { - short* buffer = new short [3*tW]; - for (int i=row_from; idata + ix, tW); - ix += 3*tW; - } - delete [] buffer; - } - else { - for (int i=row_from; i> 13; - int G = (-8017*X+15697*Y+274*Z) >> 13; - int B = (590*X-1877*Y+11517*Z) >> 13; - - /* copy RGB */ - image->data[ix++] = gamma2curve[CLIP(R)] >> 8; - image->data[ix++] = gamma2curve[CLIP(G)] >> 8; - image->data[ix++] = gamma2curve[CLIP(B)] >> 8; - } - } - } -} - -Image8* ImProcFunctions::lab2rgb (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile) { - - int tW = lab->W; - int tH = lab->H; - - if (cx<0) cx = 0; - if (cy<0) cy = 0; - if (cx+cw>tW) cw = tW-cx; - if (cy+ch>tH) ch = tH-cy; - - Image8* image = new Image8 (cw, ch); - - int X, Y, Z; - int ix = 0; - unsigned short** nL = lab->L; - short** na = lab->a; - short** nb = lab->b; - - cmsHPROFILE oprof = iccStore.getProfile (profile); - - if (oprof) { - cmsHPROFILE iprof = iccStore.getXYZProfile (); - lcmsMutex->lock (); - cmsHTRANSFORM hTransform = cmsCreateTransform (iprof, TYPE_RGB_16, oprof, TYPE_RGB_8, settings->colorimetricIntent, 0); - lcmsMutex->unlock (); - short* buffer = new short [3*cw]; - for (int i=cy; idata + ix, cw); - ix += 3*cw; - } - delete [] buffer; - cmsDeleteTransform(hTransform); - } - else { - for (int i=cy; i> 13; - int G = (-8017*X+15697*Y+274*Z) >> 13; - int B = (590*X-1877*Y+11517*Z) >> 13; - - image->data[ix++] = gamma2curve[CLIP(R)] >> 8; - image->data[ix++] = gamma2curve[CLIP(G)] >> 8; - image->data[ix++] = gamma2curve[CLIP(B)] >> 8; - } - } - } - return image; -} - -Image16* ImProcFunctions::lab2rgb16 (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile) { - - int tW = lab->W; - int tH = lab->H; - - if (cx<0) cx = 0; - if (cy<0) cy = 0; - if (cx+cw>tW) cw = tW-cx; - if (cy+ch>tH) ch = tH-cy; - - Image16* image = new Image16 (cw, ch); - - int X, Y, Z; - int ix = 0; - unsigned short** nL = lab->L; - short** na = lab->a; - short** nb = lab->b; - - cmsHPROFILE oprof = iccStore.getProfile (profile); - - if (oprof) { - for (int i=cy; ir[i-cy]; - short* ya = (short*)image->g[i-cy]; - short* za = (short*)image->b[i-cy]; - for (register int j=cx; jlock (); - cmsHTRANSFORM hTransform = cmsCreateTransform (iprof, TYPE_RGB_16_PLANAR, oprof, TYPE_RGB_16_PLANAR, settings->colorimetricIntent, 0); - lcmsMutex->unlock (); - cmsDoTransform (hTransform, image->data, image->data, image->planestride/2); - cmsDeleteTransform(hTransform); - } - else { - for (int i=cy; i> 13; - int G = (-8017*X+15697*Y+274*Z) >> 13; - int B = (590*X-1877*Y+11517*Z) >> 13; - - image->r[i-cy][j-cx] = gamma2curve[CLIP(R)]; - image->g[i-cy][j-cx] = gamma2curve[CLIP(G)]; - image->b[i-cy][j-cx] = gamma2curve[CLIP(B)]; - } - } - } - return image; -} - -void ImProcFunctions::resize (Image16* src, Image16* dst, ResizeParams params) { - - if (settings->dualThreadEnabled) { - - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::resize_), src, dst, params, 0, dst->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::resize_), src, dst, params, dst->height/2, dst->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else - resize_ (src, dst, params, 0, dst->height); -} - -void ImProcFunctions::resize_ (Image16* src, Image16* dst, ResizeParams params, int row_from, int row_to) { - if(params.method == "Downscale (Better)") { - // small-scale algorithm by Ilia - // provides much better quality on small scales - // calculates mean value over source pixels which current destination pixel covers - // works only for scales < 1 - // for scales ~1 it is analogous to bilinear - // possibly, for even less scale factors (< 0.2 possibly) boundary pixels are not needed, omitting them can give a speedup - // this algorithm is much slower on small factors than others, because it uses all pixels of the SOURCE image - // Ilia Popov ilia_popov@rambler.ru 2010 - - double delta = 1.0 / params.scale; - double k = params.scale * params.scale; - - for(int i = row_from; i < row_to; i++) { - // top and bottom boundary coordinates - double y0 = i * delta; - double y1 = (i + 1) * delta; - - int m0 = y0; - m0 = CLIPTO(m0, 0, src->height-1); - - int m1 = y1; - m1 = CLIPTO(m1, 0, src->height-1); - - // weights of boundary pixels - double wy0 = 1.0 - (y0 - m0); - double wy1 = y1 - m1; - - for(int j = 0; j < dst->width; j++) { - // left and right boundary coordinates - double x0 = j * delta; - double x1 = (j + 1) * delta; - - int n0 = x0; - n0 = CLIPTO(n0, 0, src->width-1); - int n1 = x1; - n1 = CLIPTO(n1, 0, src->width-1); - - double wx0 = 1.0 - (x0 - n0); - double wx1 = x1 - n1; - - double r = 0; - double g = 0; - double b = 0; - - // integration - // corners - r += wy0 * wx0 * src->r[m0][n0] + wy0 * wx1 * src->r[m0][n1] + wy1 * wx0 * src->r[m1][n0] + wy1 * wx1 * src->r[m1][n1]; - g += wy0 * wx0 * src->g[m0][n0] + wy0 * wx1 * src->g[m0][n1] + wy1 * wx0 * src->g[m1][n0] + wy1 * wx1 * src->g[m1][n1]; - b += wy0 * wx0 * src->b[m0][n0] + wy0 * wx1 * src->b[m0][n1] + wy1 * wx0 * src->b[m1][n0] + wy1 * wx1 * src->b[m1][n1]; - - // top and bottom boundaries - for(int n = n0 + 1; n < n1; n++) { - r += wy0 * src->r[m0][n] + wy1 * src->r[m1][n]; - g += wy0 * src->g[m0][n] + wy1 * src->g[m1][n]; - b += wy0 * src->b[m0][n] + wy1 * src->b[m1][n]; - } - - // inner rows - for(int m = m0 + 1; m < m1; m++) { - // left and right boundaries - r += wx0 * src->r[m][n0] + wx1 * src->r[m][n1]; - g += wx0 * src->g[m][n0] + wx1 * src->g[m][n1]; - b += wx0 * src->b[m][n0] + wx1 * src->b[m][n1]; - // inner pixels - for(int n = n0 + 1; n < n1; n++) { - r += src->r[m][n]; - g += src->g[m][n]; - b += src->b[m][n]; - } - } - - // overall weight is equal to the DST pixel area in SRC coordinates - r *= k; - g *= k; - b *= k; - - dst->r[i][j] = CLIP((int)r); - dst->g[i][j] = CLIP((int)g); - dst->b[i][j] = CLIP((int)b); - } - } - - return; - } - - if(params.method == "Downscale (Faster)") - { - // faster version of algo above, does not take into account border pixels, - // which are summed with non-unity weights in slow algo. So, no need - // for weights at all - // Ilia Popov ilia_popov@rambler.ru 5.04.2010 - - double delta = 1.0 / params.scale; - - int p = (int) delta; - - // if actually we are doing upscaling, behave like Nearest - if(p == 0) - p = 1; - - int q = p/2; - - // may cause problems on 32-bit systems on extremely small factors. - // In that case change 1024 to smth less - const int divider = 1024; - - // scaling factor after summation - int k = divider / (p * p); - - for(int i = row_from; i < row_to; i++) { - // y coordinate of center of destination pixel - double y = (i + 0.5) * delta; - - int m0 = (int) (y) - q; - m0 = CLIPTO(m0, 0, src->height-1); - - int m1 = m0 + p; - if(m1 > src->height) { - m1 = src->height; - m0 = m1 - p; - } - m1 = CLIPTO(m1, 0, src->height); - - for(int j = 0; j < dst->width; j++) { - // x coordinate of center of destination pixel - double x = (j + 0.5) * delta; - - int n0 = (int) (x) - q; - n0 = CLIPTO(n0, 0, src->width-1); - - int n1 = n0 + p; - if(n1 > src->width) { - n1 = src->width; - n0 = n1 - p; - } - n1 = CLIPTO(n1, 0, src->width); - - int r = 0; - int g = 0; - int b = 0; - - // integration - for(int m = m0; m < m1; m++) { - for(int n = n0; n < n1; n++) { - r += src->r[m][n]; - g += src->g[m][n]; - b += src->b[m][n]; - } - } - - dst->r[i][j] = CLIP( r * k / divider); - dst->g[i][j] = CLIP( g * k / divider); - dst->b[i][j] = CLIP( b * k / divider); - } - } - return; - } - - - if (params.method.substr(0,7)=="Bicubic") { - double Av = -0.5; - if (params.method=="Bicubic (Sharper)") - Av = -0.75; - else if (params.method=="Bicubic (Softer)") - Av = -0.25; - double wx[4], wy[4]; - for (int i=row_from; iwidth; j++) { - double Dx = j / params.scale; - int xc = (int) Dx; Dx -= (double)xc; - int xs = xc - 1; // smallest x-index used for interpolation - if (ys >= 0 && ys height-3 && xs >= 0 && xs <= src->width-3) { - // compute horizontal weights - double t1 = -Av*(Dx-1.0)*Dx; - double t2 = (3.0-2.0*Dx)*Dx*Dx; - wx[3] = t1*Dx; - wx[2] = t1*(Dx-1.0) + t2; - wx[1] = -t1*Dx + 1.0 - t2; - wx[0] = -t1*(Dx-1.0); - // compute weighted sum - int r = 0; - int g = 0; - int b = 0; -/* r = wx[0]*wy[0]*src->r[ys+0][xs+0] + wx[0]*wy[1]*src->r[ys+1][xs+0] + wx[0]*wy[2]*src->r[ys+2][xs+0] + wx[0]*wy[3]*src->r[ys+3][xs+0] + - wx[1]*wy[0]*src->r[ys+0][xs+1] + wx[1]*wy[1]*src->r[ys+1][xs+1] + wx[1]*wy[2]*src->r[ys+2][xs+1] + wx[1]*wy[3]*src->r[ys+3][xs+1] + - wx[2]*wy[0]*src->r[ys+0][xs+2] + wx[2]*wy[1]*src->r[ys+1][xs+1] + wx[2]*wy[2]*src->r[ys+2][xs+2] + wx[2]*wy[3]*src->r[ys+3][xs+2] + - wx[3]*wy[0]*src->r[ys+0][xs+3] + wx[3]*wy[1]*src->r[ys+1][xs+1] + wx[3]*wy[2]*src->r[ys+2][xs+3] + wx[3]*wy[3]*src->r[ys+3][xs+3]; - g = wx[0]*wy[0]*src->g[ys+0][xs+0] + wx[0]*wy[1]*src->g[ys+1][xs+0] + wx[0]*wy[2]*src->g[ys+2][xs+0] + wx[0]*wy[3]*src->g[ys+3][xs+0] + - wx[1]*wy[0]*src->g[ys+0][xs+1] + wx[1]*wy[1]*src->g[ys+1][xs+1] + wx[1]*wy[2]*src->g[ys+2][xs+1] + wx[1]*wy[3]*src->g[ys+3][xs+1] + - wx[2]*wy[0]*src->g[ys+0][xs+2] + wx[2]*wy[1]*src->g[ys+1][xs+1] + wx[2]*wy[2]*src->g[ys+2][xs+2] + wx[2]*wy[3]*src->g[ys+3][xs+2] + - wx[3]*wy[0]*src->g[ys+0][xs+3] + wx[3]*wy[1]*src->g[ys+1][xs+1] + wx[3]*wy[2]*src->g[ys+2][xs+3] + wx[3]*wy[3]*src->g[ys+3][xs+3]; - b = wx[0]*wy[0]*src->b[ys+0][xs+0] + wx[0]*wy[1]*src->b[ys+1][xs+0] + wx[0]*wy[2]*src->b[ys+2][xs+0] + wx[0]*wy[3]*src->b[ys+3][xs+0] + - wx[1]*wy[0]*src->b[ys+0][xs+1] + wx[1]*wy[1]*src->b[ys+1][xs+1] + wx[1]*wy[2]*src->b[ys+2][xs+1] + wx[1]*wy[3]*src->b[ys+3][xs+1] + - wx[2]*wy[0]*src->b[ys+0][xs+2] + wx[2]*wy[1]*src->b[ys+1][xs+1] + wx[2]*wy[2]*src->b[ys+2][xs+2] + wx[2]*wy[3]*src->b[ys+3][xs+2] + - wx[3]*wy[0]*src->b[ys+0][xs+3] + wx[3]*wy[1]*src->b[ys+1][xs+1] + wx[3]*wy[2]*src->b[ys+2][xs+3] + wx[3]*wy[3]*src->b[ys+3][xs+3];*/ - for (int x=0; x<4; x++) - for (int y=0; y<4; y++) { - double w = wx[x]*wy[y]; - r += w*src->r[ys+y][xs+x]; - g += w*src->g[ys+y][xs+x]; - b += w*src->b[ys+y][xs+x]; - } - dst->r[i][j] = CLIP(r); - dst->g[i][j] = CLIP(g); - dst->b[i][j] = CLIP(b); - } - else { - xc = CLIPTO(xc, 0, src->width-1); - yc = CLIPTO(yc, 0, src->height-1); - int nx = xc + 1; - if (nx>=src->width) - nx = xc; - int ny = yc + 1; - if (ny>=src->height) - ny = yc; - dst->r[i][j] = (1-Dx)*(1-Dy)*src->r[yc][xc] + (1-Dx)*Dy*src->r[ny][xc] + Dx*(1-Dy)*src->r[yc][nx] + Dx*Dy*src->r[ny][nx]; - dst->g[i][j] = (1-Dx)*(1-Dy)*src->g[yc][xc] + (1-Dx)*Dy*src->g[ny][xc] + Dx*(1-Dy)*src->g[yc][nx] + Dx*Dy*src->g[ny][nx]; - dst->b[i][j] = (1-Dx)*(1-Dy)*src->b[yc][xc] + (1-Dx)*Dy*src->b[ny][xc] + Dx*(1-Dy)*src->b[yc][nx] + Dx*Dy*src->b[ny][nx]; - } - } - } - } - else if (params.method=="Bilinear") { - for (int i=row_from; iheight-1); - double dy = i/params.scale - sy; - int ny = sy+1; - if (ny>=src->height) - ny = sy; - for (int j=0; jwidth; j++) { - int sx = j/params.scale; - sx = CLIPTO(sx, 0, src->width-1); - double dx = j/params.scale - sx; - int nx = sx+1; - if (nx>=src->width) - nx = sx; - dst->r[i][j] = (1-dx)*(1-dy)*src->r[sy][sx] + (1-dx)*dy*src->r[ny][sx] + dx*(1-dy)*src->r[sy][nx] + dx*dy*src->r[ny][nx]; - dst->g[i][j] = (1-dx)*(1-dy)*src->g[sy][sx] + (1-dx)*dy*src->g[ny][sx] + dx*(1-dy)*src->g[sy][nx] + dx*dy*src->g[ny][nx]; - dst->b[i][j] = (1-dx)*(1-dy)*src->b[sy][sx] + (1-dx)*dy*src->b[ny][sx] + dx*(1-dy)*src->b[sy][nx] + dx*dy*src->b[ny][nx]; - } - } - } - else { - for (int i=row_from; iheight-1); - for (int j=0; jwidth; j++) { - int sx = j/params.scale; - sx = CLIPTO(sx, 0, src->width-1); - dst->r[i][j] = src->r[sy][sx]; - dst->g[i][j] = src->g[sy][sx]; - dst->b[i][j] = src->b[sy][sx]; - } - } - } -} - void ImProcFunctions::getAutoExp (unsigned int* histogram, int histcompr, double expcomp, double clip, double& br, int& bl) { double sum = 0; @@ -2242,18 +481,6 @@ void ImProcFunctions::getAutoExp (unsigned int* histogram, int histcompr, doubl br = log(65535.0 / (awg-bl)) / log(2.0); if (br<0) br = 0; - -//printf ("br=%g, bl=%d, %g\n", br, bl, expcomp); - -/* - if (shc #include #include +#include namespace rtengine { using namespace procparams; -class LabImage { - private: - bool fromImage; - - public: - int W, H; - unsigned short** L; - short** a; - short** b; - - LabImage (int w, int h); - LabImage (Image16* im); - ~LabImage (); -}; - class ImProcFunctions { - protected: - struct STemp { - int cx, cy, sx, sy, oW, oH; - }; - cmsHTRANSFORM monitorTransform; - - void transform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); - void simpltransform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); - void vignetting_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); - void transform_sep_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); - void rgbProc_ (Image16* working, LabImage* lab, const ProcParams* params, int* tonecurve, SHMap* shmap, int row_from, int row_to); - void lab2rgb_ (LabImage* lab, Image8* image, int row_from, int row_to); - void colorCurve_ (LabImage* lold, LabImage* lnew, const ProcParams* params, int row_from, int row_to, double* cmultiplier); - void sharpenHaloCtrl (LabImage* lab, const ProcParams* params, unsigned short** blurmap, unsigned short** base, int W, int row_from, int row_to); - void firstAnalysis_ (Image16* original, Glib::ustring wprofile, unsigned int* histogram, int* chroma_radius, int row_from, int row_to); - void resize_ (Image16* src, Image16* dst, ResizeParams params, int row_from, int row_to); - void damping_ (float** aI, unsigned short** aO, float damping, int W, int rowfrom, int rowto); - - public: - static int* cacheL; static int* cachea; static int* cacheb; static int* xcache; static int* ycache; static int* zcache; - + static unsigned short gamma2curve[65536]; + + struct STemp { + int cx, cy, sx, sy, oW, oH; + }; + cmsHTRANSFORM monitorTransform; + int chroma_scale; int chroma_radius; + const ProcParams* params; + double scale; + bool multiThread; + + void transform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); + void simpltransform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); + void vignetting_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); + void transform_sep_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to); + void sharpenHaloCtrl (LabImage* lab, unsigned short** blurmap, unsigned short** base, int W, int H); + void firstAnalysis_ (Image16* original, Glib::ustring wprofile, unsigned int* histogram, int* chroma_radius, int row_from, int row_to); + void dcdamping (float** aI, unsigned short** aO, float damping, int W, int H); + + public: + double lumimul[3]; - static unsigned short gamma2curve[65536]; - static void initCache (); - ImProcFunctions () : monitorTransform(NULL) {} - void release (); - + ImProcFunctions (const ProcParams* iparams, bool imultiThread=true) + : monitorTransform(NULL), params(iparams), scale(1), multiThread(imultiThread) {} + ~ImProcFunctions (); + + void setScale (double iscale); void firstAnalysis (Image16* working, const ProcParams* params, unsigned int* vhist16, double gamma); - void rgbProc (Image16* working, LabImage* lab, const ProcParams* params, int* tonecurve, SHMap* shmap); + void rgbProc (Image16* working, LabImage* lab, int* tonecurve, SHMap* shmap); void luminanceCurve (LabImage* lold, LabImage* lnew, int* curve, int row_from, int row_to); - void colorCurve (LabImage* lold, LabImage* lnew, const ProcParams* params); - void sharpening (LabImage* lab, const ProcParams* params, double scale, unsigned short** buffer); - void lumadenoise (LabImage* lab, const ProcParams* params, double scale, int** buffer); - void colordenoise (LabImage* lab, const ProcParams* params, double scale, int** buffer); + void colorCurve (LabImage* lold, LabImage* lnew); + void sharpening (LabImage* lab, unsigned short** buffer); + void lumadenoise (LabImage* lab, int** buffer); + void colordenoise (LabImage* lab, int** buffer); void transform (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int sx, int sy, int oW, int oH); void simpltransform (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int sx, int sy, int oW, int oH); void vignetting (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int oW, int oH); void lab2rgb (LabImage* lab, Image8* image); - void resize (Image16* src, Image16* dst, ResizeParams params); + void resize (Image16* src, Image16* dst); - bool transCoord (const ProcParams* params, int W, int H, int x, int y, int w, int h, int& xv, int& yv, int& wv, int& hv); - bool transCoord (const ProcParams* params, int W, int H, std::vector &src, std::vector &red, std::vector &green, std::vector &blue); - void deconvsharpening(LabImage* lab, const ProcParams* params, double scale, unsigned short** buffer); + void deconvsharpening(LabImage* lab, unsigned short** buffer); Image8* lab2rgb (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile); Image16* lab2rgb16 (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile); + static bool transCoord (const ProcParams* params, int W, int H, int x, int y, int w, int h, int& xv, int& yv, int& wv, int& hv); + static bool transCoord (const ProcParams* params, int W, int H, std::vector &src, std::vector &red, std::vector &green, std::vector &blue); static void getAutoExp (unsigned int* histogram, int histcompr, double expcomp, double clip, double& br, int& bl); }; }; diff --git a/rtengine/iplab2rgb.cc b/rtengine/iplab2rgb.cc new file mode 100644 index 000000000..c7425e434 --- /dev/null +++ b/rtengine/iplab2rgb.cc @@ -0,0 +1,259 @@ +/* + * 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 +#include +#include +#include + +namespace rtengine { + +#undef CLIP +#undef CLIPTO +#undef CMAXVAL + +#define CMAXVAL 0xffff +#define CLIP(a) ((a)>0?((a)(b)?((a)<(c)?(a):(c)):(b)) + +extern const Settings* settings; + +void ImProcFunctions::lab2rgb (LabImage* lab, Image8* image) { + + if (monitorTransform) { + int ix = 0; + short* buffer = new short [3*lab->W]; + for (int i=0; iH; i++) { + unsigned short* rL = lab->L[i]; + short* ra = lab->a[i]; + short* rb = lab->b[i]; + int iy = 0; + for (int j=0; jW; j++) { + + int y_ = rL[j]; + int x_ = rL[j]+10486+ra[j]*152/chroma_scale+141556; + int z_ = rL[j]+10486-rb[j]*380/chroma_scale+369619; + + x_ = CLIPTO(x_,0,369820); + y_ = CLIPTO(y_,0,825745); + + y_ = ycache[y_]; + x_ = xcache[x_]; + z_ = zcache[z_]; + + buffer[iy++] = CLIP(x_); + buffer[iy++] = CLIP(y_); + buffer[iy++] = CLIP(z_); + } + cmsDoTransform (monitorTransform, buffer, image->data + ix, lab->W); + ix += 3*lab->W; + } + delete [] buffer; + } + else { + #pragma omp parallel for if (multiThread) + for (int i=0; iH; i++) { + unsigned short* rL = lab->L[i]; + short* ra = lab->a[i]; + short* rb = lab->b[i]; + int ix = 3*i*lab->W; + for (int j=0; jW; j++) { + + int y_ = rL[j]; + int x_ = rL[j]+10486+ra[j]*152/chroma_scale+141556; + int z_ = rL[j]+10486-rb[j]*380/chroma_scale+369619; + + x_ = CLIPTO(x_,0,369820); + y_ = CLIPTO(y_,0,825745); + + y_ = ycache[y_]; + x_ = xcache[x_]; + z_ = zcache[z_]; + + /* XYZ-D50 to RGB */ + int R = (25689*x_-13261*y_-4022*z_) >> 13; + int G = (-8017*x_+15697*y_+274*z_) >> 13; + int B = (590*x_-1877*y_+11517*z_) >> 13; + + /* copy RGB */ + image->data[ix++] = gamma2curve[CLIP(R)] >> 8; + image->data[ix++] = gamma2curve[CLIP(G)] >> 8; + image->data[ix++] = gamma2curve[CLIP(B)] >> 8; + } + } + } +} + +Image8* ImProcFunctions::lab2rgb (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile) { + + 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); + + cmsHPROFILE oprof = iccStore.getProfile (profile); + + if (oprof) { + cmsHPROFILE iprof = iccStore.getXYZProfile (); + lcmsMutex->lock (); + cmsHTRANSFORM hTransform = cmsCreateTransform (iprof, TYPE_RGB_16, oprof, TYPE_RGB_8, settings->colorimetricIntent, 0); + lcmsMutex->unlock (); + int ix = 0; + short* buffer = new short [3*cw]; + for (int i=cy; iL[i]; + short* ra = lab->a[i]; + short* rb = lab->b[i]; + int iy = 0; + for (int j=cx; jdata + ix, cw); + ix += 3*cw; + } + delete [] buffer; + cmsDeleteTransform(hTransform); + } + else { + #pragma omp parallel for if (multiThread) + for (int i=cy; iL[i]; + short* ra = lab->a[i]; + short* rb = lab->b[i]; + int ix = 3*i*cw; + for (int j=cx; j> 13; + int G = (-8017*x_+15697*y_+274*z_) >> 13; + int B = (590*x_-1877*y_+11517*z_) >> 13; + + image->data[ix++] = gamma2curve[CLIP(R)] >> 8; + image->data[ix++] = gamma2curve[CLIP(G)] >> 8; + image->data[ix++] = gamma2curve[CLIP(B)] >> 8; + } + } + } + return image; +} + +Image16* ImProcFunctions::lab2rgb16 (LabImage* lab, int cx, int cy, int cw, int ch, Glib::ustring profile) { + + 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; + + Image16* image = new Image16 (cw, ch); + + cmsHPROFILE oprof = iccStore.getProfile (profile); + + if (oprof) { + #pragma omp parallel for if (multiThread) + for (int i=cy; iL[i]; + short* ra = lab->a[i]; + short* rb = lab->b[i]; + short* xa = (short*)image->r[i-cy]; + short* ya = (short*)image->g[i-cy]; + short* za = (short*)image->b[i-cy]; + for (int j=cx; jlock (); + cmsHTRANSFORM hTransform = cmsCreateTransform (iprof, TYPE_RGB_16_PLANAR, oprof, TYPE_RGB_16_PLANAR, settings->colorimetricIntent, 0); + lcmsMutex->unlock (); + cmsDoTransform (hTransform, image->data, image->data, image->planestride/2); + cmsDeleteTransform(hTransform); + } + else { + #pragma omp parallel for if (multiThread) + for (int i=cy; iL[i]; + short* ra = lab->a[i]; + short* rb = lab->b[i]; + for (int j=cx; j> 13; + int G = (-8017*x_+15697*y_+274*z_) >> 13; + int B = (590*x_-1877*y_+11517*z_) >> 13; + + image->r[i-cy][j-cx] = gamma2curve[CLIP(R)]; + image->g[i-cy][j-cx] = gamma2curve[CLIP(G)]; + image->b[i-cy][j-cx] = gamma2curve[CLIP(B)]; + } + } + } + return image; +} + +} diff --git a/rtengine/ipresize.cc b/rtengine/ipresize.cc new file mode 100644 index 000000000..c70479f57 --- /dev/null +++ b/rtengine/ipresize.cc @@ -0,0 +1,293 @@ +/* + * 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 +#include +#include + +namespace rtengine { + +#undef CLIP +#undef CLIPTO +#undef CMAXVAL + +#define CMAXVAL 0xffff +#define CLIP(a) ((a)>0?((a)(b)?((a)<(c)?(a):(c)):(b)) + +void ImProcFunctions::resize (Image16* src, Image16* dst) { + + if(params->resize.method == "Downscale (Better)") { + // small-scale algorithm by Ilia + // provides much better quality on small scales + // calculates mean value over source pixels which current destination pixel covers + // works only for scales < 1 + // for scales ~1 it is analogous to bilinear + // possibly, for even less scale factors (< 0.2 possibly) boundary pixels are not needed, omitting them can give a speedup + // this algorithm is much slower on small factors than others, because it uses all pixels of the SOURCE image + // Ilia Popov ilia_popov@rambler.ru 2010 + + double delta = 1.0 / params->resize.scale; + double k = params->resize.scale * params->resize.scale; + + #pragma omp parallel for if (multiThread) + for(int i = 0; i < dst->height; i++) { + // top and bottom boundary coordinates + double y0 = i * delta; + double y1 = (i + 1) * delta; + + int m0 = y0; + m0 = CLIPTO(m0, 0, src->height-1); + + int m1 = y1; + m1 = CLIPTO(m1, 0, src->height-1); + + // weights of boundary pixels + double wy0 = 1.0 - (y0 - m0); + double wy1 = y1 - m1; + + for(int j = 0; j < dst->width; j++) { + // left and right boundary coordinates + double x0 = j * delta; + double x1 = (j + 1) * delta; + + int n0 = x0; + n0 = CLIPTO(n0, 0, src->width-1); + int n1 = x1; + n1 = CLIPTO(n1, 0, src->width-1); + + double wx0 = 1.0 - (x0 - n0); + double wx1 = x1 - n1; + + double r = 0; + double g = 0; + double b = 0; + + // integration + // corners + r += wy0 * wx0 * src->r[m0][n0] + wy0 * wx1 * src->r[m0][n1] + wy1 * wx0 * src->r[m1][n0] + wy1 * wx1 * src->r[m1][n1]; + g += wy0 * wx0 * src->g[m0][n0] + wy0 * wx1 * src->g[m0][n1] + wy1 * wx0 * src->g[m1][n0] + wy1 * wx1 * src->g[m1][n1]; + b += wy0 * wx0 * src->b[m0][n0] + wy0 * wx1 * src->b[m0][n1] + wy1 * wx0 * src->b[m1][n0] + wy1 * wx1 * src->b[m1][n1]; + + // top and bottom boundaries + for(int n = n0 + 1; n < n1; n++) { + r += wy0 * src->r[m0][n] + wy1 * src->r[m1][n]; + g += wy0 * src->g[m0][n] + wy1 * src->g[m1][n]; + b += wy0 * src->b[m0][n] + wy1 * src->b[m1][n]; + } + + // inner rows + for(int m = m0 + 1; m < m1; m++) { + // left and right boundaries + r += wx0 * src->r[m][n0] + wx1 * src->r[m][n1]; + g += wx0 * src->g[m][n0] + wx1 * src->g[m][n1]; + b += wx0 * src->b[m][n0] + wx1 * src->b[m][n1]; + // inner pixels + for(int n = n0 + 1; n < n1; n++) { + r += src->r[m][n]; + g += src->g[m][n]; + b += src->b[m][n]; + } + } + // overall weight is equal to the DST pixel area in SRC coordinates + r *= k; + g *= k; + b *= k; + + dst->r[i][j] = CLIP((int)r); + dst->g[i][j] = CLIP((int)g); + dst->b[i][j] = CLIP((int)b); + } + } + return; + } + + if(params->resize.method == "Downscale (Faster)") + { + // faster version of algo above, does not take into account border pixels, + // which are summed with non-unity weights in slow algo. So, no need + // for weights at all + // Ilia Popov ilia_popov@rambler.ru 5.04.2010 + + double delta = 1.0 / params->resize.scale; + + int p = (int) delta; + + // if actually we are doing upscaling, behave like Nearest + if(p == 0) + p = 1; + + int q = p/2; + + // may cause problems on 32-bit systems on extremely small factors. + // In that case change 1024 to smth less + const int divider = 1024; + + // scaling factor after summation + int k = divider / (p * p); + + #pragma omp parallel for if (multiThread) + for(int i = 0; i < dst->height; i++) { + // y coordinate of center of destination pixel + double y = (i + 0.5) * delta; + + int m0 = (int) (y) - q; + m0 = CLIPTO(m0, 0, src->height-1); + + int m1 = m0 + p; + if(m1 > src->height) { + m1 = src->height; + m0 = m1 - p; + } + m1 = CLIPTO(m1, 0, src->height); + + for(int j = 0; j < dst->width; j++) { + // x coordinate of center of destination pixel + double x = (j + 0.5) * delta; + + int n0 = (int) (x) - q; + n0 = CLIPTO(n0, 0, src->width-1); + + int n1 = n0 + p; + if(n1 > src->width) { + n1 = src->width; + n0 = n1 - p; + } + n1 = CLIPTO(n1, 0, src->width); + + int r = 0; + int g = 0; + int b = 0; + + // integration + for(int m = m0; m < m1; m++) { + for(int n = n0; n < n1; n++) { + r += src->r[m][n]; + g += src->g[m][n]; + b += src->b[m][n]; + } + } + dst->r[i][j] = CLIP( r * k / divider); + dst->g[i][j] = CLIP( g * k / divider); + dst->b[i][j] = CLIP( b * k / divider); + } + } + return; + } + if (params->resize.method.substr(0,7)=="Bicubic") { + double Av = -0.5; + if (params->resize.method=="Bicubic (Sharper)") + Av = -0.75; + else if (params->resize.method=="Bicubic (Softer)") + Av = -0.25; + #pragma omp parallel for if (multiThread) + for (int i=0; iheight; i++) { + double wx[4], wy[4]; + double Dy = i / params->resize.scale; + int yc = (int) Dy; Dy -= (double)yc; + int ys = yc - 1; // smallest y-index used for interpolation + // compute vertical weights + double t1y = -Av*(Dy-1.0)*Dy; + double t2y = (3.0-2.0*Dy)*Dy*Dy; + wy[3] = t1y*Dy; + wy[2] = t1y*(Dy-1.0) + t2y; + wy[1] = -t1y*Dy + 1.0 - t2y; + wy[0] = -t1y*(Dy-1.0); + for (int j=0; jwidth; j++) { + double Dx = j / params->resize.scale; + int xc = (int) Dx; Dx -= (double)xc; + int xs = xc - 1; // smallest x-index used for interpolation + if (ys >= 0 && ys height-3 && xs >= 0 && xs <= src->width-3) { + // compute horizontal weights + double t1 = -Av*(Dx-1.0)*Dx; + double t2 = (3.0-2.0*Dx)*Dx*Dx; + wx[3] = t1*Dx; + wx[2] = t1*(Dx-1.0) + t2; + wx[1] = -t1*Dx + 1.0 - t2; + wx[0] = -t1*(Dx-1.0); + // compute weighted sum + int r = 0; + int g = 0; + int b = 0; + for (int x=0; x<4; x++) + for (int y=0; y<4; y++) { + double w = wx[x]*wy[y]; + r += w*src->r[ys+y][xs+x]; + g += w*src->g[ys+y][xs+x]; + b += w*src->b[ys+y][xs+x]; + } + dst->r[i][j] = CLIP(r); + dst->g[i][j] = CLIP(g); + dst->b[i][j] = CLIP(b); + } + else { + xc = CLIPTO(xc, 0, src->width-1); + yc = CLIPTO(yc, 0, src->height-1); + int nx = xc + 1; + if (nx>=src->width) + nx = xc; + int ny = yc + 1; + if (ny>=src->height) + ny = yc; + dst->r[i][j] = (1-Dx)*(1-Dy)*src->r[yc][xc] + (1-Dx)*Dy*src->r[ny][xc] + Dx*(1-Dy)*src->r[yc][nx] + Dx*Dy*src->r[ny][nx]; + dst->g[i][j] = (1-Dx)*(1-Dy)*src->g[yc][xc] + (1-Dx)*Dy*src->g[ny][xc] + Dx*(1-Dy)*src->g[yc][nx] + Dx*Dy*src->g[ny][nx]; + dst->b[i][j] = (1-Dx)*(1-Dy)*src->b[yc][xc] + (1-Dx)*Dy*src->b[ny][xc] + Dx*(1-Dy)*src->b[yc][nx] + Dx*Dy*src->b[ny][nx]; + } + } + } + } + else if (params->resize.method=="Bilinear") { + #pragma omp parallel for if (multiThread) + for (int i=0; iheight; i++) { + int sy = i/params->resize.scale; + sy = CLIPTO(sy, 0, src->height-1); + double dy = i/params->resize.scale - sy; + int ny = sy+1; + if (ny>=src->height) + ny = sy; + for (int j=0; jwidth; j++) { + int sx = j/params->resize.scale; + sx = CLIPTO(sx, 0, src->width-1); + double dx = j/params->resize.scale - sx; + int nx = sx+1; + if (nx>=src->width) + nx = sx; + dst->r[i][j] = (1-dx)*(1-dy)*src->r[sy][sx] + (1-dx)*dy*src->r[ny][sx] + dx*(1-dy)*src->r[sy][nx] + dx*dy*src->r[ny][nx]; + dst->g[i][j] = (1-dx)*(1-dy)*src->g[sy][sx] + (1-dx)*dy*src->g[ny][sx] + dx*(1-dy)*src->g[sy][nx] + dx*dy*src->g[ny][nx]; + dst->b[i][j] = (1-dx)*(1-dy)*src->b[sy][sx] + (1-dx)*dy*src->b[ny][sx] + dx*(1-dy)*src->b[sy][nx] + dx*dy*src->b[ny][nx]; + } + } + } + else { + #pragma omp parallel for if (multiThread) + for (int i=0; iheight; i++) { + int sy = i/params->resize.scale; + sy = CLIPTO(sy, 0, src->height-1); + for (int j=0; jwidth; j++) { + int sx = j/params->resize.scale; + sx = CLIPTO(sx, 0, src->width-1); + dst->r[i][j] = src->r[sy][sx]; + dst->g[i][j] = src->g[sy][sx]; + dst->b[i][j] = src->b[sy][sx]; + } + } + } +} + +} diff --git a/rtengine/ipsharpen.cc b/rtengine/ipsharpen.cc new file mode 100644 index 000000000..7cc4ddfaa --- /dev/null +++ b/rtengine/ipsharpen.cc @@ -0,0 +1,201 @@ +/* + * 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 +#include +#include +#include +#include + +namespace rtengine { + +#undef CLIP +#undef CMAXVAL +#undef ABS + +#define CMAXVAL 0xffff +#define CLIP(a) ((a)>0?((a)sharpening.enabled==false || params->sharpening.deconvamount<1) + return; + + int W = lab->W, H = lab->H; + + float** tmpI = new float*[H]; + for (int i=0; iL[i][j]; + } + + float** tmp = (float**)b2; + + AlignedBuffer* buffer = new AlignedBuffer (MAX(W,H)*omp_get_max_threads()); + + float damping = params->sharpening.deconvdamping / 5.0; + bool needdamp = params->sharpening.deconvdamping > 0; + for (int k=0; ksharpening.deconviter; k++) { + + // apply blur function (gaussian blur) + gaussHorizontal (tmpI, tmp, buffer, W, H, params->sharpening.deconvradius / scale, multiThread); + gaussVertical (tmp, tmp, buffer, W, H, params->sharpening.deconvradius / scale, multiThread); + + if (!needdamp) { + #pragma omp parallel for if (multiThread) + for (int i=0; i0) + tmp[i][j] = (float)lab->L[i][j] / tmp[i][j]; + } + else + dcdamping (tmp, lab->L, damping, W, H); + + gaussHorizontal (tmp, tmp, buffer, W, H, params->sharpening.deconvradius / scale, multiThread); + gaussVertical (tmp, tmp, buffer, W, H, params->sharpening.deconvradius / scale, multiThread); + + #pragma omp parallel for if (multiThread) + for (int i=0; iL[i][j] = lab->L[i][j]*(100-params->sharpening.deconvamount) / 100 + (int)CLIP(tmpI[i][j])*params->sharpening.deconvamount / 100; + + for (int i=0; isharpening.method=="rld") { + deconvsharpening (lab, b2); + return; + } + + if (params->sharpening.enabled==false || params->sharpening.amount<1 || lab->W<8 || lab->H<8) + return; + + int W = lab->W, H = lab->H; + unsigned short** b3; + + AlignedBuffer* buffer = new AlignedBuffer (MAX(W,H)*omp_get_max_threads()); + + if (params->sharpening.edgesonly==false) { + + gaussHorizontal (lab->L, b2, buffer, W, H, params->sharpening.radius / scale, multiThread); + gaussVertical (b2, b2, buffer, W, H, params->sharpening.radius / scale, multiThread); + } + else { + b3 = new unsigned short*[H]; + for (int i=0; i (lab->L, (unsigned short**)b3, b2, W, H, params->sharpening.edges_radius / scale, params->sharpening.edges_tolerance, multiThread); + gaussHorizontal (b3, b2, buffer, W, H, params->sharpening.radius / scale, multiThread); + gaussVertical (b2, b2, buffer, W, H, params->sharpening.radius / scale, multiThread); + } + delete buffer; + + unsigned short** base = lab->L; + if (params->sharpening.edgesonly) + base = b3; + + if (params->sharpening.halocontrol==false) { + #pragma omp parallel for if (multiThread) + for (int i=0; iparams->sharpening.threshold) { + int val = lab->L[i][j] + params->sharpening.amount * diff / 100; + lab->L[i][j] = CLIP(val); + } + } + } + else + sharpenHaloCtrl (lab, b2, base, W, H); + + if (params->sharpening.edgesonly) { + for (int i=0; isharpening.halocontrol_amount); + unsigned short** nL = base; + #pragma omp parallel for if (multiThread) + for (int i=2; i params->sharpening.threshold) { + // compute maximum/minimum in a delta environment + np1 = 2*(nL[i-2][j] + nL[i-2][j+1] + nL[i-2][j+2] + nL[i-1][j] + nL[i-1][j+1] + nL[i-1][j+2] + nL[i][j] + nL[i][j+1] + nL[i][j+2]) / 27 + nL[i-1][j+1] / 3; + np2 = 2*(nL[i-1][j] + nL[i-1][j+1] + nL[i-1][j+2] + nL[i][j] + nL[i][j+1] + nL[i][j+2] + nL[i+1][j] + nL[i+1][j+1] + nL[i+1][j+2]) / 27 + nL[i][j+1] / 3; + np3 = 2*(nL[i][j] + nL[i][j+1] + nL[i][j+2] + nL[i+1][j] + nL[i+1][j+1] + nL[i+1][j+2] + nL[i+2][j] + nL[i+2][j+1] + nL[i+2][j+2]) / 27 + nL[i+1][j+1] / 3; + MINMAX3(np1,np2,np3,maxn,minn); + MAX3(max1,max2,maxn,max); + MIN3(min1,min2,minn,min); + max1 = max2; max2 = maxn; + min1 = min2; min2 = minn; + if (max < lab->L[i][j]) + max = lab->L[i][j]; + if (min > lab->L[i][j]) + min = lab->L[i][j]; + int val = lab->L[i][j] + params->sharpening.amount * diff / 100; + int newL = CLIP(val); + // applying halo control + if (newL > max) + newL = max + (newL-max) * scale / 10000; + else if (newLL[i][j] = newL; + } + } + } +} + +} diff --git a/rtengine/iptransform.cc b/rtengine/iptransform.cc new file mode 100644 index 000000000..f50c5b6ba --- /dev/null +++ b/rtengine/iptransform.cc @@ -0,0 +1,821 @@ +/* + * 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 +#include + +namespace rtengine { + +#undef CMAXVAL +#undef MAX +#undef MIN +#undef CLIP +#undef CLIPTOC + +#define CMAXVAL 0xffff +#define MAX(a,b) ((a)<(b)?(b):(a)) +#define MIN(a,b) ((a)>(b)?(b):(a)) +#define CLIP(a) ((a)>0?((a)=(b)?((a)<=(c)?(a):((c),d=true)):((b),d=true)) + +extern const Settings* settings; + +void ImProcFunctions::vignetting_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { + + int oW = sizes.oW; + int oH = sizes.oH; + int cx = sizes.cx; + int cy = sizes.cy; + + double w2 = (double) oW / 2.0 - 0.5; + double h2 = (double) oH / 2.0 - 0.5; + + double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2; + + double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; + double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; + + double mul = (1.0-v) / tanh(b); + + int val; + for (int y=row_from; ywidth; x++) { + double x_d = (double) (x + cx) - w2 ; + double r = sqrt(x_d*x_d + y_d*y_d); + double vign = v + mul * tanh (b*(maxRadius-r) / maxRadius); + val = original->r[y][x] / vign; + transformed->r[y][x] = CLIP(val); + val = original->g[y][x] / vign; + transformed->g[y][x] = CLIP(val); + val = original->b[y][x] / vign; + transformed->b[y][x] = CLIP(val); + } + } +} + +void ImProcFunctions::vignetting (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int oW, int oH) { + + STemp sizes; + sizes.cx = cx; + sizes.cy = cy; + sizes.oW = oW; + sizes.oH = oH; + + if (settings->dualThreadEnabled) { + Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::vignetting_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::vignetting_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + thread1->join (); + thread2->join (); + } + else + vignetting_ (original, transformed, params, sizes, 0, transformed->height); +} + +#include "cubint.cc" +void ImProcFunctions::transform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { + + int oW = sizes.oW; + int oH = sizes.oH; + int cx = sizes.cx; + int cy = sizes.cy; + int sx = sizes.sx; + int sy = sizes.sy; + + double w2 = (double) oW / 2.0 - 0.5; + double h2 = (double) oH / 2.0 - 0.5; + + double cost = cos(params->rotate.degree * 3.14/180.0); + double sint = sin(params->rotate.degree * 3.14/180.0); + + double max_x = (double) (sx + original->width - 1); + double max_y = (double) (sy + original->height - 1); + double min_x = (double) sx; + double min_y = (double) sy; + + const int n2 = 2; + const int n = 4; + + int mix = original->width - 1; // maximum x-index src + int miy = original->height - 1;// maximum y-index src + int mix2 = mix +1 - n; + int miy2 = miy +1 - n; + + double scale = (oW>oH) ? (double)oW / 2.0 : (double)oH / 2.0 ; + double radius = sqrt( (double)( oW*oW + oH*oH ) ); + radius /= (oWdistortion.amount; + + double d = 1.0 - a; + + // magnify image to keep size + double rotmagn = 1.0; + if (params->rotate.fill) { + double beta = atan((double)MIN(oH,oW)/MAX(oW,oH)); + rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); + } + // 1. check upper and lower border + double d1 = rotmagn - a*h2/scale; + double d2 = rotmagn - a*w2/scale; + double d3 = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; + d = MIN(d,MIN(d1,MIN(d2,d3))); + + // auxilary variables for vignetting + double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2 / scale; + + double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; + double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; + + double mul = (1.0-v) / tanh(b); + + // main cycle + double eps = 1e-10; + bool calc_r=( (fabs(a)>eps) || (fabs(1.0-v)>eps) ); + bool do_vign = (fabs(1.0-v)>eps); + + for (int y=row_from; ywidth; x++) { + double x_d = (double) (x + cx) - w2 ; + + double r=0.0; + double s = d;//10000.0; + if (calc_r) + { + r=(sqrt(x_d*x_d + y_d*y_d)) / scale; + if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); + + // Convert only valid pixels + if (valid) { + // Extract integer and fractions of source screen coordinates + int xc = (int) (Dx); Dx -= (double)xc; + int yc = (int) (Dy); Dy -= (double)yc; + int ys = yc +1 - n2 - sy; // smallest y-index used for interpolation + int xs = xc +1 - n2 - sx; // smallest x-index used for interpolation + + double vignmul = 1.0; + if (do_vign) vignmul /= (v + mul * tanh (b*(maxRadius-s*r) / maxRadius)); + + if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2) // all interpolation pixels inside image + cubint (original, xs, ys, Dx, Dy, &(transformed->r[y][x]), &(transformed->g[y][x]), &(transformed->b[y][x]), vignmul); + else { // edge pixels + int y1 = (yc>0) ? yc : 0; + if (y1>miy) y1 = miy; + int y2 = (yc0) ? xc : 0; + if (x1>mix) x1 = mix; + int x2 = (xcr[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->r[y1][x2]*Dx*(1.0-Dy) + original->r[y2][x1]*(1.0-Dx)*Dy + original->r[y2][x2]*Dx*Dy); + int g = vignmul*(original->g[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->g[y1][x2]*Dx*(1.0-Dy) + original->g[y2][x1]*(1.0-Dx)*Dy + original->g[y2][x2]*Dx*Dy); + int b = vignmul*(original->b[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->b[y1][x2]*Dx*(1.0-Dy) + original->b[y2][x1]*(1.0-Dx)*Dy + original->b[y2][x2]*Dx*Dy); + transformed->r[y][x] = CLIP(r); + transformed->g[y][x] = CLIP(g); + transformed->b[y][x] = CLIP(b); + } + } + else { + // not valid (source pixel x,y not inside source image, etc.) + transformed->r[y][x] = 0; + transformed->g[y][x] = 0; + transformed->b[y][x] = 0; + } + } + } +} + +void ImProcFunctions::simpltransform_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { + + int oW = sizes.oW; + int oH = sizes.oH; + int cx = sizes.cx; + int cy = sizes.cy; + int sx = sizes.sx; + int sy = sizes.sy; + + double w2 = (double) oW / 2.0 - 0.5; + double h2 = (double) oH / 2.0 - 0.5; + + double cost = cos(params->rotate.degree * 3.14/180.0); + double sint = sin(params->rotate.degree * 3.14/180.0); + + double max_x = (double) (sx + original->width - 1); + double max_y = (double) (sy + original->height - 1); + double min_x = (double) sx; + double min_y = (double) sy; + + const int n2 = 2; + const int n = 2; + + int mix = original->width - 1; // maximum x-index src + int miy = original->height - 1;// maximum y-index src + int mix2 = mix +1 - n; + int miy2 = miy +1 - n; + + double scale = (oW>oH) ? (double)oW / 2.0 : (double)oH / 2.0 ; + double radius = sqrt( (double)( oW*oW + oH*oH ) ); + radius /= (oWdistortion.amount; + + double d = 1.0 - a; + + + // magnify image to keep size + double rotmagn = 1.0; + if (params->rotate.fill) { + double beta = atan((double)MIN(oH,oW)/MAX(oW,oH)); + rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); + } + // 1. check upper and lower border + double d1r = rotmagn - a*h2/scale - params->cacorrection.red; + double d2r = rotmagn - a*w2/scale - params->cacorrection.red; + double d3r = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.red; + double dr = MIN(d,MIN(d1r,MIN(d2r,d3r))); + double d1b = rotmagn - a*h2/scale - params->cacorrection.blue; + double d2b = rotmagn - a*w2/scale - params->cacorrection.blue; + double d3b = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.blue; + double db = MIN(d,MIN(d1b,MIN(d2b,d3b))); + double d1g = rotmagn - a*h2/scale; + double d2g = rotmagn - a*w2/scale; + double d3g = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; + double dg = MIN(d,MIN(d1g,MIN(d2g,d3g))); + + d = MIN(dg,MIN(dr,db)); + + // auxilary variables for vignetting + double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2 / scale; + + double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; + double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; + + double mul = (1.0-v) / tanh(b); + + // main cycle + double eps = 1e-10; + bool calc_r=( (fabs(a)>eps) || (fabs(1.0-v)>eps) ); + bool do_vign = (fabs(1.0-v)>eps); + + for (int y=row_from; ywidth; x++) { + double x_d = (double) (x + cx) - w2 ; + + double r=0.0; + double s = d;//10000.0; + if (calc_r) + { + r=(sqrt(x_d*x_d + y_d*y_d)) / scale; + if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); + + // Convert only valid pixels + if (valid) { + // Extract integer and fractions of source screen coordinates + int xc = (int) (Dx); Dx -= (double)xc; + int yc = (int) (Dy); Dy -= (double)yc; + int ys = yc +1 - n2 - sy; // smallest y-index used for interpolation + int xs = xc +1 - n2 - sx; // smallest x-index used for interpolation + + double vignmul = 1.0; + if (do_vign) vignmul /= (v + mul * tanh (b*(maxRadius-s*r) / maxRadius)); + + if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2 && yc < miy-1) { // all interpolation pixels inside image + + int r = vignmul*(original->r[yc][xc]*(1.0-Dx)*(1.0-Dy) + original->r[yc][xc+1]*Dx*(1.0-Dy) + original->r[yc+1][xc]*(1.0-Dx)*Dy + original->r[yc+1][xc+1]*Dx*Dy); + int g = vignmul*(original->g[yc][xc]*(1.0-Dx)*(1.0-Dy) + original->g[yc][xc+1]*Dx*(1.0-Dy) + original->g[yc+1][xc]*(1.0-Dx)*Dy + original->g[yc+1][xc+1]*Dx*Dy); + int b = vignmul*(original->b[yc][xc]*(1.0-Dx)*(1.0-Dy) + original->b[yc][xc+1]*Dx*(1.0-Dy) + original->b[yc+1][xc]*(1.0-Dx)*Dy + original->b[yc+1][xc+1]*Dx*Dy); + transformed->r[y][x] = CLIP(r); + transformed->g[y][x] = CLIP(g); + transformed->b[y][x] = CLIP(b); + } + else { // edge pixels + int y1 = (yc>0) ? yc : 0; + if (y1>miy) y1 = miy; + int y2 = (yc0) ? xc : 0; + if (x1>mix) x1 = mix; + int x2 = (xcr[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->r[y1][x2]*Dx*(1.0-Dy) + original->r[y2][x1]*(1.0-Dx)*Dy + original->r[y2][x2]*Dx*Dy); + int g = vignmul*(original->g[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->g[y1][x2]*Dx*(1.0-Dy) + original->g[y2][x1]*(1.0-Dx)*Dy + original->g[y2][x2]*Dx*Dy); + int b = vignmul*(original->b[y1][x1]*(1.0-Dx)*(1.0-Dy) + original->b[y1][x2]*Dx*(1.0-Dy) + original->b[y2][x1]*(1.0-Dx)*Dy + original->b[y2][x2]*Dx*Dy); + transformed->r[y][x] = CLIP(r); + transformed->g[y][x] = CLIP(g); + transformed->b[y][x] = CLIP(b); + } + } + else { + // not valid (source pixel x,y not inside source image, etc.) + transformed->r[y][x] = 0; + transformed->g[y][x] = 0; + transformed->b[y][x] = 0; + } + } + } +} + + +#include "cubintch.cc" +void ImProcFunctions::transform_sep_ (Image16* original, Image16* transformed, const ProcParams* params, STemp sizes, int row_from, int row_to) { + + int oW = sizes.oW; + int oH = sizes.oH; + int cx = sizes.cx; + int cy = sizes.cy; + int sx = sizes.sx; + int sy = sizes.sy; + + double w2 = (double) oW / 2.0 - 0.5; + double h2 = (double) oH / 2.0 - 0.5; + + double cost = cos(params->rotate.degree * 3.14/180.0); + double sint = sin(params->rotate.degree * 3.14/180.0); + + double max_x = (double) (sx + original->width - 1); + double max_y = (double) (sy + original->height - 1); + double min_x = (double) sx; + double min_y = (double) sy; + + const int n2 = 2; + const int n = 4; + + int mix = original->width - 1; // maximum x-index src + int miy = original->height - 1;// maximum y-index src + int mix2 = mix +1 - n; + int miy2 = miy +1 - n; + + double scale = (oW>oH) ? (double)oW / 2.0 : (double)oH / 2.0 ; + double radius = sqrt( (double)( oW*oW + oH*oH ) ); + radius /= (oWdistortion.amount; + double d = 1.0 - a; + + double cdist[3]; + cdist[0] = params->cacorrection.red; + cdist[1] = 0.0; + cdist[2] = params->cacorrection.blue; + + // magnify image to keep size + double rotmagn = 1.0; + if (params->rotate.fill) { + double beta = atan((double)MIN(oH,oW)/MAX(oW,oH)); + rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); + } + // 1. check upper and lower border + double d1r = rotmagn - a*h2/scale - params->cacorrection.red; + double d2r = rotmagn - a*w2/scale - params->cacorrection.red; + double d3r = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.red; + double dr = MIN(d,MIN(d1r,MIN(d2r,d3r))); + double d1b = rotmagn - a*h2/scale - params->cacorrection.blue; + double d2b = rotmagn - a*w2/scale - params->cacorrection.blue; + double d3b = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.blue; + double db = MIN(d,MIN(d1b,MIN(d2b,d3b))); + double d1g = rotmagn - a*h2/scale; + double d2g = rotmagn - a*w2/scale; + double d3g = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; + double dg = MIN(d,MIN(d1g,MIN(d2g,d3g))); + + d = MIN(dg,MIN(dr,db)); + + unsigned short** chorig[3]; + chorig[0] = original->r; + chorig[1] = original->g; + chorig[2] = original->b; + + unsigned short** chtrans[3]; + chtrans[0] = transformed->r; + chtrans[1] = transformed->g; + chtrans[2] = transformed->b; + + + // auxilary variables for vignetting + double maxRadius = sqrt( (double)( oW*oW + oH*oH ) ) / 2 / scale; + + double v = 1.0 - params->vignetting.amount * 3.0 / 400.0; + double b = 1.0 + params->vignetting.radius * 7.0 / 100.0; + + double mul = (1.0-v) / tanh(b); + + // main cycle + double eps = 1e-10; + for (int y=row_from; ywidth; x++) { + double x_d = (double) (x + cx) - w2 ; + + double r = (sqrt(x_d*x_d + y_d*y_d)) / scale; + double s = 10000.0; + if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); + + // Convert only valid pixels + if (valid) { + // Extract integer and fractions of source screen coordinates + int xc = (int) (Dx); Dx -= (double)xc; + int yc = (int) (Dy); Dy -= (double)yc; + int ys = yc +1 - n2 - sy; // smallest y-index used for interpolation + int xs = xc +1 - n2 - sx; // smallest x-index used for interpolation + + if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2) // all interpolation pixels inside image + cubintch (chorig[c], xs, ys, Dx, Dy, &(chtrans[c][y][x]), vignmul); + else {// edge pixels, linear interpolation + int y1 = (yc>0) ? yc : 0; + if (y1>miy) y1 = miy; + int y2 = (yc0) ? xc : 0; + if (x1>mix) x1 = mix; + int x2 = (xc &src, std::vector &red, std::vector &green, std::vector &blue) { + + bool clipresize = true; + bool clipped = false; + + red.clear (); + green.clear (); + blue.clear (); + bool needstransform = 0;// fabs(params->rotate.degree)>1e-15 || fabs(params->distortion.amount)>1e-15 || fabs(params->cacorrection.red)>1e-15 || fabs(params->cacorrection.blue)>1e-15; + if (!needstransform) { + if (clipresize) { + // Apply resizing + if (fabs(params->resize.scale-1.0)>=1e-7) { + for (int i=0; iresize.scale, src[i].y / params->resize.scale)); + green.push_back (Coord2D (src[i].x / params->resize.scale, src[i].y / params->resize.scale)); + blue.push_back (Coord2D (src[i].x / params->resize.scale, src[i].y / params->resize.scale)); + } + for (int i=0; iresize.scale; + double rH = H*params->resize.scale; + double w2 = (double) rW / 2.0 - 0.5; + double h2 = (double) rH / 2.0 - 0.5; + double cost = cos(params->rotate.degree * 3.14/180.0); + double sint = sin(params->rotate.degree * 3.14/180.0); + + double scale = (rW>rH) ? rW / 2.0 : rH / 2.0 ; + double radius = sqrt ((double)(rW*rW + rH*rH )); + radius /= (rWdistortion.amount; + double d = 1.0 - a; + + // magnify image to keep size + double rotmagn = 1.0; + if (params->rotate.fill) { + double beta = atan(MIN(rH,rW)/MAX(rW,rH)); + rotmagn = sin(beta) / sin(fabs(params->rotate.degree) * 3.14/180.0 + beta); + } + if (params->cacorrection.red==0 && params->cacorrection.blue==0) { + // 1. check upper and lower border + double d1 = rotmagn - a*h2/scale; + double d2 = rotmagn - a*w2/scale; + double d3 = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; + d = MIN(d,MIN(d1,MIN(d2,d3))); + + for (int i=0; icacorrection.red; + cdist[1] = 0.0; + cdist[2] = params->cacorrection.blue; + + // 1. check upper and lower border + double d1r = rotmagn - a*h2/scale - params->cacorrection.red; + double d2r = rotmagn - a*w2/scale - params->cacorrection.red; + double d3r = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.red; + double dr = MIN(d,MIN(d1r,MIN(d2r,d3r))); + double d1b = rotmagn - a*h2/scale - params->cacorrection.blue; + double d2b = rotmagn - a*w2/scale - params->cacorrection.blue; + double d3b = rotmagn - a*sqrt(h2*h2+w2*w2) / scale - params->cacorrection.blue; + double db = MIN(d,MIN(d1b,MIN(d2b,d3b))); + double d1g = rotmagn - a*h2/scale; + double d2g = rotmagn - a*w2/scale; + double d3g = rotmagn - a*sqrt(h2*h2+w2*w2) / scale; + double dg = MIN(d,MIN(d1g,MIN(d2g,d3g))); + + d = MIN(dg,MIN(dr,db)); + + for (int i=0; iresize.scale-1.0)>=1e-7) { + for (int i=0; iresize.scale; + red[i].y /= params->resize.scale; + green[i].x /= params->resize.scale; + green[i].y /= params->resize.scale; + blue[i].x /= params->resize.scale; + blue[i].y /= params->resize.scale; + } + } + for (int i=0; i corners (8); + corners[0].set (x1, y1); + corners[1].set (x1, y2); + corners[2].set (x2, y2); + corners[3].set (x2, y1); + corners[4].set ((x1+x2)/2, y1); + corners[5].set ((x1+x2)/2, y2); + corners[6].set (x1, (y1+y2)/2); + corners[7].set (x2, (y1+y2)/2); + + std::vector r, g, b; + + bool result = transCoord (params, W, H, corners, r, g, b); + + std::vector transCorners; + transCorners.insert (transCorners.end(), r.begin(), r.end()); + transCorners.insert (transCorners.end(), g.begin(), g.end()); + transCorners.insert (transCorners.end(), b.begin(), b.end()); + + double x1d = transCorners[0].x; + for (int i=1; ix2d) + x2d = transCorners[i].x; + int x2v = (int)ceil(x2d); + + double y2d = transCorners[0].y; + for (int i=1; iy2d) + y2d = transCorners[i].y; + int y2v = (int)ceil(y2d); + + xv = x1v; + yv = y1v; + wv = x2v - x1v + 1; + hv = y2v - y1v + 1; + + return result; +} + +void ImProcFunctions::transform (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int sx, int sy, int oW, int oH) { + + STemp sizes; + sizes.cx = 0;//cx; + sizes.cy = 0;//cy; + sizes.oW = oW; + sizes.oH = oH; + sizes.sx = 0;//sx; + sizes.sy = 0;//sy; + + if (params->cacorrection.red==0 && params->cacorrection.blue==0) { + if (settings->dualThreadEnabled) { + Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + thread1->join (); + thread2->join (); + } + else + transform_ (original, transformed, params, sizes, 0, transformed->height); + } + else { + if (settings->dualThreadEnabled) { + Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_sep_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::transform_sep_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + thread1->join (); + thread2->join (); + } + else + transform_sep_ (original, transformed, params, sizes, 0, transformed->height); + } +} + +void ImProcFunctions::simpltransform (Image16* original, Image16* transformed, const ProcParams* params, int cx, int cy, int sx, int sy, int oW, int oH) { + + STemp sizes; + sizes.cx = 0;//cx; + sizes.cy = 0;//cy; + sizes.oW = oW; + sizes.oH = oH; + sizes.sx = 0;//sx; + sizes.sy = 0;//sy; + + if (settings->dualThreadEnabled) { + Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::simpltransform_), original, transformed, params, sizes, 0, transformed->height/2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::mem_fun(*this, &ImProcFunctions::simpltransform_), original, transformed, params, sizes, transformed->height/2, transformed->height), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); + thread1->join (); + thread2->join (); + } + else + simpltransform_ (original, transformed, params, sizes, 0, transformed->height); +} +/*void ImProcFunctions::transform (Image16* original, Image16* transformed, const ProcParams* params, int ox, int oy) { + + if (!transformed) + return; + + int oW = W, oH = H, tW = W, tH = H; + + double w2 = (double) tW / 2.0 - 0.5; + double h2 = (double) tH / 2.0 - 0.5; + double sw2 = (double) oW / 2.0 - 0.5; + double sh2 = (double) oH / 2.0 - 0.5; + + double cost = cos(params->rotate_fine * 3.14/180.0); + double sint = sin(params->rotate_fine * 3.14/180.0); + + double max_x = (double) oW; + double max_y = (double) oH; + double min_x = 0.0; + double min_y = 0.0; + + const int n2 = 2; + const int n = 4; + + int mix = oW - 1; // maximum x-index src + int miy = oH - 1;// maximum y-index src + int mix2 = mix +1 - n; + int miy2 = miy +1 - n; + + double scale = (tW>tH) ? (double)tW / 2.0 : (double)tH / 2.0 ; + double radius = sqrt( (double)( tW*tW + tH*tH ) ); + radius /= (tWlens_distortion; + + for (int y=0; yheight; y++) { + double y_d = (double) y + oy - h2 ; + + for (int x=0; xwidth; x++) { + double x_d = (double) x + ox - w2 ; + + double r = (sqrt(x_d*x_d + y_d*y_d)) / scale; + double s = 10000.0; + if (r= max_x) || (Dy >= max_y) || (Dx < min_x) || (Dy < min_y)); + + // Convert only valid pixels + if (valid) { + // Extract integer and fractions of source screen coordinates + int xc = (int) floor (Dx) ; Dx -= (double)xc; + int yc = (int) floor (Dy) ; Dy -= (double)yc; + int ys = yc +1 - n2 ; // smallest y-index used for interpolation + int xs = xc +1 - n2 ; // smallest x-index used for interpolation + + unsigned short sr[2][2], sg[2][2], sb[2][2]; + + if (ys >= 0 && ys <= miy2 && xs >= 0 && xs <= mix2) // all interpolation pixels inside image + cubint (original, xs, ys, Dx, Dy, &(transformed->r[y][x]), &(transformed->g[y][x]), &(transformed->b[y][x])); + else { // edge pixels + transformed->r[y][x] = 0; + transformed->g[y][x] = 0; + transformed->b[y][x] = 0; + } + } + else { + // not valid (source pixel x,y not inside source image, etc.) + transformed->r[y][x] = 0; + transformed->g[y][x] = 0; + transformed->b[y][x] = 0; + } + } + } +}*/ + +} + diff --git a/rtengine/labimage.cc b/rtengine/labimage.cc new file mode 100644 index 000000000..34a4cec73 --- /dev/null +++ b/rtengine/labimage.cc @@ -0,0 +1,42 @@ +#include +namespace rtengine { + +LabImage::LabImage (int w, int h) : W(w), H(h), fromImage(false) { + + L = new unsigned short*[H]; + for (int i=0; iwidth; + H = im->height; + L = im->r; + a = (short**) im->g; + b = (short**) im->b; + fromImage = true; +} + +LabImage::~LabImage () { + + if (!fromImage) { + for (int i=0; i. */ -#include +#ifndef _LABIMAGE_H_ +#define _LABIMAGE_H_ -void bilateral_unsigned (unsigned short** src, unsigned short** dst, unsigned short** buffer, Dim dim, double sigma, double sens) { +#include - bilateral (src, dst, buffer, dim, sigma, sens); +namespace rtengine { + +class LabImage { + private: + bool fromImage; + + public: + int W, H; + unsigned short** L; + short** a; + short** b; + + LabImage (int w, int h); + LabImage (Image16* im); + ~LabImage (); +}; } - -void bilateral_signed (short** src, short** dst, short** buffer, Dim dim, double sigma, double sens) { - - bilateral (src, dst, buffer, dim, sigma, sens); -} - -void bilateral_box_unsigned (unsigned short** src, unsigned short** dst, int W, int H, int sigmar, double sigmas, bilateralparams row) { - - bilateral (src, dst, W, H, sigmar, sigmas, row.row_from, row.row_to); -} - +#endif diff --git a/rtengine/rtthumbnail.cc b/rtengine/rtthumbnail.cc index 7dce60e07..7ec958787 100644 --- a/rtengine/rtthumbnail.cc +++ b/rtengine/rtthumbnail.cc @@ -206,7 +206,7 @@ IImage8* Thumbnail::processImage (const procparams::ProcParams& params, int rhei gmi = 1024.0 * gm / mul_lum; bmi = 1024.0 * bm / mul_lum; } - // resize to requested with and perform coarse transformation + // resize to requested width and perform coarse transformation int rwidth; if (params.coarse.rotate==90 || params.coarse.rotate==270) { rwidth = rheight; @@ -269,7 +269,9 @@ IImage8* Thumbnail::processImage (const procparams::ProcParams& params, int rhei int fw = baseImg->width; int fh = baseImg->height; - ImProcFunctions ipf; + ImProcFunctions ipf (¶ms, false); + ipf.setScale (sqrt(double(fw*fw+fh*fh))/sqrt(double(thumbImg->width*thumbImg->width+thumbImg->height*thumbImg->height))*scale); + unsigned int* hist16 = new unsigned int [65536]; ipf.firstAnalysis (baseImg, ¶ms, hist16, isRaw ? 2.2 : 0.0); @@ -299,7 +301,7 @@ IImage8* Thumbnail::processImage (const procparams::ProcParams& params, int rhei for (int i=0; iupdate (baseImg, buffer, shradius, ipf.lumimul, params.sh.hq); @@ -321,7 +323,7 @@ IImage8* Thumbnail::processImage (const procparams::ProcParams& params, int rhei CurveFactory::complexCurve (br, bl/65535.0, params.toneCurve.hlcompr, params.toneCurve.shcompr, params.toneCurve.brightness, params.toneCurve.contrast, logDefGain, isRaw ? 2.2 : 0, true, params.toneCurve.curve, hist16, curve, NULL, 16); LabImage* labView = new LabImage (baseImg); - ipf.rgbProc (baseImg, labView, ¶ms, curve, shmap); + ipf.rgbProc (baseImg, labView, curve, shmap); if (shmap) delete shmap; @@ -341,12 +343,11 @@ IImage8* Thumbnail::processImage (const procparams::ProcParams& params, int rhei delete [] hist16; // color processing - ipf.colorCurve (labView, labView, ¶ms); + ipf.colorCurve (labView, labView); // obtain final image Image8* readyImg = new Image8 (fw, fh); ipf.lab2rgb (labView, readyImg); - ipf.release (); delete baseImg; // calculate scale @@ -453,8 +454,7 @@ void Thumbnail::getSpotWB (const procparams::ProcParams& params, int xp, int yp, fw = thumbImg->height; fh = thumbImg->width; } - ImProcFunctions ipf; - ipf.transCoord (¶ms, fw, fh, points, red, green, blue); + ImProcFunctions::transCoord (¶ms, fw, fh, points, red, green, blue); int tr = TR_NONE; if (params.coarse.rotate==90) tr |= TR_R90; if (params.coarse.rotate==180) tr |= TR_R180; diff --git a/rtengine/shmap.cc b/rtengine/shmap.cc index b9531ac98..308bf0823 100644 --- a/rtengine/shmap.cc +++ b/rtengine/shmap.cc @@ -29,7 +29,7 @@ namespace rtengine { extern const Settings* settings; -SHMap::SHMap (int w, int h) : W(w), H(h) { +SHMap::SHMap (int w, int h, bool multiThread) : W(w), H(h), multiThread(multiThread) { map = new unsigned short*[H]; for (int i=0; i* buffer = new AlignedBuffer (MAX(W,H)*omp_get_max_threads()); - AlignedBuffer* buffer1 = new AlignedBuffer (MAX(W,H)*5); - AlignedBuffer* buffer2 = new AlignedBuffer (MAX(W,H)*5); + gaussHorizontal (map, map, buffer, W, H, radius, multiThread); + gaussVertical (map, map, buffer, W, H, radius, multiThread); - // blur - if (settings->dualThreadEnabled) { - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_unsigned), map, map, buffer1, W, 0, H/2, radius), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussHorizontal_unsigned), map, map, buffer2, W, H/2, H, radius), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_unsigned), map, map, buffer1, H, 0, W/2, radius), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(gaussVertical_unsigned), map, map, buffer2, H, W/2, W, radius), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else { - gaussHorizontal_unsigned (map, map, buffer1, W, 0, H, radius); - gaussVertical_unsigned (map, map, buffer1, H, 0, W, radius); - } - - delete buffer1; - delete buffer2; + delete buffer; } else { - if (settings->dualThreadEnabled) { - bilateralparams r1, r2; - r1.row_from = 0; - r1.row_to = H/2; - r2.row_from = H/2; - r2.row_to = H; - Glib::Thread *thread1 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_box_unsigned), map, buffer, W, H, 8000, radius, r1), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - Glib::Thread *thread2 = Glib::Thread::create(sigc::bind(sigc::ptr_fun(bilateral_box_unsigned), map, buffer, W, H, 8000, radius, r2), 0, true, true, Glib::THREAD_PRIORITY_NORMAL); - thread1->join (); - thread2->join (); - } - else { - bilateralparams r1; - r1.row_from = 0; - r1.row_to = H; - bilateral_box_unsigned (map, buffer, W, H, 8000, radius, r1); - } + #pragma omp parallel if (multiThread) + { + int tid = omp_get_thread_num(); + int nthreads = omp_get_num_threads(); + int blk = H/nthreads; + + if (tid (map, buffer, W, H, 8000, radius, tid*blk, (tid+1)*blk); + else + bilateral (map, buffer, W, H, 8000, radius, tid*blk, H); + } + // anti-alias filtering the result for (int i=0; i0 && j>0 && igetFullSize (fw, fh, tr); - ImProcFunctions ipf; + ImProcFunctions ipf (¶ms, true); Image16* baseImg; PreviewProps pp (0, 0, fw, fh, 1); @@ -83,7 +83,7 @@ IImage16* processImage (ProcessingJob* pjob, int& errorCode, ProgressListener* p fw *= params.resize.scale; fh *= params.resize.scale; baseImg = new Image16 (fw, fh); - ipf.resize (oorig, baseImg, params.resize); + ipf.resize (oorig, baseImg); delete oorig; } if (pl) @@ -117,7 +117,7 @@ IImage16* processImage (ProcessingJob* pjob, int& errorCode, ProgressListener* p SHMap* shmap = NULL; if (params.sh.enabled) { - shmap = new SHMap (fw, fh); + shmap = new SHMap (fw, fh, true); double radius = sqrt (double(fw*fw+fh*fh)) / 2.0; double shradius = radius / 1800.0 * params.sh.radius; shmap->update (baseImg, (unsigned short**)buffer, shradius, ipf.lumimul, params.sh.hq); @@ -138,7 +138,7 @@ IImage16* processImage (ProcessingJob* pjob, int& errorCode, ProgressListener* p CurveFactory::complexCurve (br, bl/65535.0, params.toneCurve.hlcompr, params.toneCurve.shcompr, params.toneCurve.brightness, params.toneCurve.contrast, imgsrc->getDefGain(), imgsrc->getGamma(), true, params.toneCurve.curve, hist16, curve, NULL); LabImage* labView = new LabImage (baseImg); - ipf.rgbProc (baseImg, labView, ¶ms, curve, shmap); + ipf.rgbProc (baseImg, labView, curve, shmap); if (shmap) delete shmap; @@ -155,15 +155,15 @@ IImage16* processImage (ProcessingJob* pjob, int& errorCode, ProgressListener* p // luminance processing CurveFactory::complexCurve (0.0, 0.0, 0.0, 0.0, params.lumaCurve.brightness, params.lumaCurve.contrast, 0.0, 0.0, false, params.lumaCurve.curve, hist16, curve, NULL); ipf.luminanceCurve (labView, labView, curve, 0, fh); - ipf.lumadenoise (labView, ¶ms, 1, buffer); - ipf.sharpening (labView, ¶ms, 1, (unsigned short**)buffer); + ipf.lumadenoise (labView, buffer); + ipf.sharpening (labView, (unsigned short**)buffer); delete [] curve; delete [] hist16; // color processing - ipf.colorCurve (labView, labView, ¶ms); - ipf.colordenoise (labView, ¶ms, 1, buffer); + ipf.colorCurve (labView, labView); + ipf.colordenoise (labView, buffer); for (int i=0; isetProgress (1.0);