Visão Computacional AlinhamentodeImagens · SergioL.Netto ProgramadeEngenhariaElétrica-COPPE/UFRJ Laboratório de Sinais, Multimídia e Telecomunicações [email protected] Julhode2017

Visão Computacional

Alinhamento de ImagensEduardo A. B. da Silva

Programa de Engenharia Elétrica - COPPE/UFRJLaboratório de Sinais, Multimídia e Telecomunicações

[email protected]

Sergio L. NettoPrograma de Engenharia Elétrica - COPPE/UFRJ

Laboratório de Sinais, Multimídia e Telecomunicaçõ[email protected]

Julho de 2017

(SMT – COPPE/UFRJ) UFRJ Julho de 2017 1 / 90

Sumário

Sumário

1 Alinhamento de ImagensTransformações 2DTranslaçãoTransformações AfinsHomografiasRANSACIPC


Alinhamento de Imagens Transformações 2D



Transformações 2D

Dado um conjunto de pontos casados em duas imagens A e B, como calcular atranformação T de A para B?



Relembrando: Transformações




Alinhamento de Imagens Translação

Caso simples: Translação

⇒ Como calcular (xt , yt)?



Deslocamento do i-ésimo ponto = (x ′i − xi , y ′i − yi)

(xt , yt) =(1n

n∑i=1

x ′i − xi ,1n

n∑i=1

y ′i − yi

)



De outra maneira:

xi + xt = x ′iyi + yt = y ′i

Sistema de equações lineares→ Quantas variáveis?

→ Quantas equações (por casamento)?



Mais equações do que incógnitas!→ Encontre a solução dos mínimos quadrados.



Mínimos Quadrados:

Para cada ponto (xi , yi)

xi + xt = x ′iyi + yt = y ′i

defina os resíduos como

rxi (xt) = (xi + xt)− x ′iryi (yt) = (yi + yt)− y ′i

Objetivo: minimizar a soma quadrática dos resíduos:

C(xt , yt) =n∑

i=1(rxi (xt)2 + ryi (yt)2)



Na forma matricial:

1 00 11 00 1... ...1 00 1

[xtyt

]=

x ′1 − x1y ′1 − y1x ′2 − x2y ′2 − y2

...x ′n − xny ′n − yn

2n×2A

2×1t =

2n×1b



Mínimos quadrados:

At = b

• Encontre t que minimize

||At − b||2

• Solução:

AT At = AT b

t = (AT A)−1AT b



Mínimos quadrados como modelagem de dadosEx: regressão linear



Regressão Linear

C(m, b) =n∑

i=1|yi − (mxi + b)|2



Regressão Linear x1 1x2 1... ...

xn 1

[mb

]=

y1y2...

yn


Alinhamento de Imagens Transformações Afins

Caso um pouco mais geral: Transformação Afim

Se conheço correspondências, posso encontrar o melhor modelo paraa transformação que mapeia essas correspondências:



Em coordenadas homogêneas:x ′y ′1

=

a b cd e f0 0 1

xy1

⇒ Quantas incógnitas?

⇒ Quantas equações por par de características?

⇒ Quantos pares de características precisamos?



Transformação Afim

• Resíduos:

rxi (a, b, c , d , e, f ) = (axi + byi + c)− x ′i

ryi (a, b, c , d , e, f ) = (dxi + eyi + f )− y ′i

• Custo:

C(a, b, c , d , e, f ) =n∑

i=1(rxi (a, b, c , d , e, f )2 + ryi (a, b, c , d , e, f )2)



Na forma matricial:

x1 y1 1 0 0 00 0 0 x1 y1 1x2 y2 1 0 0 00 0 0 x2 y2 1... ... ... ... ... ...

xn yn 1 0 0 00 0 0 xn yn 1

abcdef

=

x ′1y ′1x ′2y ′2...x ′ny ′n

2n×6A

6×1t =

2n×1b


Alinhamento de Imagens Homografias

HomografiasTransformação entre duas vistas de uma mesma superfície plana;

Ou entre imagens de duas câmeras que compartilham o mesmo centro;



Para retificar uma imagem:• Descubra a homografia H dados p e p′;

• Resolva equações da forma wp′ = Hp

� Incógnitas: w e coeficientes de H

� Homografias são definidas a menos de um fator de escala;

� Quantos pontos são necessários?



Solução para Homografiasx ′iy ′i1

=

h00 h01 h02h10 h11 h12h20 h21 h22

xiyi1

É necessário determinar 8 parâmetros, os 9 da matriz menos um correspondente aoescalamento que é irrelevante para coordenadas homogêneas:

x ′i = h00xi + h01yi + h02

h20xi + h21yi + h22

y ′i = h10xi + h11yi + h12

h20xi + h21yi + h22

x ′i (h20xi + h21yi + h22) = h00xi + h01yi + h02

y ′i (h20xi + h21yi + h22) = h10xi + h11yi + h12



Solução para Homografias:

x ′i (h20xi + h21yi + h22) = h00xi + h01yi + h02

y ′i (h20xi + h21yi + h22) = h10xi + h11yi + h12

[xi yi 1 0 0 0 −x ′i xi −x ′i yi −x ′i0 0 0 xi yi 1 −y ′i xi −y ′i yi −y ′i

]

h00h01h02h10h11h12h20h21h22

=[00

]



Solução para Homografias:

x1 y1 1 0 0 0 −x ′1x1 −x ′1y1 −x ′10 0 0 x1 y1 1 −y ′1x1 −y ′1y1 −y ′1... ... ... ... ... ... ... ... ...xn yn 1 0 0 0 −x ′nxn −x ′nyn −x ′n0 0 0 xn yn 1 −y ′nxn −y ′nyn −y ′n

h00h01h02h10h11h12h20h21h22

=

00...00

2n×9A

9×1h =

2n×10

Defina um problema de mínimos quadrados: minimize ||Ah − 0||2

• Como h só é definida a menos de um fator de escala, resolva para o vetor unitário h;

• Solução: h = autovetor de AT A com menos autovalor;

• Precisa de 4 ou mais pontos;(SMT – COPPE/UFRJ) UFRJ Julho de 2017 25 / 90


Exemplo: orientação de robô A geografia espacial da estrada é encontrada através da

retificação das formas pintadas pela homografia apropriada;



Algoritmo DLT (Direct Linear Transformation):Algoritmo linear para determinar H dado um conjunto de quatro correspondências de pontos 2Dpara 2D, xi ↔ x ′i . A transformação é dada pela equação x ′i = Hxi, ou x ′i × Hxi = 0.

Se a j-ésima linha da matriz H é denotada por hjT , então

Hxi =

h1T xih2T xih3T xi

Escrevendo x ′i = (x′i, y′i, w′i)T , o produto vetorial pode ser explicitado como

x ′i × Hxi =

y′ih3T xi − w′ih2T xi

w′ih1T xi − x′ih

3T xix′ih

2T xi − y′ih1T xi



Como hjT xi = xTi h

j para j = 1, ..., 3,

0T −w′ixTi y′ixT

i

w′ixTi 0T −x′ixT

i

−y′ixTi x′ixT

i 0T

h1

h2

h3

= 0

Estas equações são da forma Aih = 0, onde Ai é uma matriz 3× 9 e h é um vetor de tamanho 9composto pelas entradas da matriz H.

h =

h1

h2

h3

, H =

h1 h2 h3h4 h5 h6h7 h8 h9

Com hi sendo o i-ésimo elemento de h.



Como apenas duas das linhas são linearmente independentes,

[0T −w′ixT

i y′ixTi


i

]h1

h2

h3

= 0

Então, Aih = 0 é uma equação linear com h desconhecido, sendo Ai agora uma matriz 2× 9.

(x′i, y′i, w′i)T representam as coordenadas homogêneas do ponto x ′i . Escolhendo-se então w′i = 1,

(x′i, y′i) são as coordenadas medidas na imagem.

Se houver mais de 4 pares de pontos (solução super-determinada): não existe solução paraAh = 0 além da solução zero (devido a ruído). Deve-se, então, minimizar a norma ||Ah||.

Usar solução de mínimos quadrados para minimização.

É necessário que, dos 4 pontos em cada imagem, não haja 3 colineares. Caso contrário, não háuma solução única (configuração degenerada).



Solução de Mínimos Quadrados:Encontrar x que minimize ||Ax || sujeito a ||x || = 1.

Sendo A = UDV T , minimizar ||UDV T x ||. Entretanto, ||UDV T x || = ||DV T x || e ||x || = ||V T x ||.

Deve-se, então, minimizar ||DV T x || sujeito à condição ||V T x || = 1.

Escreve-se y = V T x e o problema se torna minimizar ||Dy || sujeito a ||y = 1||.

D é uma matriz diagonal com suas entradas em ordem decrescente.

A solução para este problema é y = (0, 0, ..., 0, 1)T , tendo uma entrada não zero, 1 na últimaposição.

Finalmente, x = Vy é simplesmente a última coluna de V .





SVD (Singular Value Decomposition):Dada uma matriz quadrada A, SVD é uma fatorizaçào de A tal que A = UDV T , onde U e V sãomatrizes ortogonais e D é uma matiz diagonal com valores não negativos. A decomposição podeser feita de modo que os valores na diagonal de D estejam em ordem decrescente.

SVD também existe para matrizes A não quadradas. Sendo A uma matriz m × n com m ≥ n, Ué uma matriz m × n com colunas ortogonais, A é uma matriz diagonal n × n e V é uma matrizortogonal n × n.

Como U possui colunas ortogonais, UT U = In×n. Além disso, U apresenta a seguintepropriedade: ||Ux || = ||x ||.

Em geral, UUT não é a identidade, a menos que m = n.

Para m < n, extende-se A adicionando-se linhas de zeros até tornar a matriz quadrada.

Os valores da diagonal de D são os valores singulares de A.



SVD (Singular Value Decomposition):Relação entre valores singulares e autovalores da matriz:

A = UDV T

AT A = VDUT UDV T = VD2V T

Como V é ortogonal, V T = V−1. Então

AT A = VD2V−1,

que é a equação que define os autovalores, indicando que os valores de D2 são os autovalores deAT A e que as colunas de V são os autovetores de AT A.

Os valores singulares de A, então, são as raízes quadradas dos autovalores de AT A.

Nota-se que AT A é simétrica e positiva semidefinida, então seus autovalores são reais e nãonegativos. Consequentemente, valores singulares são reais e não negativos.





Funções de CustoDistância Algébrica:

O algoritmo DLT minimiza a norma ||Ah||. O vetor ε = Ah é chamado de vetor residual, e é anorma deste vetor que é minimizada.

O vetor εi é o vetor de erro algébrico associado às correspondências de pontos xi ↔ xi′ e àmatriz H. A norma deste vetor é escalar e é chamada de distância alébrica:

dalg (x ′i ,Hxi)2 = ||εi||2 =∣∣∣∣∣∣∣∣[ 0T −w′ixT

i y′ixTi


i

]∣∣∣∣∣∣∣∣2Para quaisquer dois vetores x1 e x2:

dalg (x1, x2)2 = a21 + a2

2, onde a = (a1, a2, a3)T = x1 × x2

Dado um conjunto de correspondências, o valor ε = Ah é o vetor de erro algébrico para oconjunto completo:

∑i

dalg (x ′i ,Hxi)2 =∑i

||εi||2 = ||Ah||2 = ||ε||2



Distância Geométrica:

O vetor x representa as coordenadas na imagem medidas; x representa os valores estimados dospontos; e x os verdadeiros valores dos pontos.

Erro numa única imagem:

d(x , y) representa a distância Euclidiana entre os pontos não homogêneos representados por x ey . Então o erro de transferência para o conjunto de correspondências é

∑i

d(x ′i ,Hxi)2

A homografia estimada H é a que minimiza o erro acima.

Erro de transferência simétrico:

∑i

d(xi,H−1x ′i)2 + d(x ′i ,Hxi)2

Novamente, a homografia estimada H é a que minimiza o erro acima.



Erro de Reprojeção - ambas as imagens:

Objetivo: encontrar a homografia H e pares de pontos perfeitamente casados xi e x ′i a função deerro total

∑i

d(xi, xi)2 + d(x ′i , x ′i)2 sujeita a x ′i = Hxi ∀i

Minimizar essa função de custo envolve determinar tanto H quanto um conjunto decorrespondências subsidiárias {xi} e {x ′i}.

Esse modelo estima, por exemplo, a situação em que as correspondências medidas xi ↔ x ′i vêmde imagens de pontos num plano do mundo.

Desejamos estimar um ponto Xi de xi ↔ x ′i no plano do mundo que é então reprojetado paraestimar perfeitamente a correspondência xi ↔ x ′i .





Erro de Sampson:

Minimizá-lo requer a estimação simultânea tanto da matris de homografia quanto dos pontos xie x ′i (problema de estimação não linear).

O vetor X que minimiza o erro geométrico ||X − X ||2 é o ponto mais próximo na variedade VH àmedida X . Este ponto não pode ser extimado diretamente, exceto via iteração, por conta danatureza não linear da variedade VH . A ideia, então, é estimar a aproximação de primeira ordempara o ponto X assumindo que a função de custo é bem aproximada linearmente na vizinhançado ponto estimado.



A estimação da homografia entre dois planos pode ser pensada como o ajuste de uma superfíciede pontos em um espaço 4D, R4. Cada par de pontos de imagem x e x ′ define um único pontoX em uma medida no espaço R4, formado pela concatenação das coordenadas não homogêneasde x e x ′.

Para uma dada homografia H específica, as correspondências entre imagens x ↔ x ′ quesatisfazem x ′ × Hx = 0 definem uma variedade algébrica VH em R4 que é a interseção de duashipersuperfícies quádricas.

No caso, uma variedade é o conjunto simultâneo zero de um ou mais polinômios multivariadosdefinidos em RN .

Dados os pontos Xi = (xi, yi, x′i, y′i)T em R4, estimar a homografia é encontrar a variedade VHque aproximadamente passa pelos pontos Xi.

Para cada ponto Xi, Xi = (xi, yi, xi′, yi′)T é o ponto mais próximo de Xi que se encontra navariedade VH . Então:

||Xi − Xi||2 = (xi − xi)2 + (yi − yi)2 + (xi′ − xi′)2 + (yi′ − yi′)2 = d(xi, xi)2 + d(x′i, xi′)2



Para uma dada homografia H, qualquer ponto X = (x, y, x′, y′)T que se encontre em VH satisfazAh = 0.

Para enfatizar a dependência de X , escreveremos CH(X ) = 0, onde CH(X ) é um vetorbidimensional.

A aproximação de primeira ordem da função de custo por expansão em série de Taylor é

CH(X + δx) = CH(X ) + ∂CH (X)∂X δx

Se escrevermos δx = X − X e desejarmos que X esteja na variedade VH de modo queCH(X ) = 0, então o resultado é CH(X ) + ∂CH (X)

∂X δx = 0, que passará a ser escrito comoJδx = −ε, onde J é a matriz de derivadas parciais e ε é o custo CH(X ) assoxiado a X . Oproblema de minimização então se torna encontrar o menor δx que satisfaça a equação.



Encontrar o vetor δx que minimiza ||δx || sujeito a Jδx = −ε.

A forma padrão de se resolver esse tipo de problema é se usar multiplicadores de Lagrange:introduz-se um vetor λ de multiplicadores de Lagrange, e o problema se reduz a encontrar aextrema de δx

T δx − 2λT (Jδx + ε).

Derivando-se com respeito a δx e igualando-se a zero, tem-se:

2δxT − 2λT J = 0T

de onde se obtém que δx = JT λ. A derivada com respeito a λ dá Jδx + ε = 0, que é a restriçãooriginal. Substituindo para δx tem-se

JJT λ = −ε

que pode ser resolvido para λ dando λ = −(JJT )−1ε. Então, finalmente,

δx = −JT (JJT )−1ε

e X = X + δx . A norma ||δx ||2 é o erro de Sampson:

||δx ||2 = δxT δx = εT (JJT )−1ε



Funções de Custo Estatísticas e Estimação de MáximaVerossimilhança (Maximum Likelihood Estimation):Na ausência de erros de medida: xi′ = Hxi

Assume-se erros de medição obedecendo a uma distribuição Gaussiana isotrópica de média zero:x = x+∆x, com ∆x obedecendo a uma distribuição Gaussiana de variância σ2.

Assumindo-se que o erro de cada medida é independente, a função densidade de probabilidade(PDF) de cada ponto x medido é

Pr(x) =( 1

2πσ2

)e−d(x,x)2

2σ2



Erro numa única imagem:

A probabilidade de se obter o conjunto de correspondências {xi ↔ x ′i} é o produto de suas PDFsindividuais, já que assume-se que os erros de cada ponto são independentes:

Pr({x ′i |H}) =∏i

(1

2πσ2

)e−d(x′

i,Hxi)2

2σ2

Pr({x ′i |H}) é a probabilidade de se obter as medidas {x ′i} dado que a verdadeira homografia é H.



A verossimilhança logarítmica do conjunto de correspondências é

logPr({x ′i |H}) = − 12σ2

∑i

d(x ′i ,Hxi)2 + constante

A estimativa de máxima verossimilhança (MLE) da homografia, H, maximiza essaverossimilhança logarítmica, isto é, minimiza∑

i

d(x ′i ,Hxi)2

. Então, nota-se que a MLE equivale a minimizar a função de erro geométrico.



Erro em ambas as imagens:

Se as verdadeiras correspondências são {xi ↔ Hxi = xi′}, a PDF dos dados com ruído é

Pr({xi, x ′i}|H, {xi}) =∏i

(1

2πσ2

)e−(d(xi,xi)2+d(x′

i,Hxi)2)

2σ2

Deve-se buscar medidas na imagem "corrigidas"que se comportariam como as medidasverdadeiras (Hx). Então, a MLE da transformação projetiva H e das correspondências {xi ↔ x ′i}é a homografia H e as correspondências corrigidas {xi ↔ xi′} que minimizam∑

i

d(xi, xi)2 + d(x ′i , xi′)2

com xi′ = Hxi.

Neste caso, a MLE é idêntica a minimizar a função de erro de reprojeção.



Distância de Mahalanobis:

No caso geral Gaussiano, pode-se assumir um vetor de medidas X satisfazendo uma função dedistribuição Gaussiana com matriz de covariância Σ. Estes casos são equivalentes a uma matrizde covariância que é múltipla da identidade.

Maximizar a verossimilhança logarítmica ;e então equivalente a minimizar a distância deMahalanobis:

||X − X ||2Σ = (X − X )T Σ−1(X − X )

Havendo erros independentes em ambas as imagens, a função de custo é

||X − X ||2Σ + ||X ′ − X ′||2Σ′

onde Σ e Σ′ são as matrizes de covariância das medidas nas duas imagens.

Assumindo-se que os erros para todos os pontos xi e x ′i são independentes, com matrizes decovariância individuais Σi e Σ′i, respectivamente, então a expressão assima se expande para∑

||xi − xi||2Σi +∑||x ′i − xi′||2Σ′

i



DLT Normalizado:

Os resultados do algoritmo DLT não são invariantes a transformações de similaridade nasimagens.

Alguns sistemas de coordenadas são, então, de alguma forma, melhores do que outros para secomputar uma homografia 2D.

Antes da aplicação do algoritmo DLT, então, os dados devem ser normalizados, através detranslação e escalamento das coordenadas da imagem.

A normalização, então, melhora dramaticamente a acurácia dos resultados e torna o algoritmoinvariante com respeito a escolhas arbitrárias de escala e origem de coordenadas. Assim, nãodeve ser considerada opcional.



Escalamento Isotrópico:

1 Transladar os pontos de modo que seu centróide fique na origem;2 Escalar os pontos de modo que sua distância média da origem seja igual a

√2;

3 Aplicar essa transfor a cada uma das imagens independentemente.

No escalamento não isotrópico, após a translação, o escalamento é feito de modo que os doisprincipais momentos do conjunto de pontos sejam ambos iguais à unidade. Assim, o conjunto depontos formará uma nuvem de pontos circular aproximadamente simétrica de raio 1 ao redor daorigem.

Resultados experimentais sugerem que o esforço extra necessário para o escalamento nãoisotrópico não leva a resultados significativamente superiores aos obtidos com o escalamentoisotrópico.





Métodos de Minimização Iterativos:

Minimizar funções de custo geométricas requer o uso de técnicas iterativas. Estas técnicasgeralmente consistem de cinco passos:

1 Função de Custo: escolhe-se uma função de custo como base para a minimização.2 Parametrização: a transformação a ser computada é expressa em termos de um númerofinito de parâmetros.

3 Especificação da Função: uma função que expresse o custo em termos dos parâmetrosdefinidos deve ser especificada.

4 Inicialização: computa-se uma estimativa inicial adequada dos parâmetros. Isto geralmenteé feito com o uso de um algoritmo linear, como o DLT.

5 Iteração: partindo-se da solução inicial, os parâmetros são iterativamente refinados com oobjetivo de minimizar a função de custo.



Sobre parametrização: a estratégia geral que guia a parametrização é se selecionar um conjuntode parâmetros que cubra o espaço completo sobre o qual se está minimizando, permitindo aomesmo tempo que se compute a função de custo de uma maneira conveniente.

Não há a necessidade de se usar um número mínimo de parâmetros. Inclusive, se descobriuempiricamente que a superfície da função de custo é mais complicada quando se empregaparâmetros mínimos, o que aumenta a possibilidade de se ficar preso num mínimo local.

O uso de parâmetros mínimos também pode acabar por restringir a transformação a uma classeparticular, o que não é desejável.



Sobre especificação da função: em estimação iterativa através do ajuste de parâmetros, umasuperfície modelo S é localmente parametrizada, e estes parâmetros variam para minimizar adistância ao ponto medido. Mais especificamente,

1 Existe um vetor de medidas X ∈ RN com matriz de covariância Σ.2 Um conjunto de parâmetros é representado como um vetor P ∈ RM .3 Um mapeamento f : RM → RN é definido. O alcance deste mapeamento é, ao menoslocalmente, a superfície modelo S em RN representando o conjunto de medidas permitíveis.

4 A função de custo a ser minimizada é a distância de Mahalanobis ao quadrado:

||X − f(P)||2Σ = (X − f(P))tΣ−1(X − f(P))

Em efeito, se está tentando encontrar um conjunto de parâmetros P de modo que f(P) seja omais próximo possível de X com respeito à distância de Mahalanobis.



Erro numa única imagem: fixa-se as coordenadas dos pontos xi na primeira imagem e se variaH para se minimizar a função de custo ∑

i

d(x ′i ,Hxi)2

O vetor de medidas X é formado pelas 2n coordenadas não homogêneas dos pontos x ′i .

Pode-se escolher como parâmetros o vetor h de valores da matriz de homografia H. Assim

f : h 7→ (Hx1,Hx2, ...,Hxn)

onde se entende que Hxi indica as coordenadas não homogêneas.

Verifica-se então que ||X − f(h)||2 é igual à distância geométrica para o erro numa única imagem.



Erro de transferência simétrico: no caso da função de custo simétrica∑i

d(xi,H−1x ′i)2 + d(x ′i ,Hxi)2

escolhe-se como o vetor de medidas X o vetor 4n formado pelas coordenadas não homogêneasdos pontos xi e dos pontos x ′i .

Novamente, o vetor de parâmetros é o vetor h de valores de H, e a função f é definida por

f : h 7→ (H−1x ′1, ...,H−1x ′n,Hx1, ...,Hxn)

Verifica-se que ||X − f(h)||2 é igual à distância geométrica para o erro de transferência simétrico.



Erro de reprojeção: no caso do erro de reprojeção∑i

d(xi, xi)2 + d(x ′i , xi′)2 sujeito a xi′ = Hxi ∀i

o vetor de parâmetros é P = (h, x1, ..., xn).

O vetor de medidas contém as coordenadas não homogêneas de todos os pontos xi e x ′i .

As coordenadas de xi′ não são requeridas porque estão relacionadas aos demais parâmetros porxi = Hxi.

A função f é definida por

f : (h, x1, ..., xn) 7→ (x1, x1′, ..., xn, xn′)

onde xi = Hxi.

Verifica-se que ||X − f(P)||2, com X sendo um vetor 4n, é igual à função de custo do erro dereprojeção.

A aproximação de Sampson para o erro de reprojeção permite que este seja minimizado comsomente 9 parâmetros (os valores da homografia H), ao invés de 2n+ 9 parâmetros.



Iteração de Newton: seja X = f (P), onde X é um vetor de medidas e P é um vetor deparâmetros nos espaços Euclidianos RN e RM , respectivamente.

Um vetor de valores medidos X aproximando os verdadeiros valores X é dado, e desejamosencontrar o vetor P que melhor satisfaz a relação funcional, mais precisamente X = f (P)− ε,para a qual ||ε|| é minimizado.

Inicia-se a o processo com um valor estimado P0 e se refina esta estimativa assumindo-se que afunção f é localmente linear.

Assume-se que a função f é aproximada em P0 por f (P0 + ∆) = f (P0) + J∆, onde J é omapeamento linear representado pela matriz Jacobiana J = ∂f

∂P .

Buscamos o ponto f (P1) com P1 = P0 + ∆ que minimizaf (P1)− X = f (P0) + J∆− X = ε0 + J∆, com ε0 = f (P0)− X .

É necessário, então, se minimizar ||ε0 + J∆|| sobre ∆, que é um problema de minimização linear.



O vetor ∆ é obtido resolvendo-se as equações normais

JT J∆ = −JT ε0

Ou usando-se a pseudoinversa ∆ = −J+ε0.

O vetor solução P é obtido começando-se com uma estimativa P0 e computando-se sucessivasaproximações de acordo com a fórmula

Pi+1 = Pi + ∆i

onde ∆i é a solução para o problema linear de mínimos quadrados

J∆i = −εi

A matriz J é o Jacobiano ∂f∂P avaliado em Pi e εi = f (Pi)− X .

Espera-se que o algoritmo convirja para a desejada solução de mínimos quadrados P. Porém, elepode convergir para um mínimo local, ou mesmo não convergir.

O comportamento do algoritmo depende muito fortemente da estimativa inicial P0.



Iteração Levenberg-Marquardt:

Esta iteração é uma ligeira variação da de Newton. No caso, as equações normaisJT J∆ = −JT ε são substituídas pelas equações normais aumentadas (JT J + λI)∆ = −JT ε, paraum valor λ que varia de iteração para iteração.

Um valor inicial típico de λ é 10−3 vezes a média dos elementos da diagonal de N = JT J .

Se o valor de ∆ obtido resolvendo-se as equações normais aumentadas leva a uma redução doerro, então o incremento é aceito e λ é dividido por um fator (tipicamente 10) antes da próximaiteração. Caso o valor de ∆ leve a um aumento do erro, então λ é multiplicado pelo mesmo fatore as equações normais aumentadas são resolvidas novamente.

Esse processo continua até que seja encontrado um valor de ∆ que leve a um erro reduzido.

Este processo de se resolver equações normais aumentadas repetidamente para diferentes valoresde λ até que um ∆ aceitável seja encontrado constitui uma iteração do algoritmo LM.





Comparações Experimentais:



Algoritmo de alinhamento de imagens:

1 Encontre pontos característicos nas imagens A e B;

2 Encontre as correspondências entre os pontos característicos;

3 Calcule a homografia entre A e B usando, por exemplo, uma solução demínimos quadrados no conjunto de pares de pontos característicos;

⇒ O que pode dar errado?


Alinhamento de Imagens RANSAC

RANSACOutliers



Considere um exemplo simples de regressão linear (mínimos quadrados)

⇒ Como resolver este problema?



Ideia

Dada um reta hipotética:

• Conte o número de pontos que “concordam” com a reta;� “Corcordar” = dentro de uma distância da reta;

� Ponto que “concorda” → inlier;

• Para todas as possíveis retas, escolha aquela com o maior número de inliers;



Contando inliers



Contando inliers



Contando inliers



Como encontrar a melhor reta?

• Diferentemente do método dos mínimos quadrados, não há solução fechada;

• Teste de hipóteses:� Tente diversas retas, escolha a melhor;

� Que retas?













Translações









RANSAC

• O limiar de inliers está relacionado com a quantidade de ruído esperada;� Tipicamente o ruído é modelado como Gaussiano com algum desvio padrão (e.g., 3 pixels);

• O número de rodadas (testes de hipóteses) está relacionado à porcentagemde outliers esperada e à probabilidade de sucesso que queremos garantir

� Suponha que haja 20% de outliers, e queremos encontrar a resposta correta com 99% deprobabilidade.

� Quantas rodadas precisamos?



Quantas rodadas?

• Se temos que escolher s pontos a cada rodada� Com uma taxa de outliers e

� E queremos a resposta correta com probabilidade p

N = log(1− p)log(1− (1− e)s)





Escolha do limiar t:

Na prática, a escolha deste limiar é feita empiricamente.

Caso se assuma que os erros de medida são Gaussianos com média zero e desvio padrão σ,pode-se computar o valor de t.

O quadrado da distância dos pontos, d2⊥, é uma soma de variáveis Gaussianas ao quadrado e

segue uma distribuição χ2m, com m graus de liberdade, onde m é igual à codimensão do modelo.

Para uma reta, a codimensão é 1 (apenas a distância perpendicular à reta é medida). Se omodelo for um ponto, a codimensão é 2, e o quadrado da distância é a soma dos erros de medidade x e y elevados ao quadrado.

A probabilidade que o valor de uma variável aleatória χ2m seja menor do que k2 é dada pela

distribuição cumulativa chi quadrada, Fm(k2) =∫ k2

0 χ2m(ξ)dξ. Para a distribuição cumulativa{

inlier d2⊥ < t2

outlier d2⊥ ≥ t2 com t2 = F−1

m (α)σ2



Usualmente α é escolhido como 0.95, e modo que haja uma probabilidade de 95% de que oponto seja um inlier. Isso significa que um inlier só será incorretamente rejeitado 5% das vezes.

Valores de t para α = 0.95 e para os modelos de interesse estão na tabela abaixo:



Escolha do limiar N:

lembrando que

N = log(1− p)log(1− (1− ε)s)



Escolha do limiar T :

Geralmente encerra-se a execução do algoritmo se o tamanho do conjunto de consenso for similarà quantidade de inliers que se acredita que haja no conjunto de dados, dada a proporçãoassumida de outliers. Ou seja, para n pontos de dados, T = (1− ε)n.





Para as medidas de distância, pode-se escolher, por exemplo, dentre as apresentadas, como oerro de reprojeção.

Quanto à seleção de amostras:

• Deve-se descartar amostras degeneradas;• A amostra deve consistir de pontos com uma boa distribuição espacial sobre a imagem, por

conta do problema de extrapolação: a homografia estimada irá mapear com acurácia aregião de onde forem obtidos os pontos usados na sua computação, mas essa acuráciageralmente se deteriora com o aumento da distância a esta região.




Alinhamento de Imagens IPC

IPC: Interactive Closest Point

⇒ Não requer casamento inicial de características;

⇒ Assume que os pontos mais próximos são correspondentes;


Alinhamento de Imagens IPC

Algoritmo:

1 Defina as correspondências: atribua cada características na imagem 1 à sua vizinha maispróxima na imagem 2.

2 Estime os parâmetros da transformação que atende o casamento entre as características;

I Qualquer algortitmo pode ser usado (mínimos quadrados, SVD, etc.);

3 Transforme os pontos da imagem 1 usando os parâmetros estimados;

4 Repita os passos 1 a 3 até que a mudança seja pequena (critério de parada);

⇒ Converge para mínimo local;


Documents

Visão Computacional AlinhamentodeImagens · SergioL.Netto ProgramadeEngenhariaElétrica-COPPE/UFRJ Laboratório de Sinais, Multimídia e Telecomunicações [email protected] Julhode2017