MATLAB 图片三角风格化(low poly)

该篇博 客内容已上传至Github 链接：https://github.com/slandarer/low_poly_poisson

三角化效果:

步骤 1.图片灰度化后进行sobel卷积检测边缘

oriPic=imread('test.jpg');

% use sobel algorithm to detect image edges
if size(oriPic,3)==3
    grayPic=rgb2gray(oriPic);
else
    grayPic=oriPic;
end
sobelPic=sobelConv2_gray(grayPic);

sobel边缘检测函数： 使用矩阵运算加快运行速度

	function sobelPic=sobelConv2_gray(oriPic)
        Hx=[-1 0 1;-2 0 2;-1 0 1];
        Hy=[1 2 1;0 0 0;-1 -2 -1];

        [rows,cols]=size(oriPic);
        exPic=uint8(zeros([rows+2,cols+2]));
        exPic(2:rows+1,2:cols+1)=oriPic;
        
        exPic(2:rows+1,1)=oriPic(:,1);
        exPic(2:rows+1,cols+2)=oriPic(:,cols);
        exPic(1,2:cols+1)=oriPic(1,:);
        exPic(rows+2,2:cols+1)=oriPic(rows,:);
        
        exPic(1,1)=oriPic(1,1);
        exPic(rows+2,1)=oriPic(rows,1);
        exPic(1,cols+2)=oriPic(1,cols);
        exPic(rows+2,cols+2)=oriPic(rows,cols);
        
        Gx=zeros([rows,cols]);Gy=Gx;
        
        for ii=1:3
            for jj=1:3
                tempPic=double(exPic(ii:rows+ii-1,jj:cols+jj-1));
                Gx=Gx+tempPic.*Hx(ii,jj);
                Gy=Gy+tempPic.*Hy(ii,jj);
            end
        end
        sobelPic=uint8(sqrt(Gx.^2+Gy.^2));
    end

原图： 卷积结果：

2.泊松圆盘采样

该算法的详解可以看看这一篇： http://devmag.org.za/2009/05/03/poisson-disk-sampling/ 使用泊松圆盘取样可以产生更加均匀的随机点： 左：随机取点 中：先构建网格，在网格交错点附近添加扰动 右：泊松圆盘取样 我们对泊松圆盘采样算法进行两次改写，分别应用于边缘处，及其他部分的采点，改写分别为：

1. 边缘采点所采集点位置只能是sobel边缘检测较亮位置
2. 其他部分采点密度不同时两点之间最小距离会随之改变，且其他部分采点以边缘点为初始集

edgePic=sobelPic;
edgePic(edgePic0);
edgePntList=[edgeY,edgeX];

% set the triangle density
redge=min(size(sobelPic))/80;
rmax=min(size(sobelPic))/20;
rmin=min(size(sobelPic))/40;

% use poisson disc sampling to select points
edgePntList=poissonEdge(edgePntList,redge);
pntList=poissonDisk(sobelPic,[rmin,rmax],30,edgePntList);

其中 redge ：边缘上两个取样点最近距离 rmax ：边缘外 亮度最低时两个采样点最近距离 rmin ：边缘外 亮度最高时两个采样点最近距离

当希望使边缘更加细致时，可以将redge后面的除数变大以获得更加密集的采样，例如：

redge=min(size(sobelPic))/120;

边缘采样函数

    function resultSet=poissonEdge(edgeList,R)
        preSet=edgeList;
        resultSet=[0 0];
        resultSet(1,:)=[];
        
        times=0;
        while times=R&dis30
                pntSet=edgeList(candidate(1:30),:);
            else
                pntSet=edgeList(candidate,:);
            end
            
            flag=0;
            for j=1:size(pntSet,1)
                pnt=pntSet(j,:);
                if size(resultSet,1)==0
                    resultSet=[resultSet;pnt];
                    preSet=[preSet;pnt];
                    flag=1;
                else
                    dis=sqrt(sum((resultSet-pnt).^2,2));
                    if all(dis>=R)
                        resultSet=[resultSet;pnt];
                        preSet=[preSet;pnt];
                        flag=1;
                    end
                end
            end
            if flag==1
                preSet(tempPos,:)=[];times=0;
            else
                times=times+1;
            end 
            disp(['edge pnt num:',num2str(size(resultSet,1))]);
        end
    end

其他区域采样函数

    function resultSet=poissonDisk(grayPic,R,K,edgePntList)
        [m,n]=size(grayPic);
        
        preSet=edgePntList;
        resultSet=[edgePntList;1,1;n,m;1,m;n,1];
        grayPic=double(255-grayPic);
        cmin=min(min(grayPic));
        cmax=max(max(grayPic));
        rMap=grayPic-cmin;
        rMap=rMap./(cmax-cmin).*(R(2)-R(1))+R(1);
        
        times=0;
        while times=1&&pnt(2)>=1&&pnt(1)=R)
                        resultSet=[resultSet;pnt];
                        preSet=[preSet;pnt];
                        flag=1;
                    end
                end
            end
            if flag==1
                preSet(tempPos,:)=[];times=0;
            else
                times=times+1;
            end 
            disp(['edge pnt num:',num2str(size(resultSet,1))]);
        end
    end


    function resultSet=poissonDisk(grayPic,R,K,edgePntList)
        [m,n]=size(grayPic);
        
        preSet=edgePntList;
        resultSet=[edgePntList;1,1;n,m;1,m;n,1];
        grayPic=double(255-grayPic);
        cmin=min(min(grayPic));
        cmax=max(max(grayPic));
        rMap=grayPic-cmin;
        rMap=rMap./(cmax-cmin).*(R(2)-R(1))+R(1);
        
        times=0;
        while times=1&&pnt(2)>=1&&pnt(1)