参考链接：
- https://blog.csdn.net/LucasXu01/article/details/90764759?utm_medium=distribute.pc_relevant.none-task-blog-2defaultbaidujs_baidulandingword~default-0-90764759-blog-125619862.235v39pc_relevant_3m_sort_dl_base3&spm=1001.2101.3001.4242.1&utm_relevant_index=3
- https://www.jianshu.com/p/12ce9032d904
完整代码下载地址：https://download.csdn.net/download/m0_37567738/88622880
代码：
function labelIm = quantizeFeats (featIm, meanFeats)
height = size(featIm, 1);
width = size(featIm, 2);
dimension = size(featIm, 3);
% 计算原始图像 至 meanFeats 的距离，生成隶属度矩阵 labelIm
labelIm = ones(height, width);
for h=1:height
for w=1:width
minDis = -1;
for k=1:size(meanFeats, 1)
dis = 0;
for d=1:dimension
dis = dis + power(featIm(h, w, d) - meanFeats(k, d), 2);
end
dis = sqrt(dis);
if minDis < 0 || dis < minDis
labelIm(h, w) = k;
minDis = sqrt(dis);
end
end
end
end
return
% 给定一个h * w *d 的矩阵featim，其中h 和w 为原始图像的高度和宽度，d 表示图像中每一个像素点所提取的特征向量的维数。
% 给定一个有k 个聚类中心点的矩阵meanFeatures(矩阵维度是k*d)，其中每个中心点都是一个d 维的行向量（矩阵的一行），
% 将输入图片中的每个像素映射到其所归属的k-means 某一个中心点去。
% 函数的返回值定义为labelim，labelim 是一个h*w 的整数矩阵，用来表明每一个像素所属的聚类中心标号(1...k)。
function [textons] = createTextons(imStack, bank, k)
[row, col] = size(imStack);
bankNum = size(bank, 3);
% 用滤波器过滤元胞数组中的n个灰度图，得到n个 bankNum 维的过滤响应样本
% 再将这n个过滤响应样本组合成一个大的样本集 textonsData
textonsData = [];
for i=1:row
for j=1:col
im = imStack{i, j};
im = im2double(im);
responses = zeros(size(im, 1), size(im, 2), bankNum);
for r=1:bankNum
responses(:,:,r)=conv2(im,double(bank(:,:,r)),'same');
end
X = reshape(responses, size(responses,1)* size(responses,2), bankNum);
Xrow = size(X, 1);
ranX = X(randperm(Xrow, ceil(Xrow/1000)),:);
textonsData = [textonsData; ranX];
end
end
% 使用 kmeans 方法聚合出含有k个纹理基元的纹理编码集
[~, textons] = kmeans(textonsData, k);
return;
% 给定一个长度为n 且包含n 个灰度值图像的元胞数组imStack，以及滤波器组bank，
% 基于所有n 个图片的过滤响应样本计算一个纹理基元编码集（如一组量化的滤波器组响应）。
% 其中bank 是一个包含d 个滤波器的m*m*d 矩阵，每个滤波器的大小为m*m，textons 是一个k*d 的矩阵，
% 其中每一行代表一个纹理特征，如一个量化滤波器组的响应。
function [featIm] = extractTextonHists(origIm, bank, textons, winSize)
origIm = im2double(origIm);
[row, col] = size(origIm);
% 图像过滤，生成滤波器组响应 responses
bankNum = size(bank, 3);
responses = zeros(row, col, bankNum);
for r=1:bankNum
responses(:,:,r)=conv2(origIm,double(bank(:,:,r)),'same');
end
% 计算 滤波器组响应 与 纹理基元编码集 的距离，并生成隶属度矩阵 feattexton
X = reshape(responses, size(responses,1)* size(responses,2), bankNum);
dis2textons = dist2(textons, X);
[~, indxtexton] = max(dis2textons);
feattexton = reshape(indxtexton, row, col);
% 图像边界处理
if winSize > 1
colNumLeft = floor((winSize-1)/2);
colNumRight = ceil((winSize-1)/2);
for i=1:colNumLeft
feattexton = [feattexton(:,1), feattexton];
feattexton = [feattexton(1,:); feattexton];
end
for i=1:colNumRight
feattexton = [feattexton, feattexton(:,size(feattexton, 2))];
feattexton = [feattexton; feattexton(size(feattexton, 1), :)];
end
else
colNumLeft = 0;
colNumRight = 0;
end
% 生成纹理柱状图 featIm
featIm = zeros(row, col, size(textons, 1));
for i=(1+colNumLeft):(size(feattexton, 1)-colNumRight)
for j=(1+colNumLeft):(size(feattexton, 2)-colNumRight)
window = feattexton((i-colNumLeft):(i+colNumRight), (j-colNumLeft):(j+colNumRight));
frequency = tabulate(window(:));
for k=1:size(frequency, 1)
textonIndex = int64(frequency(k, 1));
count = int32(frequency(k, 2));
featIm(i-colNumLeft, j-colNumLeft, textonIndex) = featIm(i-colNumLeft, j-colNumLeft, textonIndex)+count;
end
end
end
return;
% 给定一张灰度图像，一个滤波器组，一份纹理编码集，构建纹理柱状图。
% 对于每个像素，基于每个纹理在其邻近范围（定义在固定大小的winSize 内的局
% 部窗口）内出现的频率。其中，texton是k*d 的矩阵。
function [colorLabelIm, textureLabelIm] = ...
compareSegmentations(origIm, bank, textons, winSize, ...
numColorRegions, numTextureRegions)
origIm = im2double(origIm);
% 获取颜色集，生成基于颜色分割的标签矩阵
colordata = reshape(origIm, size(origIm, 1)*size(origIm, 2), size(origIm, 3));
opts = statset('Display','final','MaxIter',1000);
[~, colorCenter] = kmeans(colordata, numColorRegions, 'Options', opts);
colorLabelIm = quantizeFeats(origIm, colorCenter);
% 获取纹理柱状图，计算纹理特征，生成基于纹理分割的标签矩阵
featIm = extractTextonHists(rgb2gray(origIm), bank, textons, winSize);
featImData = reshape(featIm, size(featIm, 1)*size(featIm, 2), size(featIm, 3));
[~, textureCenter] = kmeans(featImData, numTextureRegions, 'Options', opts);
textureLabelIm = quantizeFeats(featIm, textureCenter);
return
origIm5 =imread('dress.jpg');
origIm1 =imread('circle-im-1.jpg');
origIm =imread('butterfly.jpg');
origIm3 =imread('circle-im-2.jpg');
origIm4 =imread('gumballs.jpg');
%[width,height] = size(origIm);
[height,width] = size(origIm);
subplot(1,3,1);
imshow(origIm)
title('原图像');
%bank = load('filterBank.mat');
bank = cell2mat(struct2cell(load('filterBank.mat')));
str =class(bank); %判断bank类型 disp(str);
imStack = {rgb2gray(origIm)} ; %
%imStack = {rgb2gray(origIm),rgb2gray(origIm1),rgb2gray(origIm2),rgb2gray(origIm3),rgb2gray(origIm4)};
% disp(bank); %打印过滤器
[m,n,l] =size(bank);
bankNum1 = size(bank, 1);
bankNum2 = size(bank, 2);
bankNum = size(bank, 3);
disp(bankNum);
% disp(bankNum1);
% disp(m);
% disp(n);
% disp(l);
winSize = 10;
numColorRegions =10; % 聚类数
numTextureRegions = 50 ; % 聚类数
%bank =ones(5,6,7);
%textons_x = 600;
%textons =ones(textons_x,bankNum);
k = 100; %纹理基元的纹理编码集个数
textons = createTextons(imStack, bank, k);
%height = 100 ;
%winSize =ones(height,height);
%colorLabelIm = ones(height, width) ;
%textureLabelIm =ones(height, width);
[colorLabelIm, textureLabelIm] = compareSegmentations(origIm, bank, textons, winSize, numColorRegions, numTextureRegions);
% disp(colorLabelIm);
[mc,nc,lc] =size(colorLabelIm);
[mt,nt,lt] =size(textureLabelIm);
fprintf('row is %d,col is %d,%d ',mc,nc,lc);
fprintf('row is %d,col is %d,%d ',mt,nt,lt);
subplot(1,3,2);
%imshow(uint8(colorLabelIm));% 输出为黑白图
x2=[150+width,350+2*width];
y2=[50+height,50+2*height];
imagesc(x2,y2,colorLabelIm)
colorbar
%colorbar('position',[0.32 0.01 0.3 0.5])
title('颜色分割');
subplot(1,3,3);
%imshow(uint8(textureLabelIm));% 输出为黑白图
x3=[250+2*width,450+3*width];
y3=[50+2*height,50+3*height];
imagesc(x3,y3,textureLabelIm)
%imagesc(textureLabelIm)
colorbar %定义图例位置大小
title('纹理分割');
%suptitle('图像名：gumballs.jpg 滤波器： filterBank 条件：winSize=%d nColor=%d nTexture=%d texttonsNum=%d ',winSize,numColorRegions,numTextureRegions,k);
suptitle('图像名：gumballs.jpg 滤波器： filterBank 条件：winSize= nColor= nTexture= texttonsNum= ')
% imshow(colorLabelIm,[]); % 输出为黑白图
