Italian English French German Portuguese Russian Spanish

Analisi anomalie foto Apollo con strumenti matematici

Di più
14/03/2017 17:37 #8907 da doktorenko
Ti rispondo ad un messaggio della pagina precedente:

Altra cosa non chiara è perchè cerchi di stimare la rotazione tra le due foto, se hai a disposizione le due matrici camera K1, K2 relative alle due immagini Apollo, quindi con K1,K2: R^3-->P^2 che ad ogni punto della tua ricostruzione 3D associano il corrispondente punto sulle due immagini Apollo, allora l'omografia H fra l'immagine 1 e la 2 è H12=K2 K1^-1 e non c'è bisogno di stimare la rotazione.


La rotazione non la stimo, me la fornisce il programma. Io vorrei usare la tua formula, pero` prima devo capire come sono composte le matrici Kn.
Ho provato cosi`:
F = focale in pixel
px = lato immagine in px
matrice calibrazione camera K = [ [ F, 0, px/2 ], [ 0.0, F, px/2 ], [ 0.0, 0.0, 1.0 ] ]

Kn = Rn * K
(R: matrice rotazione associata alla camera n)

La matrice di trasformazione H 1->2:
H = K2*K1^1

Funziona, ma solo per due angoli su tre.

La mia difficolta` principale e` che di questi argomenti non ho trovato molte informazioni; devo quindi procedere per intuito e tentativi.

Di più
14/03/2017 18:13 #8908 da FranZeta

doktorenko ha scritto: Ho provato cosi`:

...
Kn = Rn * K
(R: matrice rotazione associata alla camera n)

La matrice di trasformazione H 1->2:
H = K2*K1^1

Funziona, ma solo per due angoli su tre.

Così ad occhio direi che una prima correzione sarebbe questa:
Kn=K*Rn
(il prodotto fra matrici non è commutativo)

La mia difficolta` principale e` che di questi argomenti non ho trovato molte informazioni; devo quindi procedere per intuito e tentativi.

Idem! Ho studiato geometria proiettiva astratta, la computer grafica è un tantino più improntata all'aspetto pratico...

Di più
16/03/2017 08:22 - 16/03/2017 16:02 #8918 da doktorenko
Vado a illustrare la procedura per ottenere una coppia stereoscopica di immagini mediante omografia.
1) per prima cosa si individuano dei punti nell`immagine di riferimento A e quelli corrispondenti nell`immagine da raddrizzare B (punti An e Bn).
1b) i punti An devono essere "mobili" (subire la parallasse dovuta allo spostamento della camera), e giacere sullo stesso piano.
2) si stima una prima matrice omografica H: H*B=A; Cn=H*Bn
3) con la forza bruta (ma si possono usare anche altri metodi) si corregge H finche` le rette rn:(Cn,An) non convergono nello stesso punto O.



I punti individuati nelle immagini sono (li copio da geogebra, perdonate l`eccesso di precisione nei punti):

A1 (1868.11769988769, 692.883470019046)	
A2 (1990.48521852825, 775.234380818113)	
A3 (2435.06354430023, 864.774779267219)


B1 (2427.28841864788, 727.392232676494)
B2 (2543.57337446140, 848.159902761366)
B3 (3045.03522141262, 984.230141465186)

Matrice omografica:
H = [[1.0087217782, 0.0492298329, -224.6818216166],
    [-0.0048663276, 1.0336044935, -441.518876724],
    [0.0000043189, 0.0000177307, 0.9504399531]]

Pero` ho un dubbio: nel calcolare la correzione della matrice H come nel punto 2, devo tener conto anche della traslazione o solo della rotazione?
Ultima modifica: 16/03/2017 16:02 da doktorenko. Motivo: matrice omografica piu` precisa

Di più
19/03/2017 03:51 #8953 da FranZeta

doktorenko ha scritto: Pero` ho un dubbio: nel calcolare la correzione della matrice H come nel punto 2, devo tener conto anche della traslazione o solo della rotazione?

Inizio dalla fine: se intendi la traslazione della camera non puoi tenerne conto tramite un'omografia perchè ti serve un'omografia nello spazio 3D. Ma questo aspetto è secondario rispetto a quanto scrivi nella prima parte. Siccome la questione tecnica è parecchio più complessa di come l'hai messa giù risponderò dopo aver studiato bene la questione, onde evitare di risponderti con sciocchezze. Nel frattempo però potresti iniziare anche tu a pensare bene a quello che vuoi effettivamente ottenere, per esempio quando dici:

si stima una prima matrice omografica H: H*B=A

Ora finchè assumevamo lo sfondo fisso potevamo cercare un'omografia tale che A=H*B a patto che con A e B si intendano le porzioni comuni di sfondo fra le due immagini, se anche lo sfondo ha parallasse non esiste nessun punto per cui valga la relazione sopra, per definizione.

La condizione che le rette AnCn siano concorrenti esprime il fatto che le camere hanno la stessa orientazione (ma andrebbe dimostrato!), quindi in sostanza la procedura serve a collimare le due camere, non fornisce però un'omografia fra le due immagini, ma su questo torneremo in seguito.

Di più
20/03/2017 07:33 #8961 da doktorenko
Per la complessita` hai ragione, e infatti prima di inviare altri "raddrizzamenti" voglio sperimentare bene con i modelli tridimensionali. Aggiungo solo che l`immagine del messaggio precedente si riferisce ancora alla fase 1 (stima con il vecchio sistema); la fase 2 (forza bruta) funziona, pero` trova piu` combinazioni di angoli per cui le rette individuate concorrono: a me invece servirebbe una posizione unica.

Adesso passo alla stima della camera nella foto AS-11-40-5866 (quella del famigerato controluce lunare); la procedura si divide in due fasi:

1) calcolo della posizione e angolazione della camera relativa al modulo lunare centrato e in bolla:
x: 2.7478 m
y: -7.67029 m
z: 1.68243 m

ax: 103.251°
ay: 0.798°
az: 34.739°

2) fatto cio`, vincolo i movimenti del modulo alla camera, e ruoto quest`ultima per collimare il paesaggio lunare con quello del modellato:
ax: 101°
ay: 2.5°
az: -44.8°

Angoli del Sole (ricavati da Stellarium):
75.672°
90.171°

Risultato:





Sto ricostruendo il modulo (semplificato) perche` i modelli liberamente disponibili mi sembrano poco accurati nelle misure.

Di più
20/03/2017 15:41 #8964 da FranZeta

doktorenko ha scritto: la fase 2 (forza bruta) funziona, pero` trova piu` combinazioni di angoli per cui le rette individuate concorrono: a me invece servirebbe una posizione unica.

Per questo ti invitavo a chiarire quello che vuoi ottenere! Applicando omografie "casuali" puoi far convergere le rette in molti modi, e verso punti diversi. Prima di parlare meglio dei procedimenti per ricostruire la scena 3D è il caso di linkare questa esposizione dell'argomento, molto chiara e non eccessivamente tecnica (per confronto questa invece è un'esposizione tecnica), in un prossimo messaggio voglio esporre la sintesi con allegato l'esempio pratico fatto sulle due foto Apollo 15.

A pagina 248 del primo documento linkato trovi il tuo metodo delle rette concorrenti, tieni presente che tutte le rette individuate da punti corrispondenti convergono verso un unico punto, che rappresenta la direzione della traslazione fra i due scatti, o meglio il punto all'infinito corrispondente al vettore traslazione. Anche in questo caso quindi la complanarità non c'entra. Inoltre la condizione di concorrenza delle rette è sufficiente ma non necessaria: se c'è solo uno spostamento perpendicolare alla direzione delle inquadrature (quindi la traslazione è solo nelle due direzioni individuate dal piano della foto, destra-sinistra e alto-basso) le rette convergono verso un punto all'infinito nel piano dell'immagine, vale a dire che sono parallele. Quest'ultimo caso però fortunatamente non si verifica nelle foto prese in considerazione.

Di più
24/03/2017 21:36 #8980 da kamiokande

doktorenko ha scritto: Sto ricostruendo il modulo (semplificato) perche` i modelli liberamente disponibili mi sembrano poco accurati nelle misure.

Forse potrebbe esserti utile questo sito http://www.ehartwell.com/LM/SCATSystems.htm


"La stampa è morta" (Egon Spengler - Ghostbuster)

Di più
27/03/2017 18:29 - 28/03/2017 02:07 #8994 da FranZeta
Come promesso inizio a esporre una sintesi del processo matematico usato per trattare immagini stereoscopiche, lo scopo del quale è ricostruire per quanto possibile la scena 3D ripresa da due fotocamere. Il testo di riferimento è quello già linkato sopra, si tratta del capitolo di un libro che purtroppo non possiedo, era nei riferimenti bibliografici nella pagina di wikipedia sulla matrice fondamentale, che come suggerisce il nome è l'oggetto cardine di tutta la teoria. Cercherò di fare il possibile per sorvolare gli aspetti tecnici che lascio alla vostra voglia di approfondire, qui più che altro voglio applicare il metodo alla coppia di foto di Apollo 15.

Un minimo di teoria però è necessario introdurla, quindi partiamo da questa figura del documento di riferimento:



La situazione dunque è questa: due fotocamere riprendono la stessa scena da due diverse posizioni (C e C') puntate verso direzioni qualunque, seguendo la notazione del libro indichiamo con X il punto della scena 3D che viene mappato nel punto x della prima immagine e x' della seconda. Nelle due immagini esistono due punti fondamentali: gli epipoli e ed e', questi risultano essere i due punti corrispondenti alla posizione della fotocamera che scatta rispettivamente la seconda e la prima immagine. Ora vediamo il motivo della loro importanza. Supponiamo di sovrapporre le due immagini (cosa che faremo tra breve), al punto x della prima corrisponderà un qualche punto x' della seconda, questo punto deve per forza trovarsi su una retta l' che parte dall'epipolo e', l' rappresenta il raggio che congiunge la fotocamera 1 col punto X della scena 3D. Stessa cosa vale invertendo i ruoli con x, e, l (si vede la rappresentazione grafica nella figura 9.1 b).

Quindi ad ogni punto x della prima immagine (d'ora in poi semplicemente 1) corrisponde una retta l' in 2 e viceversa ad x' in 2 corrisponde l in 1. Qui entra in gioco la matrice fondamentale F che rappresenta la funzione x-->l', mentre la funzione opposta x'-->l è data da F^T (trasposta), risulta perciò:
l'=F*x
l=F^T*x'
Per ragguagli sulle coordinate da usare in queste formule consiglio di guardare nel documento, non mi pare il caso di appesantire ulteriormente. Ci basti sapere che F è una matrice 3x3 di rango 2 ed è un oggetto proiettivo proprio come le omografie, quindi con una coordinata in più dell'effettiva dimensione, questi due fatti fanno sì che i coefficienti liberi di F da 9 passino a 7 (-1 perchè è proiettiva e -1 perchè det(F)=0). La matrice F stabilisce anche la seguente relazione basilare fra punti corrispondenti x e x':
1) x'^T*F*x=0
L'equazione qui sopra esprime algebricamente quanto già detto, ossia che x sta sulla retta l e x' su l'. Ogni coppia (x,x') che soddisfi la 1) rappresenta due punti immagini dello stesso punto X nella scena reale. Siccome sviluppando la 1) si ottiene un'equazione lineare nelle 9 incognite (F)_ij dei coefficienti della matrice, e siccome solo 7 di questi sono indeterminati, la matrice F può essere ricavata stabilendo 7 corrispondenze fra punti delle due immagini, cioè si risolve il sistema:
x1'^T*F*x1=0
......
x7'^T*F*x7=0
cosa che adesso facciamo subito per la coppia di immagini Apollo 15. Iniziamo con l'incollare le due immagini una sull'altra senza alcuna modifica (io l'ho fatto in Matlab usando la funzione imfuse):



Poi individuiamo 7 coppie di punti "comodi" corrispondenti, per esempio il primo punto che ho scelto è la stellina in basso a sinistra della bandiera americana. Una nota sulle coordinate: per fare questo lavoro ho copiato l'immagine sopra in geogebra, l'ho ridimensionata a 40x40 (la misura originale in mm) e poi l'ho leggermente spostata in modo da allineare la croce centrale con l'origine degli assi, visto che non era perfettamente centrata. Infine ho preso pari pari le coordinate dei punti con due decimali e ho passato il tutto al caro Maple che ha fatto i calcoli per me, ecco il risultato:



Come si può vedere per prima cosa calcolo la matrice F, si tenga presente che essendo proiettiva se moltiplico tutto per una costante diversa da 0 si ottiene sempre la stessa matrice, poi sotto ho verificato che det(F)=0 (circa!), mentre sotto ancora ci sono altri calcoli che vado a spiegare subito. Serve solo ancora un pochino di teoria, abbiate pazienza...
...allora una volta nota F è conveniente decomporla nella somma di una matrice simmetrica S e una antisimmetrica A, esplicitamente si fa così:
F=S+A
S=(F+F^T)/2
A=(F-F^T)/2
L'espressione esplicita delle due matrici è la seconda parte del calcolo riportato sopra. Alla matrice S corrisponde una conica nel piano delle due immagini incollate, mentre ad A corrisponde il punto x_a nello stesso piano. Avendo a disposizione conica e punto si possono fare le seguenti costruzioni geometriche:



Nella figura 9.10 a) vediamo la prima costruzione: la retta l_a è la polare di x_a rispetto alla conica F_S, questa retta interseca F_S in due punti che risultano essere i nostri epipoli. Gli epipoli possono essere calcolati direttamente risolvendo i sistemi:
F*e=0
F^T*e'=0
e questa è l'ultima parte del calcolo riportato sopra, le coordinate [x,y] sono quelle di e e le [X,Y] di e', se abbiamo fatto le cose bene queste devono corrispondere a quelle dei due punti trovati con la costruzione geometrica della retta polare, cosa che in effetti avviene con buona approssimazione, ecco infatti l'immagine finale:



Oltre agli oggetti già spiegati, nella figura ci sono due rette arancioni che corrispondono alla costruzione della figura 9.10 b): dato il punto x (nel nostro caso è P2) nella prima immagine, per trovare la retta l' ad esso associata si traccia prima la retta l passante per x ed e (P2 e E nell'immagine Apollo), questa interseca la conica in un punto x_c, infine la retta l' cercata è la retta passante per x_c ed e' (E'). In pratica se abbiamo fatto le cose per bene, per ogni coppia di punti Pi,Pi' sull'immagine deve valere che le due rette PiE e Pi'E' si intersecano in un punto della conica. Qui c'è la verifica per alcuni dei punti utilizzati:



Come sempre, in questi casi, bisogna tenere presente che da una parte c'è l'iperuranio, cioè il mondo platonico delle idee matematiche astratte e perfette, dall'altra c'è il mondo reale con tutte le sue imperfezioni, nello specifico la realtà ha le fattezze delle nostre foto Apollo e le intersezioni evidenziate dai circoletti sono solo approssimativamente sulla conica, ingrandendo le parti evidenziate si troverebbero delle discrepanze. Per chi fosse interessato questo è il link per il file di geogebra e questo per fare delle eventuali modifiche ad un file gemello modificabile a piacere.

Per ora mi fermerei qui, nella prossima puntata vedremo di cercare di ricavare una coppia di possibili matrici camera per le due foto Apollo.
Ultima modifica: 28/03/2017 02:07 da FranZeta. Motivo: sistemati un paio di apostrofi invertiti

Di più
28/03/2017 09:18 - 28/03/2017 14:50 #8996 da doktorenko
@FranZeta
nell`attesa della seconda parte, ho provato a calcolare la matrice fondamentale F con il modulo OpenCV, usando questo codice:

Geometria epipolare

A parte naturalmente il nome delle immagini (sx=11601; dx=11602), ho dovuto modificare solo questa riga:
sift = cv2.SIFT()
in:
sift = cv2.xfeatures2d.SIFT_create()

Risultato (mi ha trovato un sacco di punti utili!):


Matrice F:
array([[  1.54590733e-08,   1.31312522e-06,  -2.31504371e-03],
       [ -1.29520987e-06,   7.64118313e-08,  -1.06896770e-03],
       [  2.74047832e-03,   5.50595730e-04,   1.00000000e+00]])


Aggiungo il codice (in python/numpy) della prima parte dei tuoi calcoli:
F=np.matrix([[0.0002249705489,0.02234261333,-0.1725523678],\
             [-0.01959259243,0.005407583090,-0.2555901116],\
             [0.2428871895,0.2098212272,1]])

S=(F+F.transpose())/2
A=(F-F.transpose())/2

M=np.delete(np.array(F),1,0)

a=M[:,:-1]
b=-M[:,2]

x,y=np.linalg.solve(a,b)

M=np.delete(np.array(F).transpose(),1,0)

a=M[:,:-1]
b=-M[:,2]

X,Y=np.linalg.solve(a,b)

>>> x
-10.883432648836861
>>> y
7.8326029740642245
>>> X
-12.357155656359714
>>> Y
12.254998631089576

Invece, con la matrice F che ho calcolato io, le coordinate dei punti epipolari sono:
>>> px=4175
>>> mm=55.44
>>> X/px*mm
-7.1978876554283415
>>> Y/px*mm
28.010637232681304
>>> x/px*mm
-9.5717150249253162
>>> y/px*mm
23.523676064212165

Il risultato e` parecchio diverso, forse dipende anche dal fatto che dovrei filtrare la scelta in modo da avere meno coppie di punto, ma individuate piu` precisamente.
Ultima modifica: 28/03/2017 14:50 da doktorenko. Motivo: Aggiunto codice python

Di più
28/03/2017 17:28 #8997 da doktorenko
Ho provato a costruire il grafico con la soluzione che mi da` Opencv:
F=np.matrix([[  1.57550907e-05,  -7.16872070e-05,   1.09271677e-01],
       [  7.84514021e-05,  -8.04395058e-06,  -2.95069542e-02],
       [ -1.72668082e-01,   5.65996434e-02,   1.00000000e+00]])

e= (544.7,1644)
e'= (555, 2089.5)

S=(F+F.transpose())/2

Conica:
a11,a12,a13,a12,a22,a23,a13,a23,a33=np.array(S).flatten()
print ("%.9fx^2+2 * %.9fxy+%.9fy^2+2 * %.9fx+2 * %.9fy+%.9f=0" % (a11,a12,a22,a13,a23,a33) )
>>> 0.000015755x^2+2 * 0.000003382xy+-0.000008044y^2+2 * -0.031698203x+2 * 0.013546345y+1=0

Coppie di punti:
A1(342, 2184)	B1(364, 2485)
A2(583, 951)     B2(677, 1375)
A3(910, 1488)	B3(948, 1867)
A4(844, 2629)	B4(737, 2860)

Risultato (ho cambiato le proporzioni degli assi):


Allora anche questa soluzione dovrebbe andare bene? e come si fa a passare dalla matrice S all`equazione della conica?

Di più
29/03/2017 02:14 - 29/03/2017 02:58 #8998 da FranZeta

doktorenko ha scritto: Allora anche questa soluzione dovrebbe andare bene? e come si fa a passare dalla matrice S all`equazione della conica?

No, la soluzione è una sola e visto che quella che ho trovato è verificata empiricamente con le intersezioni dei punti (seppure approssimate), la soluzione reale deve essere molto vicina.
La conica associata a una matrice simmetrica S si trova, posto X=(x,y,1) sviluppando il prodotto:
X^T*S*X=0
In geogebra c'è il comando Conica[<Numero>...<Numero>] dove bisogna inserire 6 numeri, sono, nell'ordine i coefficienti S_11,S_22,S_33,2*S_12,2*S_13,2*S_23 della matrice S.

PS Poi non so esattamente come hai trovato la conica, non sono stato lì a controllare il codice, comunque la direzione della camera dovrebbe essere il centro dell'immagine, a meno che il programma non tenga già conto dello spostamento, altrimenti è normale trovare risultati differenti. Inoltre ho fatto i conti in mm, non in pixel. Conta anche l'ordine in cui si prendono le due immagini, non per la conica ma per la matrice F.

PPS I punti Pi e Pi' devono essere corrispondenti!!! Altrimenti puoi fare convergere le rette dove vuoi (come sembrerebbe dalla tua ultima figura)!
Ultima modifica: 29/03/2017 02:58 da FranZeta.

Di più
29/03/2017 04:38 - 29/03/2017 13:05 #8999 da doktorenko


Le differenze tra i due grafici dovrebbero essere:

1) i miei calcioli sono i pixel e non in mm (10 mm = 753 px)
2) il centro dell`immagine non e` all`origine degli assi ma a O(4175/2, 4175/2)
3) l`asse y e` positivo verso il basso (convenzione delle immagini informatiche)
4) forse ho usato come riferimento l`altra immagine rispetto alla tua (ordine dei punti e-e')
5) la precisione nell`individuazione dei punti di OpenCV e ` intera.

Adesso cerchero` di usare la tua convenzione per tracciare il grafico.

p.s.
come hai fatto a calcolare la dimensione di 40 mm di lato? Perche` io la tua non la vedo tagliata: i 40 mm corrispondo al quadrato interno alla griglia di 5x5 crocette.

aggiungo anche il codice per trovare la matrice fondamentale partendo dai punti:
import numpy as np
import cv2

pts1 = np.array(\
    [[ -15.21,  7.15],\
       [ -12.95, 5.57],\
       [ -9.77,7.73],\
       [ -8.05, 7.65],\
       [ -2.88, 0.8],\
       [ 9.63, 8.17],\
       [ 12.13, 8.86]])

pts2 = np.array(\
    [[ -16.23,  11.48],\
       [ -13.75, 9.74],\
       [ -10.44, 12.08],\
       [ -8.5, 12.02],\
       [ -2.78, 4.62],\
       [ 11.51, 12.98],\
       [ 14.95, 13.98]])

F, mask = cv2.findFundamentalMat(pts1,pts2,cv2.FM_LMEDS)
F=np.matrix(F[0:3])

>>> F
matrix([[  2.24970008e-04,   2.23427304e-02,  -1.72553143e-01],
        [ -1.95926949e-02,   5.40761741e-03,  -2.55590089e-01],
        [  2.42888271e-01,   2.09820763e-01,   1.00000000e+00]])

>>> x, y
(-10.883362887762464, 7.8325956676863848)
>>> X,Y
(-12.357087076459413, 12.254990884928457)


Con i tuoi punti ottengo gli stessi risultati, anche se non so perche` il programma mi restituisce addirittura 3 coniche (nel codice ho mostrato solo la prima).
Ultima modifica: 29/03/2017 13:05 da doktorenko. Motivo: aggiunte

Di più
29/03/2017 17:12 #9000 da FranZeta

doktorenko ha scritto: Le differenze tra i due grafici dovrebbero essere:

1) i miei calcioli sono i pixel e non in mm (10 mm = 753 px)
2) il centro dell`immagine non e` all`origine degli assi ma a O(4175/2, 4175/2)
3) l`asse y e` positivo verso il basso (convenzione delle immagini informatiche)
4) forse ho usato come riferimento l`altra immagine rispetto alla tua (ordine dei punti e-e')
5) la precisione nell`individuazione dei punti di OpenCV e ` intera.


In realtà io mi riferivo alla differenza nelle matrici, il grafico dovrebbe essere pressochè identico perchè la posizione sull'immagine della conica e degli epipoli non dipende dal sistema di coordinate usato. Infatti nella tua immagine ci sono evidenti discrepanze:




Qui ho tracciato le due epilinee passanti per due punti corrispondenti sulla zampa del LEM che sono ben lontane dall'intersecarsi sulla conica.

p.s.
come hai fatto a calcolare la dimensione di 40 mm di lato? Perche` io la tua non la vedo tagliata: i 40 mm corrispondo al quadrato interno alla griglia di 5x5 crocette.

Sono andato a memoria, mi pareva di ricordare che la pellicola fosse 40x40. Questo però non modifica la sostanza dei calcoli fatti sopra, è solo un cambio di unità di misura. Già che ci siamo, qual era la dimensione effettiva della pellicola?

Di più
29/03/2017 19:01 - 29/03/2017 19:04 #9001 da doktorenko
Non ti preoccupare troppo dei miei esperimenti: sono fatti per capire se e` possibile affidare ad un algoritmo la scelta dei punti.
In questa prova ho cercato di rispettare la tua convenzione: dopo la scelta dei punti (automatica) ho trasformato le coordinate in questo modo:
px=4175.
mm=55.44

pts1=(((pts1-px/2)*[1,-1])*55.44)/px
pts2=(((pts2-px/2)*[1,-1])*55.44)/px

La matrice F calcolata dalle 12000 (!) coppie di punti individuati :
F, mask = cv2.findFundamentalMat(pts1,pts2,cv2.FM_LMEDS)

>>> F
matrix([[  7.17589983e-05,  -5.93437573e-03,   2.09125440e-02],
        [  5.89581319e-03,   3.53392050e-04,   2.15287984e-01],
        [  1.01037585e-02,  -2.04270651e-01,   1.00000000e+00]])

La conica
c: 0.29x² - 0.16x y + 1.43y² + 127.06x + 45.14y = -4096

Il grafico:



Non si riesce a vedere, ma i punti che hai individuato come errati precedentemente ora convergono (approssimativamente) sulla conica.
Come ho scritto in precedenza, bisognerebbe abbassare il numero dei punti, aumentandone pero` la precisione nella scelta; purtroppo non ho ancora capito bene su quali parametri agire per ottenere questo risultato.

La dimensione della diapositiva originale la calcolo in questo modo:
spazio  (10 mm)  tra due crocette = 753 px
dimensione immagine: 4175 px
4175 px / 753 px * 10 mm  = 55.44 mm
Ultima modifica: 29/03/2017 19:04 da doktorenko.

Di più
30/03/2017 12:06 - 30/03/2017 16:38 #9003 da doktorenko
Sto cercando di capire come ricavare i punti Pxyz partendo dai punti corrispondenti Puv e P'uv.

Procedimento:

1) ricavare la matrice essenziale E (la mat. F l`abbiamo calcolata in precedenza):
lunghezza focale f = 0.0611 m
sensore f = 0.05544 m

K=np.matrix([ [ f, 0, l/2 ], [ 0.0, f, l/2 ], [ 0.0, 0.0, 1.0 ] ])
E = K**-1* F * K

A questo punto posso usare la funzione OpenCV recoverPose:
points, R, t, mask = cv2.recoverPose(E, pts1, pts2)

R, t: matrice di rotazione e vettore traslazione di una camera rispetto all`altra
(ma non capisco perche` per calcolarle ha bisogno ancora delle coppie di punti pts1 e pts2)

Risultato:
>>> R
array([[ 0.99302249, -0.11715838,  0.01342568],
       [ 0.11557984,  0.94435246, -0.30796062],
       [ 0.0234016 ,  0.30736355,  0.95130437]])
>>> t
array([[-0.99929059],
       [-0.03747715],
       [ 0.00371336]])


Con questi dati ora si possono ottenere i punti tridimensionali mediante triangolazione: io uso la funzione OpenCV triangulatePoints:
X = cv2.triangulatePoints( P1, P2, pts1.transpose(), pts2.transpose() )
X /= X[3]

Pn = matrice di proiezione della camera n
X = punti in coordinate (x,y,z,w).


L`ultima cosa che mi serve capire e` come ottenere la matrice 4x3 di proiezione P della camera; dovrebbe essere composta da:
K: matrice calibrazione 3x3
R: matrice rotazione 3x3
t: traslazione 3x1

P= K*[R|t]
Ultima modifica: 30/03/2017 16:38 da doktorenko. Motivo: correzione

Di più
31/03/2017 02:00 - 31/03/2017 02:08 #9006 da FranZeta

doktorenko ha scritto: Sto cercando di capire come ricavare i punti Pxyz partendo dai punti corrispondenti Puv e P'uv.

Questo è l'obiettivo dichiarato di tutta la procedura, ma è anche l'ultimo dei problemi, note le matrici fondamentale e essenziale come si fa è spiegato qui . Uno dei problemi veri è che se non si conosce almeno una misura sulla scena reale, sia questa la distanza fra i due scatti o le coordinate 3D di un punto, non ci sono santi nè madonne nè opencv che possano fare il calcolo. Il vettore spostamento fra i due scatti è noto solo a meno di un fattore di proporzionalità, vale a dire che si può conscere solo la retta della spostamento, per ricavare i numeri giusti pensavo di usare le misure della larghezza del lunar rover.

Per quanto riguarda le matrici rotazione quelle possibili sono due e non c'è modo a priori di scegliere la migliore, l'unica è fare delle verifiche. Cosa che pian pianino sto facendo. Altra cosa: per come sono sistemate le fotocamere la conica di cui si parlava sopra deve essere un'iperbole.

PS A me risulta E=K'^T*F*K, se K'=K: E=K^T*F*K.
Ultima modifica: 31/03/2017 02:08 da FranZeta.

Di più
31/03/2017 10:31 - 31/03/2017 12:29 #9009 da doktorenko

FranZeta ha scritto: Questo è l'obiettivo dichiarato di tutta la procedura, ma è anche l'ultimo dei problemi, note le matrici fondamentale e essenziale come si fa è spiegato qui . Uno dei problemi veri è che se non si conosce almeno una misura sulla scena reale, sia questa la distanza fra i due scatti o le coordinate 3D di un punto, non ci sono santi nè madonne nè opencv che possano fare il calcolo.


Sarebbe, come dicono i tedeschi, zu schoen, um wahr zu sein.

Il vettore spostamento fra i due scatti è noto solo a meno di un fattore di proporzionalità, vale a dire che si può conscere solo la retta della spostamento, per ricavare i numeri giusti pensavo di usare le misure della larghezza del lunar rover.


La modifica nel mio messaggio precedente era proprio per questo motivo: prima avevo scritto "distanza", poi ho corretto con "vettore traslazione"; intendendo il vettore unitario.

Per ricavare il fattore di scala lambda , ho trovato questa equazione:
(sx2, sy2, 1) = K(R|lambda*T)(X,Y,Z)

sx,sy dovrebbero essere le coordinate sulla proiezione; X,Y,Z il corrispondente punto tridimensionale di cui si conoscono i valori.

PS A me risulta E=K'^T*F*K, se K'=K: E=K^T*F*K.

Si`, ho sbagliato a copiarla:

Correlazione tra matrice fondamentale e essenziale (9.12)

Pero` non riesco lo stesso a ricavarla correttamente; cosi` me la faccio calcolare da OpenCV dalle coppie di punti:
E=cv2.findEssentialMat(pts1,pts2,focale,(cx,cy))[0]

Questa e` l`immagine dei punti tridimensionali ricostruiti:



Potrebbe essere utile per ricavare con maggiore precisione l`altimetria nella zona ai piedi del fotografo, dove i punti sono piu` densi.

p.s.
Aggiungo il codice che ottiene il risultato apparentemente corretto dell` immagine sopra:

Calcolo matrice K (metto le misure in px):
l_px=4175.
f_mm=61.1
l_mm=55.44

f_px=f_mm/l_mm*l_px

K=np.matrix([ [ f_px, 0, l_px/2 ], [ 0.0, f_px, l_px/2 ], [ 0.0, 0.0, 1.0 ] ])

Calcolo e scomposizione matrice E in matrice R e t:
E=cv2.findEssentialMat(pts1,pts2,f,(px/2,px/2))[0]
points, R, t, mask = cv2.recoverPose(E, pts1, pts2)

Valori R e t (quelli del messaggio precedente sono errati):
>>> R
array([[ 0.99953627, -0.01286357,  0.02760033],
       [ 0.01063138,  0.99677412,  0.07955081],
       [-0.02853461, -0.07922049,  0.99644864]])

>>> t
array([[-0.49655057],
       [ 0.00879161],
       [ 0.86796327]])

Triangolazione:
R1 = np.eye(3)
R2 = R

t1=np.zeros(3)[:,np.newaxis]
t2=t

P1=np.matmul(K,np.hstack((R1,t1)))
P2=np.matmul(K,np.hstack((R2,t2)))

X = cv2.triangulatePoints( P1[:3], P2[:3], pts1.transpose(), pts2.transpose() )
X /= X[3]
Ultima modifica: 31/03/2017 12:29 da doktorenko. Motivo: aggiunta codice

Di più
31/03/2017 17:11 #9011 da FranZeta
Anche se devo ancora controllare i calcoli e fare un'animazione gif delle due immagini per vedere se effettivamente la corrispondenza è buona, a questo punto anticipo quello che avevo trovato perchè i valori sono abbastanza diversi dai tuoi. Questa è la prima parte:



Prima ci sono le due matrici di calibrazione (il centro delle mie immagini ha coordinate (0,0) quindi non c'è nessuna traslazione e le matrici sono diagonali), con queste ricavo E, mentre l'ultima parte è la fattorizzazione SVD di E che serve per trovare le sue due componenti R e t, che poi calcolo qui:



Ogni volta che compare una matrice di rotazione verifico che sia effettivamente ortogonale (A*A^T=I), l'ho fatto per U,V e R, questi sono gli altri conti che compaiono sopra.

Il vettore t trovato non mi dà preoccupazioni perchè mi sembra coerente con la scena, diciamo che corrisponde per esempio ad uno spostamento in avanti di 40 cm, a sinistra di 10 e in alto di altri 10 fra i due scatti, la rotazione invece è quella che mi sembra più consona delle due possibili, l'ultima parte del calcolo è per trovare l'angolo in gradi e l'asse di rotazione associati. Come vedi ci sono differenze non marginali con i tuoi risultati, prima di pensare a ricostruire la scena 3D sarebbe meglio verificare il tutto, quando riesco faccio una gif delle due immagini sovrapposte dopo la rotazione, mi sembra il controllo migliore, però sulla parte precedente dei calcoli (conica ed epipoli) ho pochi dubbi, fossi in te verificherei se con i tuoi dati ottieni risultati simili, per esempio se ti risultano ancora gli epipoli fuori dall'immagine c'è qualcosa che non va.

Di più
31/03/2017 18:44 - 31/03/2017 18:45 #9012 da doktorenko

FranZeta ha scritto: Come vedi ci sono differenze non marginali con i tuoi risultati, prima di pensare a ricostruire la scena 3D sarebbe meglio verificare il tutto, quando riesco faccio una gif delle due immagini sovrapposte dopo la rotazione, mi sembra il controllo migliore, però sulla parte precedente dei calcoli (conica ed epipoli) ho pochi dubbi, fossi in te verificherei se con i tuoi dati ottieni risultati simili, per esempio se ti risultano ancora gli epipoli fuori dall'immagine c'è qualcosa che non va.


Faccio solo questa puntualizzazione: i risultati che ho pubblicato sopra non li considero propriamente "miei", perche` non ho fatto altro che usare il codice di esempio di OpenCV, senza modificarne i parametri; la mia rotazione- ottenuta cioe` con il sistema della collimazione della foto con il modello ricavato dai dati altimetrici- e` quella allegata all`ultima omografia.

Lo spostamento lo sto calcolando con la "forza bruta" applicata ad un modello LM con misure piu` accurate, da usare per la comparazione con lo scatto Apollo.
Ultima modifica: 31/03/2017 18:45 da doktorenko.

Tempo creazione pagina: 0.312 secondi

Per migliorare il nostro servizio, la tua esperienza di navigazione e la fruizione pubblicitaria questo sito web utilizza i cookie (proprietari e di terze parti). Per maggiori informazioni (ad esempio su come disabilitarli) leggi la nostra Cookies Policy. Chiudendo questo banner, scorrendo questa pagina o cliccando qualunque suo elemento acconsenti all'uso dei cookie. INFO