Spojování snímků pořízených z jednoho místa

 

 

Předmět:              Počítačové vidění pro informatiku

Student:                              Michal Stránský, Jan Koutník

 

 

 

 

Teorie

 

Metoda spojování snímků vychází z modelu dírkové kamery, tedy ze středového promítání.  Obraz se promítne na plošku v kameře tak, že se to dá matematicky popsat takto: Jednotlivé body předmětů ve scéně mají světové souřadnice xw , v souřadnicích cw se nachází střed promítání kamery, a (u,v) jsou souřadnice promítnutého bodu v nějaké souřadné soustavě na plošce v kameře, pak se dá odvodit, že platí:

 

                              (r1)

 

kde K je kalibrační matice (vnitřní parametr kamery), R je rotační matice udávající jak je souřadný systém kamery natočen vzhledem k světovému. Rovnice r1 se dá upravit na tvar:

 

       (r2)

 

kde byl vektor xw rozšířen na (x,y,z,1)T, koeficient a zaveden místo zc .

Pořídíme-li se z jednoho místa snímky  ze dvou směrů a překrývají-li se záběry, můžeme najít transformaci mezi jednotlivými snímky. Promítací transformace v rovnici 2 bude pro oba záběry různá, ale souřadnice polohy středu zobrazení cw bude stejná a bez újmu na obecnosti se můžeme počátek světová souřadné soustavy položit do středu promítání, tj. cw=0 , tudíž je možno zmenšit matici o poslední sloupec na 3x3 a vektor xw v homogenních souřadnicích na obyčejný 3-vektor. Nazveme-li transformační matici před vektorem třeba písmeny H1 a H2 pro první a druhý snímek, a vyjádříme-li obě transformace zpětně, tj. pro vektor světové souřadné soustavy, můžeme pro jednotlivé body prostoru které mají obraz v obou snímcích napsat rovnost zpětných transformací (odpovídající obrazy bodů mají stejný vzor). Tato rovnost pak přímo vede k transformaci odpovídajících bodů z jednoho snímku do druhého:

 

      (r3)

 

kde H je transformační matice která se má jednoduchý vztah k předchozím uvedeným(H2-1H1), i je index společných bodů, ai je zavedený koeficient transformace, druhý index u souřadnic značí číslo snímku. Pro transformaci souřadnic z jednoho snímku do druhého stačí znát matici H, respektive libovolný její nenulový násobek, protože o korekci se postará koeficient ai . Matici H vypočítáme z minimálně čtyř bodů na každém snímku, které si sobě odpovídají (žádná trojice bodů by neměla ležet na přímce). V soustavě r3 je vidět, že z třetího řádku se vypočítá ai a použije se v prvních dvou řádcích, pro  každou dvojici odpovídajících bodů vzniknou dvě rovnice. Matice H má devět prvků, ale libovolný nenulový násobek je postačující, takže je potřeba zjistit 8 členů, tj. je potřeba 8 rovnic, které vzniknou právě ze čtyř dvojic odpovídajících bodů. V případě, že jsou zadány pouze čtyři dvojice pak k nalezení matice H stačí nalézt bázi jednodimenzionálního prostoru té soustavy osmi homogenních rovnic; řešení bude jednodimenzionální, pokud byly zvoleny čtveřice tak, že žádné tři body neležely na přímce. Rovnice v soustavě se vytvářejí tak, že pro každou dvojici odpovídajících bodů se přidají tyto rovnice:

 

          (r4)

 

kde ui,k a vi,k mají stejný význam jako výše, hjk jsou členy matice H. Pro více než čtyři body bude rovnice pro H přeurčená a matice rovnice bude mít hodnost 9 tudíž  bude mít pouze triviální řešení. Pro takový případ je nutno použít aparát umožňující snížit hodnost matice soustavy na 8 a přitom tak, že se od původní liší minimálně. Je možno použít metodu SVD která rozloží mxn matici nad reálnými čísly na součin matic UDVT  kde U je rotační matice mxm, V je rotační matice nxn a D je diagonální matice mxn s prvky na diagonále kladnými nebo nulovými a od levého horního rohu dolů seřazeny sestupně. Pro minimální odlišnost dvou matic se dá použít míra která, dvěma maticím přiřadí číslo, které je součet kvadrátů členů rozdílové matice.  Pro takto navrženou míru se dá ukázat, že minimální odlišnost matice původní od matice s hodností 8, se dá zajistit vynecháním nejmenšího diagonálního členu matice D. Ještě jednodušeji pak vyplyne řešení homogenní soustavy s takto modifikovanou maticí soustavy: báze jednodimenzionálního prostoru řešení je devátý sloupec matice V. Protože by takto navržená matice soustavy obsahovala veliká čísla a členy které mají v ní byt teoreticky nulové by se po aproximaci mohly dramaticky lišit od nuly, je nutné souřadnice ve snímcích normalizovat tak, že střední hodnoty souřadnic zadaných bodů jsou nulové a směrodatná odchylka v jednotlivých směrech je co nejblíže jednotková. Toto se dá u obou snímků zajistit transformací:

                       (r5)

 

kde:

 

 

n je počet dvojic odpovídajících bodů, u’i,k v’i,k jsou normalizované souřadnice i-té dvojice bodů v k-tém snímku. Označí-li se normalizační matice pro první resp. druhý snímek L1 resp. L2 , pak se vypočítá matice transformace normalizovaných souřadnic H’ z prvního do druhého snímku výše naznačeným způsobem  podle soustavy sestavené na základě r4 použijí-li se místo původních u a v normalizované u’ a v’. Výsledná transformační matice se spočítá jako:

 

H=L2-1H’L1 .      (r6)

 

V okamžiku kdy je k dispozici matice H se může podle vztahu r3 přepočítávat souřadnice z jednoho snímku do druhého. Při transformaci je třeba vydělit výsledné souřadnice u,v hodnotou třetí složky. Nastává problém, že rastr jednoho snímku přetransformovaný do druhého bude deformovaný a souřadnice budou padat mimo rastr obrazu, tj. je třeba převzorkování. Existují různé metody, v této práci byla použita metoda nejbližšího souseda, při které se souřadnice vzniklé po transformaci zaokrouhlí na nejbližší celou hodnotu a z těchto celočíselných souřadnic se vybere pixel.

 

 

Postup spojení

 

Celý proces spojení snímků postupuje podle následujících kroků:

  1. výběr korespondujících bodů
  2. vytvoření místa u prvního snímku pro vložení převzorkovaného druhého snímku s patřičnou transformací vybraných bodů.
  3. v případě více jak č dvojic bodů normalizace snímkových souřadnic v obou snímcích (r5)
  4. sestavení rovnic pro matici homografie H, pro každou dvojici bodů se přidají 2 řádky podle (r4)
  5. v případě čtyř bodů se nalezne matice H jako přerovnaný vektor báze řešení výše sestavené rovnice, v případě více bodů se  provede SVD dekompozice matice soustavy, a výsledná matice H se spočítá podle (r6), matice H’ je přerovnaný devátý sloupec matice V vzniklý v SVD.
  6. hranice druhého obrázku se převedou do souřadnic prvního za pomocí inverzní transformace H, vzniklé body jsou v prvním obrázku uzavřeny v nejmenším obdélníkem s celočíselnými souřadnicemi rovnoběžným se souřadnou soustavou prvního snímku.
  7. pro všechny celočíselné souřadnice tohoto obdélníku v prvním snímku se vygenerují transformací H souřadnice ve druhém snímku a vezmou se jejich nejbližší celočíselné hodnoty.
  8. vygeneruje se mapa říkající které z těchto souřadnic jsou uvnitř druhého snímku, pro vnější body se nebude provádět přenos pixelů v bodě 9
  9. první snímek se doplní o body z druhého tak, že se podle masky z bodu 8 přenesou body ze souřadnic druhého snímku vygenerovaných v bodě 7 do souřadnic v prvním snímku, které byly v bodě 7 před transformací.

 

 

Příklad aplikace

 

V následující obrázcích lze vidět částečně tento postup:

 

              

 

Obr. 1-a                                                                                         Obr. 1-b

 

Na obrázcích 1-a a 1-b byly vybrány dvojice korespondujících bodů (modré hvězdičky), v tomto případě jich bylo více než 4 tudíž byla použita metoda s SVD. Obrázek 1-a zde hraje roli snímku 2 a obr. 1-b roli snímku 1 ve výše uvedené teorii. Rastr obrázku 1-b byl rozšířen pro vložení přetransformovaného snímku 1-a. Byla vypočítaná matice H a podle bodů 6 až 9 byl vložen snímek 1-a do snímku 1-b. Výsledný obrázek se je uveden v obr. 2, pro názornost obrázek nebyl oříznut do estetičtějšího tvaru, čtenář si zajisté dovede představit jak obrázek vhodně oříznout.

 

Obr. 2:  Výsledný obrázek.

 

Pro aplikaci metody v této práci byl použit systém Matlab 5.3, zdrojový skript lze nalézt v souboru spojit.m, který je vám s tímto dokumentem k dispozici. Skript dále používá soubory amat.m který sestavuje matici soustavy a imshowb.m, tyto skripty jsou rovněž k dispozici.

 

 

load body3

im1=imread('panorama1.jpg');

im1b=uint32(double(im1(:,:,1))*2^16+double(im1(:,:,2))*2^8+double(im1(:,:,3)));

im2=imread('panorama2.jpg');

im2b=uint32(double(im2(:,:,1))*2^16+double(im2(:,:,2))*2^8+double(im2(:,:,3)));

 

u1=[b2(:,1)';b1(:,1)'];

u2=[b2(:,2)';b1(:,2)'];

 

figure(1);

hold off

imshow(im1);

hold

plot(u1(2,:),u1(1,:),'*');

figure(2);

hold off

imshow(im2);

hold

plot(u2(2,:),u2(1,:),'*');

 

vyska=size(im2,1);

sirka=size(im2,2);

 

imr=[zeros(size(im2b)) im1b];

imr=[zeros(vyska/2,sirka*2);imr;zeros(vyska/2,sirka*2) ];

 

uu1=u1+[vyska/2 ;sirka]*ones(1,size(u1,2));

 

figure(3);

hold off

imshowb(imr);

hold

plot(uu1(2,:),uu1(1,:),'*');

 

 

%normalizace

bb1=mean(uu1(1,:));

cc1=mean(uu1(2,:));

aa1=mean((uu1(1,:)-bb1).^2)^0.5;

bb2=mean(u2(1,:));

cc2=mean(u2(2,:));

aa2=mean((u2(1,:)-bb2).^2)^0.5;

L1=[1 0 -bb1;0 1 -cc1]/aa1;

L2=[1 0 -bb2;0 1 -cc2]/aa2;

uu1=L1*[uu1;ones(1,size(uu1,2))];

u2=L2*[u2;ones(1,size(u2,2))];

 

 

A=[];

for i=1:floor(size(u1,2)/4)

  A=[A;amat(uu1(:,(i-1)*4+1:i*4),u2(:,(i-1)*4+1:i*4))];

end  

%A=amat(uu1,u2);

%h=reshape(null(A),3,3)';

[U,S,V]=svd(A);

h=reshape(V(:,9),3,3)';

h=inv([L2;0 0 1])*h*[L1;0 0 1];

h1=inv(h);

 

prek=50;

rohy2=[1 1 1;vyska 1 1; 1 sirka-prek 1;vyska sirka-prek 1]';

rohy=h1*rohy2;

rohy=round(rohy(1:2,:)./([1;1]*rohy(3,:)));

ramec=[min(rohy') max(rohy')];

[y,x]=meshgrid(ramec(1):ramec(3),ramec(2):ramec(4));

indx=[y(:) x(:)]';

indx2=h*[indx;ones(1,size(indx,2))];

indx2=round(indx2(1:2,:)./([1;1]*indx2(3,:)));

uvnitr=(indx2(1,:)>=1)&(indx2(1,:)<=vyska)&(indx2(2,:)>=1)&(indx2(2,:)<=sirka);

indx=indx(1,:)+vyska*2*(indx(2,:)-1);

indx2=indx2(1,:)+vyska*(indx2(2,:)-1);

 

imr(indx(uvnitr))=im2b(indx2(uvnitr));

figure (4);

imshowb(imr);

 

Výpis skriptu spojit.m

 

              

function A=amat(u1,u2)

A=[];

for i=1:4

   x=[u2(:,i);1]*[u1(:,i)' 1];

   y=[x(3,:) 0 0 0 ; 0 0 0 x(3,:)];

   y=[y -x(1:2,:)];

   A=[A;y];

end

 

Výpis funkce amat.m

 

 

function imshowb(x)

y=uint8(zeros([size(x) 3]));

y(:,:,1)=uint8(double(x)/2^16);

y(:,:,2)=bitand( uint32(double(x)/2^8),255);

y(:,:,3)=bitand( uint32(double(x)),255);

 

imshow(y);

 

               Výpis funkce imshowb.m

 

 

Závěr

Výsledky spojení snímků se liší byla-li použita normalizace souřadnic vstupních bodů či nikoliv. Ovšem ani jedna z těchto variant nedává uspokojivý výsledek, což se dá odůvodnit buď posunem kamery nebo radiálním zkreslením optiky. V navázaném snímku se neprojevila závislost kvality spojení na vzdálenosti od kamery, tudíž přichází v úvahu radiální zkreslení.

 

Použité zdroje:

 

Přednášky předmětu Počítačové vidění pro informatiku, přednášející: Tomáš Pajdla, ČVUT FEL, 1999-2000.