步骤:1、计算累加直方图;2、计算一阶累积矩;3、计算图像的总体灰度平均值;4、计算每个灰度值对应的方差;5、找到目标阈值。
C++:
#include#includeusing namespace std;using namespace cv;//计算图像灰度直方图Mat calcgrayhist(const Mat&image){Mat histogram = Mat::zeros(Size(256, 1), CV_32SC1);//图像宽高int rows = image.rows;int cols = image.cols;for (int i = 0; i < rows; i++){for (int j = 0; j < rows; j++){int index = int(image.at(i, j));histogram.at(0, index) += 1;}}return histogram;};//OTSUint OTSU(const Mat&image,Mat &OTSU_image){int rows = image.rows;int cols = image.cols;Mat histogram = calcgrayhist(image);//归一化直方图Mat normhist;histogram.convertTo(normhist, CV_32FC1, 1.0 / (rows*cols), 0.0);//计算累加直方图和一阶累积矩Mat zeroaccumulate = Mat::zeros(Size(256, 1), CV_32FC1);Mat oneaccumulate = Mat::zeros(Size(256, 1), CV_32FC1);for (int i = 0; i < 256; i++){if (i == 0){zeroaccumulate.at(0, i) = normhist.at(0, i);oneaccumulate.at(0, i) =i * normhist.at(0, i);}else{zeroaccumulate.at(0, i) = zeroaccumulate.at(0, i-1)+normhist.at(0, i);oneaccumulate.at(0, i) = oneaccumulate.at(0, i - 1)+i* normhist.at(0, i);}}//计算间类方差Mat variance= Mat::zeros(Size(256, 1), CV_32FC1);float mean = oneaccumulate.at(0, 255);for (int i = 0; i < 255; i++){if (zeroaccumulate.at(0, i) == 0 || zeroaccumulate.at(0, i) == 1)variance.at(0, i) = 0;else{float cofficient = zeroaccumulate.at(0, i)*(1.0 -zeroaccumulate.at(0, i));variance.at(0, i) = pow(mean*zeroaccumulate.at(0, i)- oneaccumulate.at(0, i), 2.0) / cofficient;}}//找到阈值;Point maxloc;//计算矩阵中最大值minMaxLoc(variance, NULL, NULL, NULL, &maxloc);int otsuthreshold = maxloc.x;threshold(image, OTSU_image, otsuthreshold, 255, THRESH_BINARY);return otsuthreshold;};int main(){Mat img, gray_img,dst;img = imread("D:/testimage/orange.jpg");cvtColor(img, gray_img, COLOR_BGR2GRAY);OTSU(gray_img,dst);imshow("otsu image", dst);waitKey(0);return 0;};
结果:
![fdefd6158ff4e832e2386409e1e69e3e.png](https://i-blog.csdnimg.cn/blog_migrate/f230ffec581592fe9958108125961f3e.jpeg)
Otsu
Python:
import cv2 as cvimport numpy as npimport mathdef calc_grayhist(image): #图像宽高 rows,cols=image.shape grayhist=np.zeros([256],np.uint64) for i in range(rows): for j in range(cols): grayhist[image[i][j]]+=1 return grayhistdef OTSU(image): rows, cols = image.shape grayhist=calc_grayhist(image) #归一化直方图 uniformgrayhist=grayhist/float(rows*cols) #计算零阶累积矩和一阶累积矩 zeroaccumulat = np.zeros([256],np.float32) oneaccumulat = np.zeros([256], np.float32) for k in range(256): if k==0: zeroaccumulat[k]=uniformgrayhist[0] oneaccumulat[k]=k*uniformgrayhist[0] else: zeroaccumulat[k]=zeroaccumulat[k-1]+uniformgrayhist[k] oneaccumulat[k]=oneaccumulat[k-1]+k*uniformgrayhist[k] #计算间类方差 variance=np.zeros([256],np.float32) for k in range(256): if zeroaccumulat[k]==0 or zeroaccumulat[k]==1: variance[k]=0 else: variance[k]=math.pow(oneaccumulat[255]*zeroaccumulat[k]- oneaccumulat[k],2)/(zeroaccumulat[k]*(1.0-zeroaccumulat[k])) threshLoc=np.where(variance[0:255]==np.max(variance[0:255])) thresh=threshLoc[0][0] #阈值分割 threshold=np.copy(image) threshold[threshold>thresh]=255 threshold[threshold <= thresh] = 0 return thresholdif __name__=="__main__": img=cv.imread("D:/testimage/cxk.jpg") gray_dst = cv.cvtColor(img, cv.COLOR_BGR2GRAY) Otsu_dst=OTSU(gray_dst) cv.imshow("gray dst",gray_dst) cv.imshow("Otsu dst", Otsu_dst) cv.imwrite("D:/testimage/result-cxk.jpg",Otsu_dst) cv.waitKey(0) cv.destroyAllWindows()
CSDN:@WitnesS
感谢支持