parallelization and bug fixes

modules/xphoto/CMakeLists.txt

@@ -1,3 +1,3 @@
 set(the_description "Addon to basic photo module")
 ocv_warnings_disable(CMAKE_CXX_FLAGS -Wundef)
-ocv_define_module(xphoto opencv_core opencv_imgproc OPTIONAL opencv_highgui)
+ocv_define_module(xphoto opencv_core opencv_imgproc OPTIONAL opencv_photo opencv_highgui)

modules/xphoto/src/dct_image_denoising.cpp

@@ -54,28 +54,93 @@
 
 namespace cv
 {
-    void grayDctDenoising(const Mat_<float> &src, Mat_<float> &dst, const double sigma, const int psize)
+    struct grayDctDenoisingInvoker : public ParallelLoopBody
     {
-        CV_Assert( src.channels() == 1 );
+    public:
+        grayDctDenoisingInvoker(const Mat_<float> &src, std::vector < Mat_<float> > &patches, const double sigma, const int psize);
+        ~grayDctDenoisingInvoker();
 
-        Mat_<float> res( src.size(), 0.0f ),
-                    num( src.size(), 0.0f );
+        void operator() (const Range &range) const;
+
+    private:
+        const Mat_<float> &src;
+        std::vector < Mat_<float> > &patches; // image decomposition into sliding patches
 
-        double threshold = 2.0*log(psize)*sigma;
+        const int psize; // size of block to compute dct
+        const double sigma; // expected noise standard deviation
+        const double thresh; // thresholding estimate
+    };
 
-        for (int i = 0; i <= src.rows - psize; ++i)
-            for (int j = 0; j <= src.cols - psize; ++j)
+    grayDctDenoisingInvoker::grayDctDenoisingInvoker(const Mat_<float> &src, std::vector < Mat_<float> > &patches,
+                                                     const double sigma, const int psize)
+        : src(src), patches(patches), sigma(sigma), thresh(3*sigma), psize(psize) {}
+    grayDctDenoisingInvoker::~grayDctDenoisingInvoker(){}
+
+    void grayDctDenoisingInvoker::operator() (const Range &range) const
+    {
+        for (int i = range.start; i <= range.end - 1; ++i)
         {
-                Mat_<float> patch = src( Rect(j, i, psize, psize) ).clone();
+            int y = i / (src.cols - psize);
+            int x = i % (src.cols - psize);
+
+            Rect patchNum( x, y, psize, psize );
+
+            Mat_<float> patch(psize, psize);
+            src(patchNum).copyTo( patch );
 
             dct(patch, patch);
-                float * ptr = (float *) patch.data;
+            float *data = (float *) patch.data;
             for (int k = 0; k < psize*psize; ++k)
-                    if (fabs(ptr[k]) < threshold)
-                        ptr[k] = 0.0f;
-                idct(patch, patch);
+                data[k] *= fabs(data[k]) > thresh;
+            idct(patch, patches[i]);
+        }
+    }
+
+    void grayDctDenoising(const Mat_<float> &src, Mat_<float> &dst, const double sigma, const int psize)
+    {
+        CV_Assert( src.channels() == 1 );
+
+        //Mat_<float> res( src.size(), 0.0f ),
+        //            num( src.size(), 0.0f );
+        //
+        //double threshold = 3*sigma;
+        //
+        //for (int i = 0; i <= src.rows - psize; ++i)
+        //    for (int j = 0; j <= src.cols - psize; ++j)
+        //    {
+        //        Mat_<float> patch = src( Rect(j, i, psize, psize) ).clone();
+        //
+        //        dct(patch, patch);
+        //        float * ptr = (float *) patch.data;
+        //        for (int k = 0; k < psize*psize; ++k)
+        //            if (fabs(ptr[k]) < threshold)
+        //                ptr[k] = 0.0f;
+        //        idct(patch, patch);
+        //
+        //        res( Rect(j, i, psize, psize) ) += patch;
+        //        num( Rect(j, i, psize, psize) ) += Mat_<float>::ones(psize, psize);
+        //    }
+        //res /= num;
+        //
+        //res.convertTo( dst, src.type() );
+
+        int npixels = (src.rows - psize)*(src.cols - psize);
+
+        std::vector < Mat_<float> > patches;
+        for (int i = 0; i < npixels; ++i)
+            patches.push_back( Mat_<float>(psize, psize) );
+        parallel_for_( cv::Range(0, npixels),
+            grayDctDenoisingInvoker(src, patches, sigma, psize) );
+
+        Mat_<float> res( src.size(), 0.0f ),
+                    num( src.size(), 0.0f );
+
+        for (int k = 0; k < npixels; ++k)
+        {
+            int i = k / (src.cols - psize);
+            int j = k % (src.cols - psize);
 
-                res( Rect(j, i, psize, psize) ) += patch;
+            res( Rect(j, i, psize, psize) ) += patches[k];
             num( Rect(j, i, psize, psize) ) += Mat_<float>::ones(psize, psize);
         }
         res /= num;
@@ -96,9 +161,9 @@ namespace cv
         for (Mat_<Vec3f>::const_iterator it = src.begin(); it != src.end(); ++it, ++outIt)
         {
             Vec3f rgb = *it;
-            *outIt = Vec3f(M[0]*rgb[0] + M[3]*rgb[1] + M[6]*rgb[2],
-                           M[1]*rgb[0] + M[4]*rgb[1] + M[7]*rgb[2],
-                           M[2]*rgb[0] + M[5]*rgb[1] + M[8]*rgb[2]);
+            *outIt = Vec3f(M[0]*rgb[0] + M[1]*rgb[1] + M[2]*rgb[2],
+                           M[3]*rgb[0] + M[4]*rgb[1] + M[5]*rgb[2],
+                           M[6]*rgb[0] + M[7]*rgb[1] + M[8]*rgb[2]);
         }
 
         /*************************************/
@@ -106,7 +171,7 @@ namespace cv
         split(dst, mv);
 
         for (int i = 0; i < mv.size(); ++i)
-            grayDctDenoising(mv[0], mv[0], sigma, psize);
+            grayDctDenoising(mv[i], mv[i], sigma, psize);
 
         merge(mv, dst);
         /*************************************/
@@ -114,9 +179,9 @@ namespace cv
         for (Mat_<Vec3f>::iterator it = dst.begin(); it != dst.end(); ++it)
         {
             Vec3f rgb = *it;
-            *it = Vec3f(M[0]*rgb[0] + M[1]*rgb[1] + M[2]*rgb[2],
-                        M[3]*rgb[0] + M[4]*rgb[1] + M[5]*rgb[2],
-                        M[6]*rgb[0] + M[7]*rgb[1] + M[8]*rgb[2]);
+            *it = Vec3f(M[0]*rgb[0] + M[3]*rgb[1] + M[6]*rgb[2],
+                        M[1]*rgb[0] + M[4]*rgb[1] + M[7]*rgb[2],
+                        M[2]*rgb[0] + M[5]*rgb[1] + M[8]*rgb[2]);
         }
     }
 

modules/xphoto/test/dct_image_denoising.cpp

@@ -1,15 +1,18 @@
 #include "test_precomp.hpp"
 
+#define NO_COMPARISON
+
 namespace cvtest
 {
     TEST(xphoto_dctimagedenoising, regression)
     {
         cv::String dir = cvtest::TS::ptr()->get_data_path() + "dct_image_denoising/";
-        int nTests = 12;
+        int nTests = 1;
+
+        float psnrThreshold[] = {0.5};
 
-		int psize = 3.0;
-		float psnrThreshold = 40.0;
-		float sigma = 15.0;
+        int psize[] = {8};
+        double sigma[] = {9.0};
 
         for (int i = 0; i < nTests; ++i)
         {
@@ -19,16 +22,51 @@ namespace cvtest
             cv::String previousResultName = dir + cv::format( "results/%02d.png", i + 1 );
             cv::Mat previousResult = cv::imread( previousResultName, 1 );
 
+            cv::Mat sqrError = ( src - previousResult ).mul( src - previousResult );
+            cv::Scalar mse = cv::sum(sqrError) / cv::Scalar::all( sqrError.total()*sqrError.channels() );
+            double psnr = 10*log10(3*255*255/(mse[0] + mse[1] + mse[2]));
+
+
             cv::Mat currentResult, fastNlMeansResult;
-            cv::dctDenoising(src, currentResult, sigma, psize);
-            //cv::fastNlMeansDenoising(src, fastNlMeansResult);
 
-            cv::Mat sqrError = ( currentResult - previousResult )
+#ifndef NO_COMPARISON
+            double currentTime = clock() / double(CLOCKS_PER_SEC);
+#endif
+            cv::dctDenoising(src, currentResult, sigma[i], psize[i]);
+#ifndef NO_COMPARISON
+            currentTime = clock() / double(CLOCKS_PER_SEC) - currentTime;
+            std::cout << "---- dct denoising time = " << currentTime << " (sec) ----" << std::endl;
+#endif
+
+            cv::Mat sqrError1 = ( currentResult - previousResult )
                 .mul( currentResult - previousResult );
-            cv::Scalar mse = cv::sum(sqrError) / cv::Scalar::all( sqrError.total()*sqrError.channels() );
-			double psnr = 10*log10(3*255*255/(mse[0] + mse[1] + mse[2]));
+            cv::Scalar mse1 = cv::sum(sqrError1) / cv::Scalar::all( sqrError1.total()*sqrError1.channels() );
+            double psnr1 = 10*log10(3*255*255/(mse1[0] + mse1[1] + mse1[2])) - psnr;
+#ifndef NO_COMPARISON
+            std::cout << "---- dct PSNR rate = " << psnr1 << " ----" << std::endl;
+#endif
+
+#ifndef NO_COMPARISON
+            double fastNlMeansTime = clock() / double(CLOCKS_PER_SEC);
+
+            if ( src.channels() == 3 )
+                cv::fastNlMeansDenoisingColored(src, fastNlMeansResult);
+            else if ( src.channels() == 1 )
+                cv::fastNlMeansDenoising(src, fastNlMeansResult);
+
+            fastNlMeansTime = clock() / double(CLOCKS_PER_SEC) - fastNlMeansTime;
+#ifdef NO_COMPARISON
+            std::cout << "---- nonlocal means denoising time = " << fastNlMeansTime << " (sec) ----" << std::endl;
+#endif
+
+            cv::Mat sqrError2 = ( fastNlMeansResult - previousResult )
+                .mul( fastNlMeansResult - previousResult );
+            cv::Scalar mse2 = cv::sum(sqrError2) / cv::Scalar::all( sqrError2.total()*sqrError2.channels() );
+            double psnr2 = 10*log10(3*255*255/(mse2[0] + mse2[1] + mse2[2])) - psnr;
+            std::cout << "---- nonlocal means PSNR rate = " << psnr2 << " ----" << std::endl;
+#endif
 
-            EXPECT_GE( psnr, psnrThreshold );
+            EXPECT_GE( psnr1, psnrThreshold[i] );
         }
     }
 }

modules/xphoto/test/test_precomp.hpp

@@ -13,8 +13,11 @@
 #include "opencv2/imgproc.hpp"
 #include "opencv2/imgproc/types_c.h"
 #include "opencv2/highgui.hpp"
+#include "opencv2/photo.hpp"
 #include "opencv2/xphoto.hpp"
 #include "opencv2/ts.hpp"
+
+#include <ctime>
 #include <iostream>
 
 #endif