C++实现三维高斯filter------3D Gaussian filter(大量matrix验证通过，结果完全与MATLAB完全一致)

工作上遇到了一个3D Gaussian filter conversion 问题，在MATLAB中可以直接调用B = imgaussfilt3(A)完成，且可以自定义三维的kernelSize和sigma。
但是在C++中，OpenCV只能处理二维，加上Z轴并不是简单的for-loop，除非Z轴值为1。
测试过手写convolution，但是速度很慢；也测试过使用MATLAB的code generator，速度也达不到项目要求。

自己通过openMP和OpenCV完成了conversion，300张图片只需要耗时0.5不到。

vector<cv::Mat> Aline3DPAR::Apply3DGaussianFilter(vector<cv::Mat>& images, int kernelSizeX, int kernelSizeY, int kernelSizeZ, double sigmaX, double sigmaY, double sigmaZ) {
	//  X, Y Gaussian filter
	vector<cv::Mat> gaussianResXY2D(imageCounts);
#pragma omp parallel for
	for (int i = 0; i < images.size(); ++i) {
		cv::Mat filteredSlice;
		images[i].convertTo(filteredSlice, CV_32F);
		cv::GaussianBlur(filteredSlice, gaussianResXY2D[i], cv::Size(kernelSizeX, kernelSizeY), sigmaX, sigmaY, cv::BORDER_REPLICATE);
	}
	int X = images[0].cols;
	int Y = images[0].rows;
	int Z = images.size();
	cv::Mat tempImageZ = cv::Mat::zeros(X * Y, Z, CV_32FC1);

	// Convert 3D volume to 2D
#pragma omp parallel for
	for (int y = 0; y < Y; ++y) {
		for (int x = 0; x < X; ++x) {
			for (int z = 0; z < Z; ++z) {
				tempImageZ.at<float>(y * X + x, z) = gaussianResXY2D[z].at<float>(y, x);
			}
		}
	}

	// Z dimension Gaussian filter 
	cv::Mat gaussianResZ;
	cv::GaussianBlur(tempImageZ, gaussianResZ, cv::Size(kernelSizeZ, 1), sigmaZ, 0, cv::BORDER_REPLICATE);


	// Convert 2D image to 3D volume
	std::vector<cv::Mat> finalOutput(Z);
#pragma omp parallel for
	for (int z = 0; z < Z; ++z) {
		finalOutput[z] = cv::Mat(Y, X, CV_64FC1);
	}

	#pragma omp parallel for
	for (int z = 0; z < Z; ++z) {
		for (int y = 0; y < Y; ++y) {
			for (int x = 0; x < X; ++x) {
				finalOutput[z].at<double>(y, x) = gaussianResZ.at<float>(y * X + x, z);
			}
		}
	}

	return finalOutput;
}