彩色图像生成素描

colorgrad.m

%梯度边缘提取算法
function [VG, A, PPG] = colorgrad(f, T)
if (ndims(f)~=3) || (size(f,3)~=3)              %ndims是matlab中求一个数组维数的函数，size(f,3)返回f的行数
    error('Input image must be RGB');
end
sh = fspecial('sobel');                          %sobel用于边缘提取
sv = sh';                                        %转置
Rx = imfilter(double(f(:,:,1)), sh, 'replicate');%线性空间滤波函数；replicate表示图像大小经过复制外边界的值来扩展，
Ry = imfilter(double(f(:,:,1)), sv, 'replicate');%symmetric表示图像大小经过沿自身的边界进行滤镜
Gx = imfilter(double(f(:,:,2)), sh, 'replicate');%circular表示经过将图像做为二位周期函数的一个周期来扩展
Gy = imfilter(double(f(:,:,2)), sv, 'replicate');
Bx = imfilter(double(f(:,:,3)), sh, 'replicate');
By = imfilter(double(f(:,:,3)), sv, 'replicate');
 
gxx = Rx.^2 + Gx.^2 + Bx.^2;
gyy = Ry.^2 + Gy.^2 + By.^2;
gxy = Rx.*Ry + Gx.*Gy + Bx.*By;
A = 0.5*(atan(2*gxy./(gxx-gyy+eps)));            
%%atan计算反正切函数                            |||||
%%./点除                                       |||||          若是a、b是矩阵，a./b就是a、b中对应的每一个元素相除，获得一个新的矩阵；若是a、b是两个数，那么a./b就是普通的除法
%%eps是一个函数，能够返回某一个数N的最小浮点数精度|||||

G1 = 0.5*((gxx+gyy) + (gxx-gyy).*cos(2*A) + 2*gxy.*sin(2*A));  %
A = A + pi/2;
G2 = 0.5*((gxx+gyy) + (gxx-gyy).*cos(2*A) + 2*gxy.*sin(2*A));
G1 = G1.^0.5;
G2 = G2.^0.5;
VG = mat2gray(max(G1, G2));            %mat2gray实现图像矩阵的归一化操做
 
RG = sqrt(Rx.^2 + Ry.^2);
GG = sqrt(Gx.^2 + Gy.^2);
BG = sqrt(Bx.^2 + By.^2);
 
PPG = mat2gray(RG + GG + BG);
 
if nargin ==2                           %是用来判断输入变量个数的函数
    VG = (VG>T).*VG;
    PPG = (PPG>T).*PPG;
end

---

demo.m

close all;clear all;clc;
f = imread('wo.jpg');
[VG,A,PPG] = colorgrad(f);
ppg = im2uint8(PPG); %把图像数据类型转换为无符号八位整型
ppgf = 255 - ppg;
[M,N] = size(ppgf);
T=150;
ppgf1 = zeros(M,N);
for ii = 1:M
    for jj = 1:N
        if ppgf(ii,jj)<T %边缘区域
%             ppgf1(ii,jj)=0;
            ppgf1(ii,jj)=(T-ppgf1(ii,jj))/5;
        else %平滑区域
            ppgf1(ii,jj)=235/(255-T)*(ppgf(ii,jj)-T);
        end
    end
end
ppgf1 = uint8(ppgf1);
figure;
subplot(221);imshow(ppgf);
subplot(222);imshow(ppgf1);
subplot(223);imhist(ppgf);
subplot(224);imhist(ppgf1);
 
figure;imshow(ppgf1);

实验原图:

实验效果：