NIQE c++实现

根据https://github.com/Shiyuuuu/image_metrics/blob/main/niqe.py代码改编，结果精度会不同，如果可以调整的很好的话，欢迎告知

测试代码：

①多线程调用模式

//多线程调用
void test_speed()
{
	string Path = R"(E:\2024Work\2\image_result\pathologic\0clear\0_0.jpg)";
	Mat input_mat = imread(Path);
	string file_path = R"(E:\2024Work\2\ILNIQE\test\Param.txt)";
	ScanQualityEvaluator sqe;
	sqe.start_calc(6, file_path);  //启动6个线程，来计算
 
	int loop_counts = 200;
	auto st = CMyLog::getCurrentTime();
	for (int i = 0; i < loop_counts; i++)
	{
		// 图片添加到队列
		sqe.push_back(input_mat.data, input_mat.rows, input_mat.cols,1);
	}
	sqe.stop_calc();  //等待停止
	double total_times = CMyLog::getSpanTimes(st, CMyLog::getCurrentTime());
	cout << total_times / loop_counts << endl;
}

②单测模式

 
void test_1()
{
	string Path = R"(E:\2024Work\2\image_result\pathologic\1blur\0_0.jpg)";
	Mat rgbImg = imread(Path);
    //模型文件
	string file_path = R"(E:\2024Work\2\ILNIQE\test\Param.txt)";
 
	ScanQualityEvaluator::_CalculateILNIQE_Ready(file_path);
 
	ScanQualityEvaluator::_CalculateILNIQE(cut_img);
}

代码：

#pragma once
#include <opencv2/opencv.hpp>
#include "../BaseCode/BnvQueue.hpp"
 
//计算扫描质量的接口
class ScanQualityEvaluator
{
public:
	//Path 是均值和方差文件
	static double _CalculateILNIQE(const cv::Mat img, int block_size_h = 96, int block_size_w = 96);
	static int _CalculateILNIQE_Ready(std::string path);
 
	
	//type=0  raw 1 rgb
	int push_back(unsigned char* data, int rows, int cols, int type = 0);
	//启动多少个线程来处理
	int start_calc(int thread_nums,std::string path);
	double stop_calc(); //要等待线程结束  返回计算的值
 
private:
	static std::tuple<float, float, float> _estimate_aggd_param(const cv::Mat& block);
	static std::vector<float> _compute_feature(const cv::Mat& block);
 
private:
	void _Start();
	double _Stop();
	int setStatus(int status);
	int config(int thread_nums, std::string path);
private:
	int m_status;  //0 是表示启动  1表示停止
	std::mutex m_img_q_vals_lock;
	std::vector<double> m_img_q_vals;
	BnvQueue<cv::Mat> m_rgb_img_queue;
	std::vector<std::thread> m_threads;
 
	//NIQE的均值，方差矩阵
	static cv::Mat mu_pris_param;
	static cv::Mat cov_pris_param;
	static cv::Mat gaussian_window;
	static std::vector<float> gam;
	static std::vector<float> r_gam;
};

#include "ScanQualityEvaluator.hpp"
#include <vector>
#include <numeric>
#include "MyLog.h"
#include <fstream>
#include <Tool.h>
using namespace cv;
using namespace std;
 
cv::Mat ScanQualityEvaluator::mu_pris_param;
cv::Mat ScanQualityEvaluator::cov_pris_param;
cv::Mat ScanQualityEvaluator::gaussian_window;
std::vector<float> ScanQualityEvaluator::gam;
std::vector<float> ScanQualityEvaluator::r_gam;
 
std::tuple<float, float, float> ScanQualityEvaluator::_estimate_aggd_param(const Mat& block) {
	// 将图像块展平
	Mat block_flattened = block.clone().reshape(1, 1); // 展平为行向量
 
	//计算标准差
	float negative_sum = 0.0, positive_sum = 0.0;
	int negative_counts = 0, positive_counts = 0;
 
 
	float* p = block_flattened.ptr<float>(0);
	float val;
	for (int i = 0; i < block_flattened.cols; ++i)
	{
		val = *p++;
		if (val < 0)
		{
			negative_sum += val *val;
			++negative_counts;
		}
		else
		{
			positive_sum += val * val;
			++positive_counts;
		}
	}
 
 
	float left_std = sqrt(negative_sum / negative_counts);
	float right_std = sqrt(positive_sum/ positive_counts);
 
	float gammahat = left_std / right_std;
	float rhat = std::pow(mean(abs(block_flattened))[0], 2) / mean(block_flattened.mul(block_flattened))[0];
	float rhatnorm = (rhat * (std::pow(gammahat, 3) + 1) * (gammahat + 1)) / std::pow((std::pow(gammahat, 2) + 1), 2);
	size_t array_position = std::distance(r_gam.begin(), std::min_element(r_gam.begin(), r_gam.end(), [&](float a, float b) { return std::abs(a - rhatnorm) < std::abs(b - rhatnorm); }));
 
	float alpha = gam[array_position];
	float beta_l = left_std * std::sqrt(std::tgamma(1 / alpha) / std::tgamma(3 / alpha));
	float beta_r = right_std * std::sqrt(std::tgamma(1 / alpha) / std::tgamma(3 / alpha));
 
	return std::make_tuple(alpha, beta_l, beta_r);
}
 
 
//50ms
std::vector<float> ScanQualityEvaluator::_compute_feature(const Mat& block) {
	std::vector<float> feat;
	// 计算 alpha、beta_l 和 beta_r 参数
	float alpha, beta_l, beta_r;
	std::tie(alpha, beta_l, beta_r) = _estimate_aggd_param(block);
 
	// 添加 alpha 和 (beta_l + beta_r) / 2 到特征向量中
	feat.push_back(alpha);
	feat.push_back((beta_l + beta_r) / 2);
 
	// 定义偏移量
	std::vector<std::vector<int>> shifts = { {0, 1}, {1, 0}, {1, 1}, {1, -1} };
 
	// 遍历每个偏移量
	// 滚动图像块
	
	for (size_t i = 0; i < shifts.size(); ++i) {
		Mat shifted_block = block.clone();
		Mat cut_block;
		int y_shift = shifts[i][0];   //y_shift 上下行移动
		int x_shift = shifts[i][1];   //x_shift 是行内移动
		if (y_shift != 0)
		{
			if (y_shift > 0)
			{
				cut_block = shifted_block(Rect(0, shifted_block.rows - y_shift, shifted_block.cols, y_shift));
				shifted_block = shifted_block(Rect(0, 0, shifted_block.cols, shifted_block.rows - y_shift));
				vconcat(cut_block, shifted_block, shifted_block);
			}
			else
			{
				cut_block = shifted_block(Rect(0, 0, shifted_block.cols, abs(y_shift)));
				shifted_block = shifted_block(Rect(0, abs(y_shift), shifted_block.cols, shifted_block.rows - abs(y_shift)));
				vconcat(shifted_block, cut_block, shifted_block);
			}
		}
		if (x_shift != 0)
		{
			if(x_shift > 0)
			{
				cut_block = shifted_block(Rect(shifted_block.cols - x_shift, 0, x_shift, shifted_block.rows));
				shifted_block = shifted_block(Rect(0, 0, shifted_block.cols - x_shift, shifted_block.rows));
				hconcat(cut_block, shifted_block, shifted_block);
			}
			else
			{
				cut_block = shifted_block(Rect(0, 0, abs(x_shift), shifted_block.rows));
				shifted_block = shifted_block(Rect(abs(x_shift), 0, shifted_block.cols - abs(x_shift), shifted_block.rows));
				hconcat(shifted_block, cut_block, shifted_block);
			}
		}
 
		// 计算 alpha、beta_l 和 beta_r 参数
		std::tie(alpha, beta_l, beta_r) = _estimate_aggd_param(block.mul(shifted_block));
 
		// Eq. 8
		float mean = (beta_r - beta_l) * (std::tgamma(2 / alpha) / std::tgamma(1 / alpha));
		feat.push_back(alpha);
		feat.push_back(mean);
		feat.push_back(beta_l);
		feat.push_back(beta_r);
	}
 
 
	return feat;
}
 
 
double ScanQualityEvaluator::_CalculateILNIQE(const cv::Mat img, int block_size_h, int block_size_w)
{
	auto st = CMyLog::getCurrentTime();
	// 将图像转换为浮点型，并归一化到[0, 1]范围
	Mat floatImage;
	img.convertTo(floatImage, CV_32F, 1.0f / 255);
	
	// 使用 OpenCV 函数实现颜色转换和亮度调整
	Mat y_Mat;
	transform(floatImage, y_Mat, cv::Matx13f(65.481f, 128.553f, 24.966f));
	y_Mat += 16.0f;
 
 
	float* p;
	for (int y = 0; y < y_Mat.rows; ++y) {
		p = y_Mat.ptr<float>(y);
		for (int x = 0; x < y_Mat.cols; ++x) {
			*p = std::round(*p);
			++p;
		}
	}
 
 
	int h = y_Mat.rows;
	int w = y_Mat.cols;
 
	int num_block_h = std::floor(h / block_size_h);
	int num_block_w = std::floor(w / block_size_w);
 
	int cut_width = block_size_w * num_block_w;
	int cut_height = block_size_h * num_block_h;
 
	Rect rect(0, 0, cut_width, cut_height);
	y_Mat = y_Mat(rect); 
 
 
 
	vector<float> scales{ 1,2 };
	vector <Mat> distparam;
 
 
	for (auto scale : scales)
	{
		Mat mu;
		filter2D(y_Mat, mu, -1, gaussian_window, Point(-1, -1), 0, BORDER_REPLICATE);
 
		Mat mu_squared = mu.mul(mu);
 
		Mat img_squared = y_Mat.mul(y_Mat);
 
		Mat conv_result;
		filter2D(img_squared, conv_result, -1, gaussian_window, Point(-1, -1), 0, BORDER_REPLICATE);
 
		Mat sigma;
 
		sqrt(abs(conv_result - mu_squared), sigma);
 
		Mat img_nomalized = (y_Mat - mu) / (sigma + 1);
 
		vector <vector<float>> feat;
		vector <float> all_feat;
		for (int idx_w = 0; idx_w < num_block_w; ++idx_w)
		{
			for (int idx_h = 0; idx_h < num_block_h; ++idx_h)
			{
				Rect block_rect(idx_w * block_size_w / scale, idx_h * block_size_h / scale,
					block_size_w / scale, block_size_h / scale);
				Mat block = img_nomalized(block_rect);
 
				// 计算当前块的特征并存储到 feat 中
				auto vals = _compute_feature(block);
				feat.push_back(vals);
				all_feat.insert(all_feat.end(), vals.begin(), vals.end());
			}
		}
 
 
		if (scale == 1)
		{
			resize(y_Mat / 255., y_Mat, Size(), 0.5, 0.5, INTER_LINEAR);
			y_Mat = y_Mat * 255;
		}
 
		Mat featMat(feat.size(), feat[0].size(), CV_32FC1, all_feat.data());
		distparam.emplace_back(featMat.clone());
	}
 
 
	Mat distparam_mat = distparam[0];
	for (int i = 1; i < distparam.size(); ++i)
	{
		hconcat(distparam_mat, distparam[i], distparam_mat);
	}
 
	Mat mu_distparam;
	reduce(distparam_mat, mu_distparam, 0, REDUCE_AVG);
 
 
	Mat cov_distparam;
	calcCovarMatrix(distparam_mat, cov_distparam, mu_distparam, COVAR_NORMAL | COVAR_ROWS | COVAR_SCALE);
 
	cov_distparam.convertTo(cov_distparam, CV_32F);
	mu_distparam.convertTo(mu_distparam, CV_32F);
 
	// 计算 niqe quality
	Mat invcov_param = (cov_pris_param + cov_distparam) / 2;
	invert(invcov_param, invcov_param, DECOMP_SVD);
 
	Mat diff_mat = mu_pris_param - mu_distparam;
 
	Mat quality = diff_mat * invcov_param * diff_mat.t();
	double quality_value = std::sqrt(quality.at<float>(0, 0));
 
	std::cout << "Quality: " << quality_value << std::endl;
	做个规约转到2-10转为100-0%   分数越小越好
	//quality_value = quality_value < 2 ? 2 : quality_value;
	//quality_value = quality_value >10 ? 10 : quality_value;
 
	//quality_value = 1 - (quality_value - 2) / 8;
 
 
	return quality_value;
 
}
 
int ScanQualityEvaluator::_CalculateILNIQE_Ready(std::string path)
{
	if (mu_pris_param.empty() || cov_pris_param.empty())
	{
		ifstream in_file(path);
		vector<string> lines;
		string line;
		while (getline(in_file, line))
		{
			lines.push_back(line);
		}
		if (lines.size() < 37)
		{
			return -1;
		}
		in_file.close();
		vector<float> mu_pris_param_vec;
		Tool::split_f(lines[0], mu_pris_param_vec, ' ');
 
 
		vector<float> cov_pris_param_vec;
		for (int i = 1; i < 37; ++i)
		{
			vector<float> tmp;
			Tool::split_f(lines[i], tmp, ' ');
			cov_pris_param_vec.insert(cov_pris_param_vec.end(), tmp.begin(), tmp.end());
		}
 
		int param_counts = mu_pris_param_vec.size();
		mu_pris_param = Mat(1, param_counts, CV_32FC1, mu_pris_param_vec.data()).clone();
		cov_pris_param = Mat(param_counts, param_counts, CV_32FC1, cov_pris_param_vec.data()).clone();
	}
 
 
	if (gaussian_window.empty())
	{//生成高斯窗口
		gaussian_window = getGaussianKernel(7, 7.0 / 6, CV_32F);
 
		// 将高斯窗口变为二维矩阵
		gaussian_window = gaussian_window * gaussian_window.t();
 
		// 归一化高斯窗口
		gaussian_window /= sum(gaussian_window)[0];
	}
 
 
	// 计算参数
	for (float i = 0.2; i <= 10.001; i += 0.001) {
		gam.push_back(i);
	}
	std::vector<float> gam_reciprocal(gam.size());
	std::transform(gam.begin(), gam.end(), gam_reciprocal.begin(), [](float x) { return 1.0f / x; });
 
	r_gam.resize(gam.size());
	for (size_t i = 0; i < gam.size(); ++i) {
		float gamma_reciprocal_2 = std::tgamma(gam_reciprocal[i] * 2);
		float gamma_reciprocal_3 = std::tgamma(gam_reciprocal[i] * 3);
		r_gam[i] = std::pow(gamma_reciprocal_2, 2) / (std::tgamma(gam_reciprocal[i]) * gamma_reciprocal_3);
	}
 
	return 0;
}
 
 
int ScanQualityEvaluator::config(int thread_nums, std::string path)
{
	if (_CalculateILNIQE_Ready(path) != 0)
	{
		return -1;
	}
	m_rgb_img_queue.init();
	m_threads.clear();
	m_img_q_vals.clear();
	setStatus(0);
	for (int i = 0; i < thread_nums; ++i)
	{
		m_threads.push_back(thread(&ScanQualityEvaluator::_Start, this));
	}
	return 0;
}
 
//raw 转 RGB
int ScanQualityEvaluator::push_back(unsigned char* data, int rows, int cols, int type)
{
	Mat rgbImg;
	if (type == 0)
	{
		const Mat srcImg(rows, cols, CV_8UC1, data);
		cvtColor(srcImg, rgbImg, COLOR_BayerRG2RGB);
	}
	else  //rgb
	{
		rgbImg = Mat(rows, cols, CV_8UC3, data).clone();
	}
 
	//裁取到最长边1024  加速处理
	int max_size = max(rgbImg.rows, rgbImg.cols);
	double scale = 1024.0 / max_size;
	int new_rows = scale * rgbImg.rows;
	int new_cols = scale * rgbImg.cols;
 
	//在_img中心裁取new_rows*new_cols的图像
	int x = (rgbImg.cols - new_cols) / 2;
	int y = (rgbImg.rows - new_rows) / 2;
	Mat cut_img = rgbImg(Rect(x, y, new_cols, new_rows));
 
	//如果不要裁取就解开注释
	//Mat cut_img = rgbImg;
	
	m_rgb_img_queue.push_back(cut_img.clone());
	return 0;
}
 
int ScanQualityEvaluator::start_calc(int thread_nums, std::string path)
{
	return config(thread_nums, path);
}
 
double ScanQualityEvaluator::stop_calc()
{
	return _Stop();
}
 
int ScanQualityEvaluator::setStatus(int status)
{
	m_status = status;
	return m_status;
}
 
void ScanQualityEvaluator::_Start()
{
	Mat rgbImg;
	while (!(m_status == 1 && m_rgb_img_queue.empty()))
	{
		if (m_rgb_img_queue.pop(rgbImg))
		{
			double ilniqe = _CalculateILNIQE(rgbImg);
			{
				lock_guard<mutex> lk(m_img_q_vals_lock);
				m_img_q_vals.push_back(ilniqe);
			}
		}
		else
		{
			this_thread::sleep_for(10ms);
		}
	}
}
 
double ScanQualityEvaluator::_Stop()
{
	setStatus(1);
	m_rgb_img_queue.finish();
	for (auto& th : m_threads)
	{
		th.join();
	}
	m_threads.clear();
	int offset = 0;
	sort(m_img_q_vals.begin(), m_img_q_vals.end());
 
	if (m_img_q_vals.size() > 100)  //超过100个，前20%不参与运算
	{
		offset = m_img_q_vals.size() * 0.2;
	}
	double mean = std::accumulate(m_img_q_vals.begin() + offset, m_img_q_vals.end(), 0.0) / (m_img_q_vals.size() - offset);
 
	return mean;
}
 
 
 

Param.txt文件的制作(Matlab来生成的)

setDir = fullfile('E:\2024Work\2\ILNIQE\扫描质量评分\all');
imds = imageDatastore(setDir,'FileExtensions',{'.jpg'});
 
model = fitniqe(imds);
 
mean = model.Mean;
cov = model.Covariance;
 
 
 dlmwrite('Param.txt', mean, ' ');
 dlmwrite('Param.txt', cov,'delimiter', ' ', '-append');