correzione della distorsione del fisheye a livello di programmazione

AGGIORNAMENTO STATUS BOUNTY:

Ho scoperto come mappare una lente lineare , dalle coordinate di destination alle coordinate di source .

Come calcoli la distanza radiale dal centro per passare dal fisheye al rettilineo?

  • 1). In realtà faccio fatica a capovolgerlo e per mappare le coordinate della sorgente alle coordinate della destinazione. Qual è l’inverso, nel codice nello stile delle funzioni di conversione che ho postato?

  • 2). Vedo anche che la mia imperfezione è imperfetta su alcuni obiettivi, presumibilmente quelli che non sono strettamente lineari. Qual è l’equivalente a-e-da le coordinate di origine e destinazione per tali obiettivi? Ancora una volta, più codice che solo formule matematiche per favore …


Domanda come originariamente dichiarato:

Ho alcuni punti che descrivono le posizioni in una foto scattata con un objective fisheye.

Voglio convertire questi punti in coordinate rettilinee. Voglio smascherare l’immagine.

Ho trovato questa descrizione su come generare un effetto fisheye, ma non su come invertirlo.

C’è anche un post sul blog che descrive come utilizzare gli strumenti per farlo; queste immagini sono da quella:

(1) : SOURCE Collegamento foto originale

Input: immagine originale con distorsione fish-eye da correggere.

(2) : DESTINATION Link foto originale

Output: immagine corretta (tecnicamente anche con correzione prospettica, ma questo è un passaggio separato).

Come calcoli la distanza radiale dal centro per passare dal fisheye al rettilineo?

Il mio stub di funzioni ha questo aspetto:

 Point correct_fisheye(const Point& p,const Size& img) { // to polar const Point centre = {img.width/2,img.height/2}; const Point rel = {px-centre.x,py-centre.y}; const double theta = atan2(rel.y,rel.x); double R = sqrt((rel.x*rel.x)+(rel.y*rel.y)); // fisheye undistortion in here please //... change R ... // back to rectangular const Point ret = Point(centre.x+R*cos(theta),centre.y+R*sin(theta)); fprintf(stderr,"(%d,%d) in (%d,%d) = %f,%f = (%d,%d)\n",px,py,img.width,img.height,theta,R,ret.x,ret.y); return ret; } 

In alternativa, potrei in qualche modo convertire l’immagine da fisheye a rettilineo prima di trovare i punti, ma sono completamente confuso dalla documentazione di OpenCV . C’è un modo semplice per farlo in OpenCV, e funziona abbastanza bene da farlo per un feed video dal vivo?

La descrizione che si menziona afferma che la proiezione di una fotocamera pin-hole (una che non introduce la distorsione della lente) è modellata da

 R_u = f*tan(theta) 

e la proiezione delle comuni telecamere con objective fisheye (cioè distorta) è modellata da

 R_d = 2*f*sin(theta/2) 

Conoscete già R_d e theta e se conosceste la lunghezza focale della fotocamera (rappresentata da f), correggere l’immagine equivarrebbe a calcolare R_u in termini di R_d e theta. In altre parole,

 R_u = f*tan(2*asin(R_d/(2*f))) 

è la formula che stai cercando. La stima della lunghezza focale f può essere risolta calibrando la fotocamera o altri mezzi, ad esempio consentendo all’utente di fornire un feedback su come l’immagine viene corretta o utilizzando la conoscenza della scena originale.

Per risolvere lo stesso problema usando OpenCV, dovresti ottenere i parametri intrinseci della fotocamera e i coefficienti di distorsione dell’objective. Vedi, ad esempio, il capitolo 11 di Learning OpenCV (non dimenticare di controllare la correzione ). Quindi puoi usare un programma come questo (scritto con i collegamenti Python per OpenCV) per invertire la distorsione della lente:

 #!/usr/bin/python # ./undistort 0_0000.jpg 1367.451167 1367.451167 0 0 -0.246065 0.193617 -0.002004 -0.002056 import sys import cv def main(argv): if len(argv) < 10: print 'Usage: %s input-file fx fy cx cy k1 k2 p1 p2 output-file' % argv[0] sys.exit(-1) src = argv[1] fx, fy, cx, cy, k1, k2, p1, p2, output = argv[2:] intrinsics = cv.CreateMat(3, 3, cv.CV_64FC1) cv.Zero(intrinsics) intrinsics[0, 0] = float(fx) intrinsics[1, 1] = float(fy) intrinsics[2, 2] = 1.0 intrinsics[0, 2] = float(cx) intrinsics[1, 2] = float(cy) dist_coeffs = cv.CreateMat(1, 4, cv.CV_64FC1) cv.Zero(dist_coeffs) dist_coeffs[0, 0] = float(k1) dist_coeffs[0, 1] = float(k2) dist_coeffs[0, 2] = float(p1) dist_coeffs[0, 3] = float(p2) src = cv.LoadImage(src) dst = cv.CreateImage(cv.GetSize(src), src.depth, src.nChannels) mapx = cv.CreateImage(cv.GetSize(src), cv.IPL_DEPTH_32F, 1) mapy = cv.CreateImage(cv.GetSize(src), cv.IPL_DEPTH_32F, 1) cv.InitUndistortMap(intrinsics, dist_coeffs, mapx, mapy) cv.Remap(src, dst, mapx, mapy, cv.CV_INTER_LINEAR + cv.CV_WARP_FILL_OUTLIERS, cv.ScalarAll(0)) # cv.Undistort2(src, dst, intrinsics, dist_coeffs) cv.SaveImage(output, dst) if __name__ == '__main__': main(sys.argv) 

Si noti inoltre che OpenCV utilizza un modello di distorsione dell'objective molto diverso da quello della pagina Web a cui si è collegati.

(Poster originale, fornendo un’alternativa)

La seguente funzione mappa le coordinate di destinazione (rettilineo) alle coordinate sorgente (distorta fisheye). (Apprezzerei l’aiuto nel capovolgerlo)

Sono arrivato a questo punto per tentativi ed errori: non capisco fondamentalmente perché questo codice funzioni, le spiegazioni e l’accuratezza migliorata siano apprezzate !

 def dist(x,y): return sqrt(x*x+y*y) def correct_fisheye(src_size,dest_size,dx,dy,factor): """ returns a tuple of source coordinates (sx,sy) (note: values can be out of range)""" # convert dx,dy to relative coordinates rx, ry = dx-(dest_size[0]/2), dy-(dest_size[1]/2) # calc theta r = dist(rx,ry)/(dist(src_size[0],src_size[1])/factor) if 0==r: theta = 1.0 else: theta = atan(r)/r # back to absolute coordinates sx, sy = (src_size[0]/2)+theta*rx, (src_size[1]/2)+theta*ry # done return (int(round(sx)),int(round(sy))) 

Se utilizzato con un fattore di 3.0, riesce a distinguere con successo le immagini utilizzate come esempi (non ho fatto alcun tentativo di interpolazione di qualità):

Collegamento morto

(E questo è dal post del blog, per confronto 🙂

Utilizzando Panotools

Se pensi che le tue formule siano esatte, puoi calcolare una formula esatta con trig, in questo modo:

 Rin = 2 f sin(w/2) -> sin(w/2)= Rin/2f Rout= f tan(w) -> tan(w)= Rout/f (Rin/2f)^2 = [sin(w/2)]^2 = (1 - cos(w))/2 -> cos(w) = 1 - 2(Rin/2f)^2 (Rout/f)^2 = [tan(w)]^2 = 1/[cos(w)]^2 - 1 -> (Rout/f)^2 = 1/(1-2[Rin/2f]^2)^2 - 1 

Tuttavia, come dice @jmbr, la distorsione effettiva della fotocamera dipenderà dall’objective e dallo zoom. Piuttosto che fare affidamento su una formula fissa, potresti provare una espansione polinomiale:

 Rout = Rin*(1 + A*Rin^2 + B*Rin^4 + ...) 

Modificando prima i coefficienti A, quindi quelli di ordine superiore, è ansible calcolare qualsiasi funzione locale ragionevole (la forma dell’espansione sfrutta la simmetria del problema). In particolare, dovrebbe essere ansible calcolare i coefficienti iniziali per approssimare la funzione teorica di cui sopra.

Inoltre, per ottenere buoni risultati, sarà necessario utilizzare un filtro di interpolazione per generare l’immagine corretta. Finché la distorsione non è troppo grande, puoi usare il tipo di filtro che useresti per ridimensionare l’immagine in modo lineare senza troppi problemi.

Modifica: secondo la tua richiesta, il fattore di scala equivalente per la formula precedente:

 (Rout/f)^2 = 1/(1-2[Rin/2f]^2)^2 - 1 -> Rout/f = [Rin/f] * sqrt(1-[Rin/f]^2/4)/(1-[Rin/f]^2/2) 

Se tracciate la formula sopra insieme a tan (Rin / f), potete vedere che sono molto simili nella forma. Fondamentalmente, la distorsione dalla tangente diventa grave prima che il peccato (w) diventi molto diverso da w.

La formula inversa dovrebbe essere qualcosa del tipo:

 Rin/f = [Rout/f] / sqrt( sqrt(([Rout/f]^2+1) * (sqrt([Rout/f]^2+1) + 1) / 2 ) 

Ho preso ciò che ha fatto JMBR e sostanzialmente lo ho invertito. Ha preso il raggio dell’immagine distorta (Rd, ovvero la distanza in pixel dal centro dell’immagine) e ha trovato una formula per Ru, il raggio dell’immagine non distorta.

Vuoi andare dall’altra parte. Per ogni pixel non distorto (immagine elaborata), si desidera sapere qual è il pixel corrispondente nell’immagine distorta. In altre parole, dato (xu, yu) -> (xd, yd). Quindi sostituisci ogni pixel nell’immagine non distorta con il pixel corrispondente dall’immagine distorta.

Iniziando dove ha fatto JMBR, faccio il contrario, trovando Rd come una funzione di Ru. Ottengo:

 Rd = f * sqrt(2) * sqrt( 1 - 1/sqrt(r^2 +1)) 

dove f è la lunghezza focale in pixel (ti spiegherò più avanti), e r = Ru/f .

La lunghezza focale della mia fotocamera era di 2,5 mm. La dimensione di ciascun pixel sul mio CCD era di 6 um quadrati. f era quindi 2500/6 = 417 pixel. Questo può essere trovato per tentativi ed errori.

Trovare Rd consente di trovare il pixel corrispondente nell’immagine distorta usando le coordinate polari.

L’angolo di ogni pixel dal punto centrale è lo stesso:

theta = arctan( (yu-yc)/(xu-xc) ) dove xc, yc sono i punti centrali.

Poi,

 xd = Rd * cos(theta) + xc yd = Rd * sin(theta) + yc 

Assicurati di sapere in quale quadrante ci sei.

Ecco il codice C # che ho usato

  public class Analyzer { private ArrayList mFisheyeCorrect; private int mFELimit = 1500; private double mScaleFESize = 0.9; public Analyzer() { //A lookup table so we don't have to calculate Rdistorted over and over //The values will be multiplied by focal length in pixels to //get the Rdistorted mFisheyeCorrect = new ArrayList(mFELimit); //i corresponds to Rundist/focalLengthInPixels * 1000 (to get integers) for (int i = 0; i < mFELimit; i++) { double result = Math.Sqrt(1 - 1 / Math.Sqrt(1.0 + (double)i * i / 1000000.0)) * 1.4142136; mFisheyeCorrect.Add(result); } } public Bitmap RemoveFisheye(ref Bitmap aImage, double aFocalLinPixels) { Bitmap correctedImage = new Bitmap(aImage.Width, aImage.Height); //The center points of the image double xc = aImage.Width / 2.0; double yc = aImage.Height / 2.0; Boolean xpos, ypos; //Move through the pixels in the corrected image; //set to corresponding pixels in distorted image for (int i = 0; i < correctedImage.Width; i++) { for (int j = 0; j < correctedImage.Height; j++) { //which quadrant are we in? xpos = i > xc; ypos = j > yc; //Find the distance from the center double xdif = i-xc; double ydif = j-yc; //The distance squared double Rusquare = xdif * xdif + ydif * ydif; //the angle from the center double theta = Math.Atan2(ydif, xdif); //find index for lookup table int index = (int)(Math.Sqrt(Rusquare) / aFocalLinPixels * 1000); if (index >= mFELimit) index = mFELimit - 1; //calculated Rdistorted double Rd = aFocalLinPixels * (double)mFisheyeCorrect[index] /mScaleFESize; //calculate x and y distances double xdelta = Math.Abs(Rd*Math.Cos(theta)); double ydelta = Math.Abs(Rd * Math.Sin(theta)); //convert to pixel coordinates int xd = (int)(xc + (xpos ? xdelta : -xdelta)); int yd = (int)(yc + (ypos ? ydelta : -ydelta)); xd = Math.Max(0, Math.Min(xd, aImage.Width-1)); yd = Math.Max(0, Math.Min(yd, aImage.Height-1)); //set the corrected pixel value from the distorted image correctedImage.SetPixel(i, j, aImage.GetPixel(xd, yd)); } } return correctedImage; } } 

Ho trovato questo file pdf e ho dimostrato che la matematica è corretta (eccetto per la riga vd = *xd**fv+v0 which should say vd = **yd**+fv+v0 ).

http://perception.inrialpes.fr/CAVA_Dataset/Site/files/Calibration_OpenCV.pdf

Non utilizza tutti i coefficienti più recenti che OpenCV ha a disposizione, ma sono sicuro che potrebbe essere adattato abbastanza facilmente.

 double k1 = cameraIntrinsic.distortion[0]; double k2 = cameraIntrinsic.distortion[1]; double p1 = cameraIntrinsic.distortion[2]; double p2 = cameraIntrinsic.distortion[3]; double k3 = cameraIntrinsic.distortion[4]; double fu = cameraIntrinsic.focalLength[0]; double fv = cameraIntrinsic.focalLength[1]; double u0 = cameraIntrinsic.principalPoint[0]; double v0 = cameraIntrinsic.principalPoint[1]; double u, v; u = thisPoint->x; // the undistorted point v = thisPoint->y; double x = ( u - u0 )/fu; double y = ( v - v0 )/fv; double r2 = (x*x) + (y*y); double r4 = r2*r2; double cDist = 1 + (k1*r2) + (k2*r4); double xr = x*cDist; double yr = y*cDist; double a1 = 2*x*y; double a2 = r2 + (2*(x*x)); double a3 = r2 + (2*(y*y)); double dx = (a1*p1) + (a2*p2); double dy = (a3*p1) + (a1*p2); double xd = xr + dx; double yd = yr + dy; double ud = (xd*fu) + u0; double vd = (yd*fv) + v0; thisPoint->x = ud; // the distorted point thisPoint->y = vd; 

Ho implementato ciecamente le formule da qui , quindi non posso garantire che farebbe ciò di cui hai bisogno.

Usa auto_zoom per ottenere il valore per il parametro di zoom .

 def dist(x,y): return sqrt(x*x+y*y) def fisheye_to_rectilinear(src_size,dest_size,sx,sy,crop_factor,zoom): """ returns a tuple of dest coordinates (dx,dy) (note: values can be out of range) crop_factor is ratio of sphere diameter to diagonal of the source image""" # convert sx,sy to relative coordinates rx, ry = sx-(src_size[0]/2), sy-(src_size[1]/2) r = dist(rx,ry) # focal distance = radius of the sphere pi = 3.1415926535 f = dist(src_size[0],src_size[1])*factor/pi # calc theta 1) linear mapping (older Nikon) theta = r / f # calc theta 2) nonlinear mapping # theta = asin ( r / ( 2 * f ) ) * 2 # calc new radius nr = tan(theta) * zoom # back to absolute coordinates dx, dy = (dest_size[0]/2)+rx/r*nr, (dest_size[1]/2)+ry/r*nr # done return (int(round(dx)),int(round(dy))) def fisheye_auto_zoom(src_size,dest_size,crop_factor): """ calculate zoom such that left edge of source image matches left edge of dest image """ # Try to see what happens with zoom=1 dx, dy = fisheye_to_rectilinear(src_size, dest_size, 0, src_size[1]/2, crop_factor, 1) # Calculate zoom so the result is what we wanted obtained_r = dest_size[0]/2 - dx required_r = dest_size[0]/2 zoom = required_r / obtained_r return zoom