#include"eyemEdge.h"

int eyemEdgesPixel(EyemImage tpImage, double dThreshold)
{
	cv::Mat image = cv::Mat(tpImage.iHeight, tpImage.iWidth, MAKETYPE(tpImage.iDepth, tpImage.iChannels), tpImage.vpImage).clone();

	if (image.empty()) {
		return FUNC_IMAGE_NOT_EXIST;
	}
	int X = image.cols, Y = image.rows;
	uchar *upF = image.data;
	/*计算偏导*/
	double *dpFx = (double *)malloc(X*Y * sizeof(double));
	double *dpFy = (double *)malloc(X*Y * sizeof(double));
	/*梯度幅值*/
	double *dpMag = (double *)malloc(X*Y * sizeof(double));

	cv::parallel_for_(cv::Range(1, Y - 1), [&](const cv::Range& range) -> void {
		for (int y = range.start; y < range.end; y++) {
			for (int x = 1; x < X - 1; x++) {
				dpFx[(x)+(y)*X] = 0.5*((double)(upF[(x + 1) + (y)*X]) - (double)(upF[(x - 1) + (y)*X]));
				dpFy[(x)+(y)*X] = 0.5*((double)(upF[(x)+(y + 1)*X]) - (double)(upF[(x)+(y - 1)*X]));
				dpMag[(x)+(y)*X] = sqrt(dpFx[(x)+(y)*X] * dpFx[(x)+(y)*X] + dpFy[(x)+(y)*X] * dpFy[(x)+(y)*X]);
			}
		}
	});

	unsigned char *ucpLabel = (unsigned char *)malloc(X*Y * sizeof(unsigned char));
	memset(ucpLabel, 0, X*Y * sizeof(unsigned char));
	cv::Mat label = cv::Mat(Y, X, CV_8UC1, ucpLabel);

	cv::parallel_for_(cv::Range(1, Y - 1), [&](const cv::Range& range) -> void {
		for (int y = range.start; y < range.end; y++)
			for (int x = 1; x < (X - 1); x++) {
				if (dpMag[(x)+(y)*X] > dThreshold) {
					//判断边缘
					if (abs(dpFy[(x)+(y)*X]) >= abs(dpFx[(x)+(y)*X]) &&
						abs(dpFy[(x)+(y)*X]) >= abs(dpFy[(x)+(y - 1)*X]) && abs(dpFy[(x)+(y)*X]) > abs(dpFy[(x)+(y + 1)*X])) {
						ucpLabel[(x)+(y)*X] = 255;//垂直边缘-2
					}
					else if (abs(dpFx[(x)+(y)*X]) > abs(dpFy[(x)+(y)*X]) &&
						abs(dpFx[(x)+(y)*X]) >= abs(dpFx[(x - 1) + (y)*X]) && abs(dpFx[(x)+(y)*X]) > abs(dpFx[(x + 1) + (y)*X])) {
						ucpLabel[(x)+(y)*X] = 255;//水平边缘-1
					}
				}
			}
	});
	//标记连通域,根据不同连通域显示不同颜色
	cv::Mat labels;
	int ilabel = cv::connectedComponents(label, labels);
	//绘制连通域
	std::vector<cv::Vec3b> labelColor(ilabel);
	for (int i = 0; i < ilabel; i++) {
		labelColor[i] = cv::Vec3b(rand() % 256, rand() % 256, rand() % 256);
	}
	cv::Mat labelImage;
	cv::cvtColor(image, labelImage, cv::COLOR_GRAY2BGR);
	for (int y = 1; y < Y - 1; y++) {
		for (int x = 1; x < X - 1; x++) {
			int lb = labels.at<int>(y, x);
			if (lb != 0) {
				labelImage.at<cv::Vec3b>(y, x) = labelColor[lb];
			}
		}
	}
	//释放资源
	free((void *)ucpLabel);
	free((void *)dpFx);
	free((void *)dpFy);
	free((void *)dpMag);
	return FUNC_OK;
}

int eyemEdgesSubpixel(EyemImage tpImage, int iFilter, int iLow, int iHigh)
{
	cv::Mat image = cv::Mat(tpImage.iHeight, tpImage.iWidth, MAKETYPE(tpImage.iDepth, tpImage.iChannels), tpImage.vpImage).clone();
	if (image.empty()) {
		return FUNC_IMAGE_NOT_EXIST;
	}

	int X = image.cols, Y = image.rows;
	uchar *upF = image.data;
	/*计算偏导*/
	double *dpFx = (double *)malloc(X*Y * sizeof(double));
	double *dpFy = (double *)malloc(X*Y * sizeof(double));
	/*梯度幅值*/
	double *dpMag = (double *)malloc(X*Y * sizeof(double));

	cv::parallel_for_(cv::Range(1, Y - 1), [&](const cv::Range& range) -> void {
		for (int y = range.start; y < range.end; y++) {
			for (int x = 1; x < X - 1; x++) {
				dpFx[(x)+(y)*X] = 0.5*((double)(upF[(x + 1) + (y)*X]) - (double)(upF[(x - 1) + (y)*X]));
				dpFy[(x)+(y)*X] = 0.5*((double)(upF[(x)+(y + 1)*X]) - (double)(upF[(x)+(y - 1)*X]));
				dpMag[(x)+(y)*X] = sqrt(dpFx[(x)+(y)*X] * dpFx[(x)+(y)*X] + dpFy[(x)+(y)*X] * dpFy[(x)+(y)*X]);
			}
		}
	});

	unsigned char *ucpLabel = (unsigned char *)malloc(X*Y * sizeof(unsigned char));
	memset(ucpLabel, 0, X*Y * sizeof(unsigned char));
	cv::Mat label = cv::Mat(Y, X, CV_8UC1, ucpLabel);

	cv::parallel_for_(cv::Range(1, Y - 1), [&](const cv::Range& range) -> void {
		for (int y = range.start; y < range.end; y++)
			for (int x = 1; x < (X - 1); x++) {
				if (dpMag[(x)+(y)*X] > 25) {
					//判断边缘
					if (abs(dpFy[(x)+(y)*X]) >= abs(dpFx[(x)+(y)*X]) &&
						abs(dpFy[(x)+(y)*X]) >= abs(dpFy[(x)+(y - 1)*X]) && abs(dpFy[(x)+(y)*X]) > abs(dpFy[(x)+(y + 1)*X])) {
						ucpLabel[(x)+(y)*X] = 255;//垂直边缘-2
					}
					else if (abs(dpFx[(x)+(y)*X]) > abs(dpFy[(x)+(y)*X]) &&
						abs(dpFx[(x)+(y)*X]) >= abs(dpFx[(x - 1) + (y)*X]) && abs(dpFx[(x)+(y)*X]) > abs(dpFx[(x + 1) + (y)*X])) {
						ucpLabel[(x)+(y)*X] = 255;//水平边缘-1
					}
				}
			}
	});
	//所有边缘
	EyemOcsDXY edge;
	std::vector<EyemOcsDXY> edges;
	cv::parallel_for_(cv::Range(1, Y - 1), [&](const cv::Range& range) -> void {
		for (int y = range.start; y < range.end; y++)
			for (int x = 1; x < (X - 1); x++) {
				if (ucpLabel[(x)+(y)*X] != 0) {
					//判断边缘
					if (abs(dpFy[(x)+(y)*X]) >= abs(dpFx[(x)+(y)*X]) &&
						abs(dpFy[(x)+(y)*X]) >= abs(dpFy[(x)+(y - 1)*X]) && abs(dpFy[(x)+(y)*X]) > abs(dpFy[(x)+(y + 1)*X])) {
						//垂直边缘
						double a, b, c, u;
						a = dpMag[(x)+(y - 1)*X];
						b = dpMag[(x)+(y)*X];
						c = dpMag[(x)+(y + 1)*X];
						u = 0.5*(a - c) / (a - b - b + c);
						edge.dX = (float)x + 0.5f;
						edge.dY = (float)y + 0.5f + (float)u;
						edges.push_back(edge);
					}
					else if (abs(dpFx[(x)+(y)*X]) > abs(dpFy[(x)+(y)*X]) &&
						abs(dpFx[(x)+(y)*X]) >= abs(dpFx[(x - 1) + (y)*X]) && abs(dpFx[(x)+(y)*X]) > abs(dpFx[(x + 1) + (y)*X])) {
						//水平边缘
						double a, b, c, u;
						a = dpMag[(x - 1) + (y)*X];
						b = dpMag[(x)+(y)*X];
						c = dpMag[(x + 1) + (y)*X];
						u = 0.5*(a - c) / (a - b - b + c);
						edge.dX = (float)x + 0.5f + (float)u;
						edge.dY = (float)y + 0.5f;
						edges.push_back(edge);
					}
				}
			}
	});
	//释放资源
	free((void *)ucpLabel);
	free((void *)dpFx);
	free((void *)dpFy);
	free((void *)dpMag);
	return FUNC_OK;
}

int eyemSobelAmp(EyemImage tpImage, EyemImage &ImaAmp)
{
	cv::Mat image(tpImage.iHeight, tpImage.iWidth, CV_8UC1, tpImage.vpImage);

	if (image.empty()) {
		return FUNC_IMAGE_NOT_EXIST;
	}
	cv::Mat dx, dy;
	cv::spatialGradient(image, dx, dy);

	cv::Mat mag;
	cv::magnitude(dx, dy, mag);

	return FUNC_OK;
}

int eyemSkeleton(EyemImage tpImage, cv::Mat &skeleton)
{
	cv::Mat image(tpImage.iHeight, tpImage.iWidth, CV_8UC1, tpImage.vpImage);

	if (image.empty()) {
		return FUNC_IMAGE_NOT_EXIST;
	}
	return FUNC_OK;
}

int	eyemAutoCanny(EyemImage tpImage, double dSigma = 0.33)
{
	cv::Mat image(tpImage.iHeight, tpImage.iWidth, CV_8UC1, tpImage.vpImage);

	if (image.empty()) {
		return FUNC_IMAGE_NOT_EXIST;
	}
	cv::Mat F;
	cv::GaussianBlur(image, F, cv::Size(3, 3), 0, 0);
	////get the median value of the matrix
	//double v = medianMat(output);
	////generate the thresholds
	//int lower = (int)std::max(0.0, (1, 0 - sigma)*v);
	//int upper = (int)std::min(255.0, (1, 0 + sigma)*v);
	////apply canny operator
	//cv::Canny(output, output, lower, upper, 3);
	return FUNC_OK;
}