用函式畫出可愛的卡通貓
函式共 268 頁 11783 行，前後折騰了近半個月。可惜還是有些小瑕疵，不然文章標題就該叫《震驚！死宅竟用三角函式畫出可愛的老婆》了......
用函式畫的 Pusheen
畫 Pusheen 的函式
前言
在推上看到這樣一張圖：
Twitter截圖
覺著用函式畫妹子什麼的好有趣，決定自己試著實現一下。
思路
從圖中可以看到，繪製妹子的函式為三角函式項的和。不難發現，函式以 t 為引數，兩個一組分別為正弦級數和餘弦級數，首項(a0/b0)為常數，很容易看出來，函式應該和傅立葉級數有關係。
大體思路如下：
- 預處理圖片。用 MATLAB 將想要用函式繪製的圖片讀入
- 提取輪廓
- 獲取輪廓的座標
- 處理座標。獲取到的座標為 N 點長的離散序列，將其進行 N 點離散傅立葉級數展開（DFS）
- 驗證。將輸出結果的函式繪製成影象，觀察效果
實現
先畫個簡單點的。桌面上有張可愛的 Pusheen 圖片，就用它啦。
萌萌的 Pusheen 貓
1. 預處理圖片
將 Pusheen 圖片讀入MATLAB。
cd C:\Users\Desktop
pic = imread('C:\Users\Ascend\Desktop\pusheen.jpg');轉換為二值影象。
pic = im2bw(pic);
show(pic)轉換為二值影象的原因是便於提取輪廓。看下效果。
二值影象效果
2. 提取輪廓
這裡直接呼叫了 MATLAB 自帶提取輪廓的演算法。
pic = edge(pic,'sobel');
imshow(pic);效果如下圖，看起來不錯，輪廓分明線條幹淨。
輪廓
3. 獲取輪廓座標
獲取輪廓座標，並將能夠圍成封閉曲線的座標分組，則每組為一個有限長的離散序列，再將每個序列進行 DFS 展開。
獲取輪廓座標並不難，難在擬合封閉曲線這裡。想了想毫無頭緒，上網搜尋找到了 MATLAB 自帶的 imcontour 函式。精度一般，使用條件苛刻，後文再敘，先湊合用。
point = imcontour(pic);得到一個 2x126494 的矩陣 point，裡面有我們需要的座標。
4. 處理座標
處理矩陣
座標已經得到，先試著用座標畫一下影象。
width=point(2,:);
height=point(1,:);
plot(width,height)卻得到了一副非常奇怪的影象。
一團亂麻
按道理，座標描點得到的影象就應該為 Pusheen 貓的輪廓。注意到影象有幾個異常點，橫座標特別大。檢視矩陣 point 的結構。
矩陣 point 的結構
第一行為縱座標 height ，第二行為橫座標 width，每一列為一個點 (height,width)。
第一個點的數值很奇怪。橫座標為 2435 ，縱座標為 0.1，怎麼看都不是輪廓上的點。也許是某種標記？於是向後翻了2435個點，看到第 2437 列：
2437
又是一個 (0.1,2427)。再向後翻 2427 個點：
4865
第 4865 列，又是一個 (0.1,2531)。於是這個點的意義就清楚了：橫座標為 0.1 的點為特殊標記，其值為該組點的個數。依此可以將這個 2x126494 的矩陣拆成多個小矩陣，每個矩陣為一組可以擬合成封閉曲線的座標(imcontour函式認為的)。
拆分矩陣
%判斷point(1,:)第一行 0.1 的個數，從而判斷length的長度
length=0;
for i=1:size(point,2)
if point(1,i)==0.1
length=length+1;
end
end
%開始分段，拆成a1，a2，a3，...，一共length個二維陣列
for i=1:length
temp=point(2,1);
point(:,1)=[];
eval(sprintf('a%d=point(1:2,1:temp);', i)); %將第i段資料送入第i個數組
point(:,1:temp)=[];
end
%縱座標取負數，將離散資料轉為座標
for i=1:length
eval(sprintf('a%d(2,:)=-a%d(2,:);',i,i));
end離散傅立葉級數展開
離散傅立葉級數的公式如下，根據公式，寫出 MATLAB 程式碼。
DFS公式
file = fopen('outer.txt', 'w'); %I/O方式輸出結果，結果輸出到桌面的outer.txt中
for count=1:length
eval(sprintf('a=a%d;',count));
N = size(a,2);
%a->ax(虛陣列)
are=a(1,:);
aim=a(2,:);
ax=are+aim*j;
clear are aim;
%A:傅立葉正變換陣列
A=[];
for i = 0:N-1
ang = -2 * 1j * pi * i / N;
r = 0;
for j=0:N-1
r = r + exp(ang * j) * ax(j+1);
end
A(i+1)=r;
end
clear ang i j r;
fprintf(file,['t=1:',num2str(N),';\n']);
fprintf(file,'x(t)=');
for i = 0:N-1
ang = 2 * pi * i / N;
fprintf(file,'+');
fprintf(file,[num2str(real(A(i+1) / N)),'*cos(', num2str(ang), '*t)-']);
fprintf(file,[num2str(imag(A(i+1) / N)),'*sin(', num2str(ang), '*t)']);
end
fprintf(file,';\n');
fprintf(file,'y(t)=');
for i = 0:N-1
ang = 2 * pi * i / N;
fprintf(file,[num2str(imag(A(i+1) / N)),'*cos(', num2str(ang), '*t)+']);
fprintf(file,[num2str(real(A(i+1) / N)),'*sin(', num2str(ang), '*t)+']);
end
fprintf(file,'0;\n');
fprintf(file,'plot(x,y);hold on;\n');
end
fclose(file);5. 驗證
輸出的函式在桌面的 outer.txt 檔案中。
結果
將函式輸入到 MATLAB 中畫影象，結果如下。
Pusheen
驗證無誤。輸出那一堆就是能畫出可愛的 Pusheen 的函式啦 ~
其他
卡通小貓都畫出來了，打算也畫一個卡通妹子。在網上找了張艾米莉亞線稿
EMT
使用了用同樣的方法，結果如下。
混亂的艾米莉亞
一臉茫然，大概輪廓倒能看出來，可怎麼會亂成這樣子[笑cry]
仔細想了想，原因可能出在閉合曲線分組那兒。卡通貓輪廓上的蜜汁線條也可能因為此。
對卡通貓的函式進行 Debug：
step1
step2
step3
step4
step5
這些就差不多能看出原因了。imcontour 函式將鼻子、眼睛也認為是外輪廓的一部分，錯誤地將其分為一組資料。對於像 Pusheen 這種線條簡單的圖形來說，精度一般，分錯組數較少。對於艾米莉亞這種複雜的人像來說，精度遠遠不夠，錯誤就太多了。
回到從推上看到的圖片。人家是如何畫出那麼完美的人像的呢？翻了翻博主的其他推文，發現人家利用遊戲引擎 Hot Soup Processor 設計了一個畫圖板程式，在程式上手繪一組組封閉曲線，遊戲引擎會用傅立葉級數將這組封閉曲線展開為函式。這個過程不需要用程式將離散序列分組，精度極高。
在沒找到更好的分組方法前，就只能畫畫輪廓簡單的圖片啦[笑cry]
一些其他的函式曲線：
Twitter Curve：
twitterline
twitterfunction
1:
one
onefunction