Exposure fusion. Code, tests.

modules/photo/include/opencv2/photo.hpp

@@ -96,8 +96,10 @@ CV_EXPORTS_W void fastNlMeansDenoisingColoredMulti( InputArrayOfArrays srcImgs,
 
 CV_EXPORTS_W void makeHDR(InputArrayOfArrays srcImgs, const std::vector<float>& exp_times, OutputArray dst);
 
-CV_EXPORTS_W void tonemap(InputArray src, OutputArray dst, tonemap_algorithms algorithm, std::vector<float>& params = std::vector<float>());
+CV_EXPORTS_W void tonemap(InputArray src, OutputArray dst, tonemap_algorithms algorithm, 
+	                      const std::vector<float>& params = std::vector<float>());
 
+CV_EXPORTS_W void exposureFusion(InputArrayOfArrays srcImgs, OutputArray dst, float wc = 1, float ws = 1, float we = 0);
 } // cv
 
 #endif

modules/photo/src/hdr_fusion.cpp

@@ -64,14 +64,12 @@ static void generateResponce(float responce[])
     responce[0] = responce[1];
 }
 
-void makeHDR(InputArrayOfArrays _images, const std::vector<float>& _exp_times, OutputArray _dst)
+static void checkImages(std::vector<Mat>& images, bool hdr, const std::vector<float>& _exp_times = std::vector<float>())
 {
-	std::vector<Mat> images;
-    _images.getMatVector(images);
 	if(images.empty()) {
 		CV_Error(Error::StsBadArg, "Need at least one image");
 	}
-	if(images.size() != _exp_times.size()) {
+	if(hdr && images.size() != _exp_times.size()) {
 		CV_Error(Error::StsBadArg, "Number of images and number of exposure times must be equal.");
 	}
 	int width = images[0].cols;
@@ -85,8 +83,16 @@ void makeHDR(InputArrayOfArrays _images, const std::vector<float>& _exp_times, O
 			CV_Error(Error::StsBadArg, "Images must have CV_8UC3 type.");
 		}
 	}
+}
+
+void makeHDR(InputArrayOfArrays _images, const std::vector<float>& _exp_times, OutputArray _dst)
+{
+	std::vector<Mat> images;
+    _images.getMatVector(images);
+	checkImages(images, true, _exp_times);
 	_dst.create(images[0].size(), CV_32FC3);
 	Mat result = _dst.getMat();
+
 	std::vector<float> exp_times(_exp_times.size());
 	for(size_t i = 0; i < exp_times.size(); i++) {
 		exp_times[i] = log(_exp_times[i]);
@@ -122,4 +128,88 @@ void makeHDR(InputArrayOfArrays _images, const std::vector<float>& _exp_times, O
 	result = result / max;
 }
 
+void exposureFusion(InputArrayOfArrays _images, OutputArray _dst, float wc, float ws, float we)
+{
+	std::vector<Mat> images;
+    _images.getMatVector(images);
+	checkImages(images, false);
+
+	std::vector<Mat> weights(images.size());
+	Mat weight_sum = Mat::zeros(images[0].size(), CV_32FC1);
+	for(size_t im = 0; im < images.size(); im++) {
+		Mat img, gray, contrast, saturation, wellexp;
+		std::vector<Mat> channels(3);
+
+		images[im].convertTo(img, CV_32FC3, 1.0/255.0);
+		cvtColor(img, gray, COLOR_RGB2GRAY);
+		split(img, channels);
+
+		Laplacian(gray, contrast, CV_32F);
+		contrast = abs(contrast);
+
+		Mat mean = (channels[0] + channels[1] + channels[2]) / 3.0f;
+		saturation = Mat::zeros(channels[0].size(), CV_32FC1);
+		for(int i = 0; i < 3;  i++) {
+			Mat deviation = channels[i] - mean;
+			pow(deviation, 2.0, deviation);
+			saturation += deviation;
+		}
+		sqrt(saturation, saturation);
+
+		wellexp = Mat::ones(gray.size(), CV_32FC1);
+		for(int i = 0; i < 3; i++) {
+			Mat exp = channels[i] - 0.5f;
+			pow(exp, 2, exp);
+			exp = -exp / 0.08;
+			wellexp = wellexp.mul(exp);
+		}
+
+		pow(contrast, wc, contrast);
+		pow(saturation, ws, saturation);
+		pow(wellexp, we, wellexp);
+
+		weights[im] = contrast;
+		weights[im] = weights[im].mul(saturation);
+		weights[im] = weights[im].mul(wellexp);
+		weight_sum += weights[im];
+	}
+	int maxlevel = (int)(log((double)max(images[0].rows, images[0].cols)) / log(2.0)) - 1;
+	std::vector<Mat> res_pyr(maxlevel + 1);
+
+	for(size_t im = 0; im < images.size(); im++) {
+		weights[im] /= weight_sum;
+		Mat img;
+		images[im].convertTo(img, CV_32FC3, 1/255.0);
+		std::vector<Mat> img_pyr, weight_pyr;
+		buildPyramid(img, img_pyr, maxlevel);
+		buildPyramid(weights[im], weight_pyr, maxlevel);
+		for(int lvl = 0; lvl < maxlevel; lvl++) {
+			Mat up;
+			pyrUp(img_pyr[lvl + 1], up, img_pyr[lvl].size());
+			img_pyr[lvl] -= up;
+		}
+		for(int lvl = 0; lvl <= maxlevel; lvl++) {
+			std::vector<Mat> channels(3);
+			split(img_pyr[lvl], channels);
+			for(int i = 0; i < 3; i++) {
+				channels[i] = channels[i].mul(weight_pyr[lvl]);
+			}
+			merge(channels, img_pyr[lvl]);
+			if(res_pyr[lvl].empty()) {
+				res_pyr[lvl] = img_pyr[lvl];
+			} else {
+				res_pyr[lvl] += img_pyr[lvl];
+			}
+		}
+	}
+	for(int lvl = maxlevel; lvl > 0; lvl--) {
+		Mat up;
+		pyrUp(res_pyr[lvl], up, res_pyr[lvl - 1].size());
+	    res_pyr[lvl - 1] += up;
+	}
+	_dst.create(images[0].size(), CV_32FC3);
+	Mat result = _dst.getMat();
+	res_pyr[0].copyTo(result);
+}
+
 };
\ No newline at end of file

modules/photo/src/tonemap.cpp

@@ -45,146 +45,147 @@
 
 namespace cv
 {
-	static float getParam(std::vector<float>& params, size_t i, float defval) 
-	{
-		if(params.size() > i) {
-			return params[i];
-		} else {
-			return defval;
-		}
-		
-	}
-	static void DragoMap(Mat& src_img, Mat &dst_img, std::vector<float>& params)
-	{
-		float bias_value = getParam(params, 1, 0.85f);
-		Mat gray_img;
-		cvtColor(src_img, gray_img, COLOR_RGB2GRAY);
-		Mat log_img;
-		log(gray_img, log_img);
-		float mean = exp((float)sum(log_img)[0] / log_img.total());
-		gray_img /= mean;
-		log_img.release();
-
-		double max;
-		minMaxLoc(gray_img, NULL, &max);
-
-		Mat map;
-		log(gray_img + 1.0f, map);
-		Mat div;
-		pow(gray_img / (float)max, log(bias_value) / log(0.5f), div);
-		log(2.0f + 8.0f * div, div);
-		map = map.mul(1.0f / div);
-		map = map.mul(1.0f / gray_img);
-		div.release();
-		gray_img.release();
-
-		std::vector<Mat> channels(3);
-		split(src_img, channels);
-		for(int i = 0; i < 3; i++) {
-			channels[i] = channels[i].mul(map);
-		}
-		map.release();
-		merge(channels, dst_img);
-	}
-
-	static void ReinhardDevlinMap(Mat& src_img, Mat &dst_img, std::vector<float>& params)
-	{
-		float intensity   = getParam(params, 1, 0.0f);
-		float color_adapt = getParam(params, 2, 0.0f);
-		float light_adapt = getParam(params, 3, 1.0f);
-
-		Mat gray_img;
-		cvtColor(src_img, gray_img, COLOR_RGB2GRAY);
-		Mat log_img;
-		log(gray_img, log_img);
-
-		float log_mean = (float)sum(log_img)[0] / log_img.total();
-		double log_min, log_max;
-		minMaxLoc(log_img, &log_min, &log_max);
-		log_img.release();
-
-		double key = (float)((log_max - log_mean) / (log_max - log_min));
-		float map_key = 0.3f + 0.7f * pow((float)key, 1.4f);
-		intensity = exp(-intensity);
-		Scalar chan_mean = mean(src_img);
-		float gray_mean = (float)mean(gray_img)[0];
-
-		std::vector<Mat> channels(3);
-		split(src_img, channels);
-
-		for(int i = 0; i < 3; i++) {
-			float global = color_adapt * (float)chan_mean[i] + (1.0f - color_adapt) * gray_mean;
-			Mat adapt = color_adapt * channels[i] + (1.0f - color_adapt) * gray_img;
-			adapt = light_adapt * adapt + (1.0f - light_adapt) * global;
-			pow(intensity * adapt, map_key, adapt);
-			channels[i] = channels[i].mul(1.0f / (adapt + channels[i]));		
-		}
-		gray_img.release();
-		merge(channels, dst_img);
-	}
-
-	static void DurandMap(Mat& src_img, Mat& dst_img, std::vector<float>& params)
-	{
-		float contrast   = getParam(params, 1, 4.0f);
-		float sigma_color = getParam(params, 2, 2.0f);
-		float sigma_space = getParam(params, 3, 2.0f);
-
-		Mat gray_img;
-		cvtColor(src_img, gray_img, COLOR_RGB2GRAY);
-		Mat log_img;
-		log(gray_img, log_img);
-		Mat map_img;
-		bilateralFilter(log_img, map_img, -1, sigma_color, sigma_space);
-		
-		double min, max;
-		minMaxLoc(map_img, &min, &max);
-		float scale = contrast / (float)(max - min);
-
-		exp(map_img * (scale - 1.0f) + log_img, map_img);
-		log_img.release();
-		map_img = map_img.mul(1.0f / gray_img);
-		gray_img.release();
-
-		std::vector<Mat> channels(3);
-		split(src_img, channels);
-		for(int i = 0; i < 3; i++) {
-			channels[i] = channels[i].mul(map_img);
-		}
-		merge(channels, dst_img);
-	}
-
-	void tonemap(InputArray _src, OutputArray _dst, tonemap_algorithms algorithm, 
-				 std::vector<float>& params)
-	{
-		typedef void (*tonemap_func)(Mat&, Mat&, std::vector<float>&);
-		const unsigned param_count[TONEMAP_COUNT] = {0, 1, 3, 3};
-		tonemap_func functions[TONEMAP_COUNT] = {
-			NULL, DragoMap, ReinhardDevlinMap, DurandMap};
-
-		Mat src = _src.getMat();
-		if(src.empty()) {
-			CV_Error(Error::StsBadArg, "Empty input image");
-		}
-		if(algorithm < 0 || algorithm >= TONEMAP_COUNT) {
-			CV_Error(Error::StsBadArg, "Wrong algorithm index");
-		}
-
-		_dst.create(src.size(), CV_32FC3);
-		Mat dst = _dst.getMat();
-		src.copyTo(dst);
-
-		double min, max;
-		minMaxLoc(dst, &min, &max);
-		if(max - min < 1e-10f) {
-			return;
-		}
-		dst = (dst - min) / (max - min);
-		if(functions[algorithm]) {
-			functions[algorithm](dst, dst, params);
-		}
-		minMaxLoc(dst, &min, &max);
-		dst = (dst - min) / (max - min);
-		float gamma = getParam(params, 0, 1.0f);		
-		pow(dst, 1.0f / gamma, dst);			
-	}
+
+static float getParam(const std::vector<float>& params, size_t i, float defval) 
+{
+    if(params.size() > i) {
+        return params[i];
+    } else {
+        return defval;
+    }
+        
+}
+
+static void DragoMap(Mat& src_img, Mat &dst_img, const std::vector<float>& params)
+{
+    float bias_value = getParam(params, 1, 0.85f);
+    Mat gray_img;
+    cvtColor(src_img, gray_img, COLOR_RGB2GRAY);
+    Mat log_img;
+    log(gray_img, log_img);
+    float mean = exp((float)sum(log_img)[0] / log_img.total());
+    gray_img /= mean;
+    log_img.release();
+
+    double max;
+    minMaxLoc(gray_img, NULL, &max);
+
+    Mat map;
+    log(gray_img + 1.0f, map);
+    Mat div;
+    pow(gray_img / (float)max, log(bias_value) / log(0.5f), div);
+    log(2.0f + 8.0f * div, div);
+    map = map.mul(1.0f / div);
+    map = map.mul(1.0f / gray_img);
+    div.release();
+    gray_img.release();
+
+    std::vector<Mat> channels(3);
+    split(src_img, channels);
+    for(int i = 0; i < 3; i++) {
+        channels[i] = channels[i].mul(map);
+    }
+    map.release();
+    merge(channels, dst_img);
+}
+
+static void ReinhardDevlinMap(Mat& src_img, Mat &dst_img, const std::vector<float>& params)
+{
+    float intensity   = getParam(params, 1, 0.0f);
+    float color_adapt = getParam(params, 2, 0.0f);
+    float light_adapt = getParam(params, 3, 1.0f);
+
+    Mat gray_img;
+    cvtColor(src_img, gray_img, COLOR_RGB2GRAY);
+    Mat log_img;
+    log(gray_img, log_img);
+
+    float log_mean = (float)sum(log_img)[0] / log_img.total();
+    double log_min, log_max;
+    minMaxLoc(log_img, &log_min, &log_max);
+    log_img.release();
+
+    double key = (float)((log_max - log_mean) / (log_max - log_min));
+    float map_key = 0.3f + 0.7f * pow((float)key, 1.4f);
+    intensity = exp(-intensity);
+    Scalar chan_mean = mean(src_img);
+    float gray_mean = (float)mean(gray_img)[0];
+
+    std::vector<Mat> channels(3);
+    split(src_img, channels);
+
+    for(int i = 0; i < 3; i++) {
+        float global = color_adapt * (float)chan_mean[i] + (1.0f - color_adapt) * gray_mean;
+        Mat adapt = color_adapt * channels[i] + (1.0f - color_adapt) * gray_img;
+        adapt = light_adapt * adapt + (1.0f - light_adapt) * global;
+        pow(intensity * adapt, map_key, adapt);
+        channels[i] = channels[i].mul(1.0f / (adapt + channels[i]));		
+    }
+    gray_img.release();
+    merge(channels, dst_img);
+}
+
+static void DurandMap(Mat& src_img, Mat& dst_img, const std::vector<float>& params)
+{
+    float contrast   = getParam(params, 1, 4.0f);
+    float sigma_color = getParam(params, 2, 2.0f);
+    float sigma_space = getParam(params, 3, 2.0f);
+
+    Mat gray_img;
+    cvtColor(src_img, gray_img, COLOR_RGB2GRAY);
+    Mat log_img;
+    log(gray_img, log_img);
+    Mat map_img;
+    bilateralFilter(log_img, map_img, -1, sigma_color, sigma_space);
+        
+    double min, max;
+    minMaxLoc(map_img, &min, &max);
+    float scale = contrast / (float)(max - min);
+
+    exp(map_img * (scale - 1.0f) + log_img, map_img);
+    log_img.release();
+    map_img = map_img.mul(1.0f / gray_img);
+    gray_img.release();
+
+    std::vector<Mat> channels(3);
+    split(src_img, channels);
+    for(int i = 0; i < 3; i++) {
+        channels[i] = channels[i].mul(map_img);
+    }
+    merge(channels, dst_img);
+}
+
+void tonemap(InputArray _src, OutputArray _dst, tonemap_algorithms algorithm, 
+             const std::vector<float>& params)
+{
+    typedef void (*tonemap_func)(Mat&, Mat&, const std::vector<float>&);
+    tonemap_func functions[TONEMAP_COUNT] = {
+        NULL, DragoMap, ReinhardDevlinMap, DurandMap};
+
+    Mat src = _src.getMat();
+    if(src.empty()) {
+        CV_Error(Error::StsBadArg, "Empty input image");
+    }
+    if(algorithm < 0 || algorithm >= TONEMAP_COUNT) {
+        CV_Error(Error::StsBadArg, "Wrong algorithm index");
+    }
+
+    _dst.create(src.size(), CV_32FC3);
+    Mat dst = _dst.getMat();
+    src.copyTo(dst);
+
+    double min, max;
+    minMaxLoc(dst, &min, &max);
+    if(max - min < 1e-10f) {
+        return;
+    }
+    dst = (dst - min) / (max - min);
+    if(functions[algorithm]) {
+        functions[algorithm](dst, dst, params);
+    }
+    minMaxLoc(dst, &min, &max);
+    dst = (dst - min) / (max - min);
+    float gamma = getParam(params, 0, 1.0f);		
+    pow(dst, 1.0f / gamma, dst);			
+}
 }
\ No newline at end of file

modules/photo/test/test_hdr.cpp

@@ -47,7 +47,7 @@
 using namespace cv;
 using namespace std;
 
-TEST(Photo_MakeHdr, regression)
+TEST(Photo_HdrFusion, regression)
 {
 	string folder = string(cvtest::TS::ptr()->get_data_path()) + "hdr/";
 	
@@ -75,6 +75,14 @@ TEST(Photo_MakeHdr, regression)
 	double max = 1.0;
 	minMaxLoc(abs(result - expected), NULL, &max);
 	ASSERT_TRUE(max < 0.01);
+
+	expected_path = folder + "grand_canal_exp_fusion.png";
+	expected = imread(expected_path);
+	ASSERT_FALSE(expected.empty()) << "Could not load input image " << expected_path;
+	exposureFusion(images, result);
+	result.convertTo(result, CV_8UC3, 255);
+	minMaxLoc(abs(result - expected), NULL, &max);
+	ASSERT_FALSE(max > 0);
 }
 
 TEST(Photo_Tonemap, regression)