Auto Threshold (Huang2) C++

Java -> C++ 로 변경

 

Code

int Huang2(float* histo, int n_length) {
	// Implements Huang's fuzzy thresholding method 
	// Uses Shannon's entropy function (one can also use Yager's entropy function) 
	// Huang L.-K. and Wang M.-J.J. (1995) "Image Thresholding by Minimizing  
	// the Measures of Fuzziness" Pattern Recognition, 28(1): 41-51
	// Reimplemented (to handle 16-bit efficiently) by Johannes Schindelin Jan 31, 2011

	// find first and last non-empty bin
	int first = 1, last;
	for (first = 0; first < n_length && histo[first] == 0; first++)
		; // do nothing
	for (last = n_length - 1; last > first && histo[last] == 0; last--)
		; // do nothing
	if (first == last)
		return 0;

	// calculate the cumulative density and the weighted cumulative density
	double* S = nullptr;
	double* W = nullptr;
	double* Smu = nullptr;

	S = new double[last + 1];
	W = new double[last + 1];
	
	int start = max(1, first);
	S[start - 1] = histo[start - 1];
	W[start - 1] = histo[start - 1];

	for (int i = start; i <= last; i++)
	{
		S[i] = S[i - 1] + histo[i];
		W[i] = W[i - 1] + i * histo[i];
	}

	// precalculate the summands of the entropy given the absolute difference x - mu (integral)
	double C = last - first;
	int Smu_length = last + 1 - first;
	Smu = new double[Smu_length];
	
	for (int i = 1; i < Smu_length; i++) 
	{
		double mu = 1 / (1 + abs(i) / C);
		Smu[i] = -mu * log(mu) - (1 - mu) * log(1 - mu);
	}

	// calculate the threshold
	int bestThreshold = 0;
	double bestEntropy = DBL_MAX;

	for (int threshold = first; threshold <= last; threshold++) 
	{
		double entropy = 0;
		int mu = (int)round(W[threshold] / S[threshold]);

		for (int i = first; i <= threshold; i++)
			entropy += Smu[abs(i - mu)] * histo[i];

		mu = (int)round((W[last] - W[threshold]) / (S[last] - S[threshold]));

		for (int i = threshold + 1; i <= last; i++)
			entropy += Smu[abs(i - mu)] * histo[i];

		if (bestEntropy > entropy) {
			bestEntropy = entropy;
			bestThreshold = threshold;
		}
	}

	if (S)
		delete[] S;
	
	if (W)
		delete[] W;

	if (Smu)
		delete[] Smu;

	return bestThreshold;
}

 

사용 예시

		std::string str_image_path;
		cv::Mat raw_img = cv::imread(str_image_path, -1);
		cv::Mat raw_img_8bit;
		
		// 8Bit 1Channel 만 가능...
		raw_img.convertTo(raw_img_8bit, CV_8UC1, 256.0 / 65536.0);

		int histSize = 256;    
		float range[] = { 0, 255 };
		const float* histRange = { range };
		cv::Mat hist;
		cv::calcHist(&raw_img_8bit, 1, 0, cv::Mat(), hist, 1, &histSize, &histRange);

		int n_length = raw_img.rows * raw_img.cols;
		uint8_t threshold_value = Huang2((float*)hist.data, 256);
	
		cv::Mat threshold_img;
		cv::threshold(raw_img_8bit, threshold_img, threshold_value, 255, cv::THRESH_BINARY);

 

imageJ 관련 문서 : https://imagej.net/plugins/auto-threshold

github link : https://github.com/fiji/Auto_Threshold/blob/master/src/main/java/fiji/threshold/Auto_Threshold.java