MATLABでマウスインタラクティブGUI ～ 画像を馴染ませてコピー ～ Poisson Image Editing
まえがき
MATLABは、ちょっとしたマウス・インタラクティブなGUIも簡単に作れます。（要 Image Processing Toolbox）
今回は「Poisson Image Editing」を題材に、その使い方を解説してみようかと思います。
実装する機能としては以下の通りです。
2つの画像、src（コピー元）とdst（コピー先）を表示
src 側で範囲指定した画像を、dst 側に「馴染ませて」表示
src 側は位置とサイズ、dst 側は表示位置をマウスで選択可能
dst 側のマウスクリックで、「単なるコピペ／馴染ませる」を切り替え
操作例はこちらをご覧ください。
Poisson Image Editing
最近また、「Computational photography」という言葉を聞くようになりました。
「見た目そのまま」ではなく、コンピューターの助けを借りた、人の手ではできないような「レタッチ」処理技術全般を指す言葉と理解しています。
私が最初にその言葉を目にしたのは、多分2000年初頭の SIGGRAPH 。
あの頃は SIGGRAPH に限らず、SLAM、Light Field Photography、Coded Aperture 等々、面白い「Computational photography」技術が続々と発表されていました。
「Poisson Image Editing」もその中の一つ、Microsoft UKが2003年に発表した論文です。（前年にHDR応用も出していますが）
その処理内容から、「Gradient-domain image processing」とも呼ばれます。
歴史的には、2002年、Adobe Photshop7 に healing Brush なるものが実装されました。
それまでのコピーペーストではなく、周りと馴染むようにコピー元を修正して貼り付けることができる技術です。
詳細は不明ですが、Adobeは「熱力学応用の、（MSの）より高次の処理」みたいなことを言っていたと思います。
おそらく基本的には今回のと似たような、偏微分方程式を解く技術だとは思われますが。
このような機能は Photoshop だけでなく色々なアプリに実装されているので、最近の人は特に不思議にも思わないかもしれませんが、当時は衝撃的でした。
考え方自体は比較的シンプルで、src 画像のエッジ情報を保存したまま、dst 画像との境界で双方の輝度差が 0 となるよう、最小変更の src 修正画像を求める、というものです。
P. Pérez, M. Gangnet, A. Blake. Poisson Image Editing.
ACM Transactions on Graphics (SIGGRAPH'03), 22(3):313-318, 2003.
数学的には、Dirichlet境界条件でポアソン方程式を解くことになります。
最終的に Ax=b の x を求めることになるのですがここで問題が・・。
画像なので上式のそれぞれが巨大な行列となりますし、普通にやろうとすると簡単には解けません。
スパース行列
「巨大な行列」自体（とその基本演算）に対しては、スパース行列を使います。スパース行列はゼロ要素を取り除いた表現形式で、メモリの節約と高速化につながります。
詳しくはMATLABのドキュメント 「スパース行列の作成」 をご覧ください。
スパース ラプラス演算子行列も使うのですが、これは「クロネッカーのテンソル積」で求めることができます。
MATLABのドキュメント にも書いてあるそのままです。
「巨大な行列の逆行列」に対しては、DSTを使った高速演算法等も存在しますが、人類にはMATLABがあります。
そう、x = A\b で解けてしまいます！
しかも、\ は逆行列そのものを求める手法をとっていないので、そこそこ高速に求まります（ガウスの消去法を使っているようです）。
はい、解決。(^-^)
インタラクティブ ROI
さて、今回のテーマはここからが本題です。
画面上にインタラクティブなオブジェクト（ROI）を作成するには、drawrectangle() を使います。
他にも、フリーハンド、円、多角形等、様々な形状の ROI を作成することができます。こちらもMATLABのドキュメント「ROI 作成の概要」をご参照ください。
この ROI オブジェクトに addlistner() によってイベントリスナー（とその callback ）を作成することによって、マウス操作に対応した処理が行えます。
roiS = drawrectangle(ax1,'Position', [136 180 127 127],'FixedAspectRatio',true);
roiS.Label = 'src';addlistener(roiS,'ROIMoved',@(src,evt)allevents(src,evt,dstI,dstROI,srcImg,srcROI,srcI,dstImgAx));上記の様に ROI にラベルを付けておくと、callback 内でイベント別以外にそのラベル別でも処理を行うことができます。
動作例
まずは未処理状態の画像を。
左： src（コピー元） 右：dst（コピー先）
左が src（コピー元）で、右が dst（コピー先）です。
ここにそれぞれの ROI を設定し、src 側の ROI を dst 画像に馴染ませる処理をしてペーストします。
src の矩形はマウスで位置とサイズを変えられます（処理簡素化のためアスペクト比は 1:1 固定にしています）。
位置やサイズをマウスで変えてクリックを離すと新たな src 画像での処理が行われます。
dst 側の circle マーカーを移動させるとそこを中心にペーストされます。
セレクト状態（「Ctrl+左クリック」で false、「クリック」で true）も見られますので今回は、dst 側の circle を「Ctrl+左クリック」で「単なるコピペ」、「クリック」で「馴染ませる」に切り替えています。
（マーカー色が赤に変わる）
dst 側 circleを「右クリック -> 円の削除」 で、もう dst 側は動かせなくなりますがマーカーを消すこともできます。
MATLABコード
サンプル画像付きMATLABコード
ソースコード全体
function GDIP
% Poisson image editing, also called Gradient domain image processing
clear srcROI_old dstOI_old
% set source image
srcI = imread('ayaka18_007.jpg');
f1 = figure(1);
f1.MenuBar='none';
f1.Pointer = 'hand';
f1.InnerPosition = [10 20 size(srcI,2)*2 size(srcI,1)];
f1.NumberTitle = 'off';
f1.Name = 'Poisson image editing';
ax1 = subplot(1,2,1);
ax1.Units = 'Normalize';
ax1.Position = [0 0 0.5 1];
imshow(srcI, Parent=ax1);
% draw rectangle
roiS = drawrectangle(ax1,'Position', [136 180 127 127],'FixedAspectRatio',true);
roiS.Label = 'src';
roiS.LabelVisible = 'hover';
% get ROI
srcROI = roiS.Position;
srcROI(3) = srcROI(3);
srcROI(4) = srcROI(4);
srcImg = imcrop(srcI,srcROI);
% set destination image
dstI = imread('ayaka18_004.jpg');
ax2 = subplot(1,2,2);
ax2.Units = 'Normalize';
ax2.Position = [0.5 0 0.5 1];
dstImgAx = imshow(dstI,Parent=ax2);
% draw circle
roiD = drawcircle(ax2,'Center',[size(dstI,2)/2, size(dstI,1)/2],'Radius',10, ...
'InteractionsAllowed','translate', 'Color','r','SelectedColor' ,'g',...
'MarkerSize',1);
roiD.Label = 'dst';
roiD.LabelVisible = 'off';
% get ROI center
cnt = roiD.Center;
dstROI(1) = cnt(1) - fix(srcROI(3)/2);
dstROI(2) = cnt(2) - fix(srcROI(4)/2);
dstROI(3) = srcROI(3);
dstROI(4) = srcROI(4);
% copy paste
pasteI = dstI;
pasteI(dstROI(2):dstROI(2)+dstROI(4), dstROI(1):dstROI(1)+dstROI(3),:) = srcImg;
dstImgAx.CData = pasteI; % updata data only
% set listeners (source position/size, destination position, GDIE ON/OFF)
addlistener(roiS,'ROIMoved',@(src,evt)allevents(src,evt,dstI,dstROI,srcImg,srcROI,srcI,dstImgAx));
addlistener(roiD,'ROIMoved',@(src,evt)allevents(src,evt,dstI,dstROI,srcImg,srcROI,srcI,dstImgAx));
addlistener(roiD,'ROIClicked',@(src,evt)allevents(src,evt,dstI,dstROI,srcImg,srcROI,srcI,dstImgAx));
end
% listner callback
function allevents(src,evt,dstI,dstROI,srcImg,srcROI,srcI,ax2)
persistent srcROI_old dstROI_old gdie_on
if isempty(srcROI_old)
srcROI_old = srcROI;
dstROI_old = dstROI;
gdie_on = false;
end
srcLabel = src.Label;
evname = evt.EventName;
switch(srcLabel)
case 'src'
% disp('src moved')
p = fix(evt.CurrentPosition);
srcROI = p;
srcROI(3) = srcROI(3);
srcROI(4) = srcROI(4);
srcROI_old = srcROI;
srcImg = imcrop(srcI,srcROI);
dstROI = dstROI_old; % recall
dstROI(3) = srcROI(3);
dstROI(4) = srcROI(4);
dstROI_old = dstROI; % update
case 'dst'
% disp('dst moved')
if strcmp(evname,'ROIClicked')
srcROI = srcROI_old; % recall
srcImg = imcrop(srcI,srcROI);
dstROI = dstROI_old; % recall
gdie_on = evt.CurrentSelected;
else % 'ROIMoved'
p = fix(evt.CurrentCenter);
srcROI = srcROI_old; % recall
srcImg = imcrop(srcI,srcROI);
dstROI(1) = p(1) - fix(srcROI(3)/2);
dstROI(2) = p(2) - fix(srcROI(4)/2);
dstROI(3) = srcROI(3);
dstROI(4) = srcROI(4);
dstROI_old = dstROI;
end
end
% check over the border (only right/bottom)
[maxy, maxx, ~] = size(dstI);
y1 = dstROI(2);
x1 = dstROI(1);
y2 = dstROI(2) + dstROI(4);
ov = y2 - maxy;
if ov > 0
y1 = y1 - ov; % move to top
y2 = y1 + dstROI(4);
dstROI(2) = y1;
end
x2 = dstROI(1) + dstROI(3);
ov = x2 - maxx;
if ov > 0
x1 = x1 - ov; % move to left
x2 = x1 + dstROI(3);
dstROI(1) = x1;
end
if gdie_on
dstImg = imcrop(dstI,dstROI);
cmpImg = zeros(size(srcImg));
cmpImg(:,:,1) = gdie(srcImg(:,:,1),dstImg(:,:,1));
cmpImg(:,:,2) = gdie(srcImg(:,:,2),dstImg(:,:,2));
cmpImg(:,:,3) = gdie(srcImg(:,:,3),dstImg(:,:,3));
cmpImg = max(cmpImg,0);
cmpImg = min(cmpImg,1);
pasteI = dstI;
pasteI(y1:y2, x1:x2,:) = uint8(cmpImg*255);
else % copy paste
pasteI = dstI;
pasteI(y1:y2, x1:x2,:) = srcImg;
end
ax2.CData = pasteI; % updata data only keeping listner
end
% Poisson image editing (Gradient domain image processing)
function outI = gdie(src, dst)
src = im2double(src);
dst = im2double(dst);
r = length(src); % take max
c = r;
% Find the Laplacian matrix using Kronecker's tensor product
I = speye(c,r);
E = sparse(2:c,1:r-1,1,c,r);
D = E + E' - 2*I;
A = -kron(D,I)-kron(I,D);
dOmg = zeros(r,c);
dOmg(:,1) = 1; dOmg(:,c) = 1;
dOmg(1,:) = 1; dOmg(r,:) = 1;
dOmg = logical(dOmg);
% Poisson solution -> Laplace equation with Dirichlet boundary conditions
Lap = [0 -1 0; -1 4 -1; 0 -1 0];
srcGradient = conv2(src,Lap,'same');
divf = srcGradient;
dst_img = dst;
divf = divf(:);
d2f = dst_img(dOmg);
d2A = A(:,dOmg);
g = d2A * d2f;
d1b = divf - g;
d1A = A(~dOmg,~dOmg);
d1b = d1b(~dOmg,:);
d1x = d1A\d1b;
outI = zeros(size(divf,1),1);
outI(~dOmg,1) = d1x;
outI(dOmg,1) = d2f;
outI = reshape(outI,r,c);
endその他 Tips
・src/dst 2つの画像を一つの figure に隙間なく表示し、別々に addlistener するために以下のような設定にしています。
ax1 = subplot(1,2,1);
ax1.Units = 'Normalize';
ax1.Position = [0 0 0.5 1];ax2 = subplot(1,2,2);
ax2.Units = 'Normalize';
ax2.Position = [0.5 0 0.5 1];tiledlayout は設定できるプロパティが少ないので axes で・・。
もっとスマートなやり方があれば教えてください！
・新たな画像を imshow() してしまうと listener も消えてしまうので、データのみ （ CData ）を書き換えています。
ax2.CData = pasteI;
・静止画のデータ型はuint8（[0, 255]）です。
このような演算を行う場合、最初にim2double()（[0, 1]へのスケーリングもしてくれます）したくなるところですが、それをやると処理が遅くなるので、メイン処理まではuint8のままにしています。
・imcrop(x, y, w, h) の結果は、(w+1, h+1) のサイズになります。
・イメージ系関数の引数は (x, y, w, h) ですが、行列して扱う場合は r（行）/ c（列）の順番になるので (y, x, h, w) で指定します。
・先にソース画像全体の gradient を求めてしまっても良いのですが、それ自体は軽いのであまり速くはならないと思います。
あとがき
境界がエッジをまたいでいたり、src/dst 画像の組み合わせによっては馴染まないことも多々ありますので、色々試してみてください。今回の主題はこちらではないので、とりあえず深掘りはやめました。(^-^;
Poisson Image Editing には、これ以外にもHDR画像生成や色情報を保存した白黒化等、色々応用先も広い技術ですのでご興味のある方は掘ってみてください。
それでは、素敵なインタラクティブ・ライフを！
モデル：綾夏
この記事が気に入ったらサポートをしてみませんか？