铅笔特效算法

铅笔特效方面的算法很多，像photoshop这种牛逼软件具有这个功能就不足为奇了。偶尔闲来无事，收集一些这方面的资料，觉得蛮有意思的。最近看到一篇香港中文大学贾佳亚教授的paper，当然他只是挂个名，因为他不是第一作者，应该是他学生的作品。废话不多说，这篇paper的题目是“Combining Sketch and Tone for Pencil Drawing Production”。看到这篇paper的实现效果，简直让人惊呆了。这是我见过的最好的铅笔特效算法，没有之一。文章的原理也很简单。希望感兴趣的朋友去下来研究研究。你值得拥有！

下面贴出自己的代码 ，甚为粗糙，还望看官莫要见笑。

clear all;
close all;
clc;
ks=10;
dirNum=8;
[FileName,PathName]=uigetfile({'*.jpg';'*.png';'*.bmp'},'Open an Image File');
Img=imread([PathName FileName]);
figure,imshow(Img);
[H, W, sc] = size(Img);
    if sc == 3
        I = rgb2gray(Img);
    end
PC = phasecongmono(I, 5, 3, 2.5, 0.55, 2.0);
figure,imshow(PC);
PC=PC.^0.65;

%% convolution kernel with horizontal direction 
    kerRef = zeros(ks*2+1);
    kerRef(ks+1,:) = 1;

%% classification 
    response = zeros(H,W,dirNum);
    Sresponse = zeros(H,W,dirNum);
    for ii = 0 : (dirNum-1)
        ker = imrotate(kerRef, ii*180/dirNum, 'bilinear', 'crop');
        response(:,:,ii+1) = conv2(PC, ker, 'same');
    end
    [~ , index] = sort(response, 3); 
    Dirs = (index(:,:,end)-1)*180/dirNum;
    maxindex=index(:,:,end);
    C = zeros(H,W,dirNum);
    for k=1:dirNum
    for i=1:H
        for j=1:W
            if(maxindex(i,j)==k)
                C(i,j,k)=PC(i,j);
                    
            end
        end
    end
    end
Cresponse = zeros(H,W,dirNum);
for jj = 0 : (dirNum-1)
        ker = imrotate(kerRef, jj*180/dirNum, 'bilinear', 'crop');

        Cresponse(:,:,jj+1) = conv2(C(:,:,jj+1), ker, 'same');
end
S=0.079*Cresponse(:,:,1)+0.079*Cresponse(:,:,2)+0.079*Cresponse(:,:,3)+0.079*Cresponse(:,:,4)+0.079*Cresponse(:,:,5)+0.079*Cresponse(:,:,6)+0.079*Cresponse(:,:,7)+0.079*Cresponse(:,:,8);
figure,imshow(1-S);
a=0:255;
p1 = 1/9*exp(-(255-a)/9);
p3=1/sqrt(2*pi*11) * exp(-(a-90).*(a-90) / (2.0*11*11));
p2= 0:255;
p2(:)=1/(225-105);
p2(1:55)=0;
p2(225:256)=0;
p = 0.76 * p1 + 0.22 * p2 + 0.02 * p3;
figure,plot(a,p);
hgram=round(p*H*W);

Inew=histeq(uint8(I),hgram);
H=fspecial('gaussian',13,2);
Inew = imfilter(Inew,H,'replicate');
Inew=double(Inew);
Sn=1-mapminmax(S);
out=Inew.*Sn*0.55;
figure,imshow(uint8(out));

需要提醒的是，本代码未能实现texture rendering的步骤，所以效果与文章还是有一定差距的。