la correction de la déformation du fisheye par programmation

La description que vous mentionnez états que la projection par un trou d'épingle de la caméra (qui n'introduit pas de distorsion de l'objectif) est modélisé par

R_u = f*tan(theta)

et la projection d'un commun objectif fisheye caméras (qui est, déformé) est modélisé par

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

Vous le savez déjà R_d et thêta et si vous connaissais l'appareil photo de la longueur focale (représenté par f) la correction de l'image serait à l'informatique R_u en termes de R_d et thêta. En d'autres termes,

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

est la formule que vous cherchez. L'estimation de la distance focale f peut être résolu par l'étalonnage de la caméra ou par d'autres moyens tels que la possibilité pour l'utilisateur de fournir de la rétroaction sur la façon dont l'image est corrigée ou à l'aide de connaissances à partir de la scène d'origine.

Afin de résoudre le même problème en utilisant OpenCV, vous devez obtenir l'appareil photo, les paramètres intrinsèques et de la distorsion de l'objectif coefficients. Voir, par exemple, le Chapitre 11 de L'Apprentissage D'OpenCV (n'oubliez pas de vérifier la correction). Ensuite, vous pouvez utiliser un programme comme celui-ci (écrit avec les liaisons Python pour OpenCV) dans le but d'inverser la distorsion de l'objectif:

#!/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)

Également noter que OpenCV utilise très différente de la distorsion de l'objectif du modèle à la une de la page web lié à.

Si vous pensez que vos formules sont exactes, vous pouvez comput une formule exacte avec la trigo, comme suit:

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

Cependant, comme @jmbr dit, la caméra de distorsion dépend de l'objectif et le zoom. Plutôt que de compter sur une formule fixe, vous voudrez peut-être essayer un polynôme d'extension:

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

En jouant d'abord Un, puis de l'ordre supérieur des coefficients, vous pouvez calculer raisonnable de fonction locale (la forme de l'expansion prend avantage de la symétrie du problème). En particulier, il devrait être possible de calculer initiale coefficients de l'approximation de l'théorique de la fonction ci-dessus.

Aussi, pour obtenir de bons résultats, vous aurez besoin d'utiliser un filtre d'interpolation pour générer votre image corrigée. Aussi longtemps que la distorsion n'est pas trop grande, vous pouvez utiliser le type de filtre que vous voudriez utiliser pour redimensionner l'image de façon linéaire sans trop de problème.

Edit: selon votre demande, l'équivalent facteur d'échelle pour la formule ci-dessus:

(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)

Si vous tracez la formule ci-dessus aux côtés de tan(Rin/f), vous pouvez voir qu'ils sont très similaires dans leur forme. Fondamentalement, la distorsion de la tangente devient sévère avant le péché(w) devient très différent de w.

La formule d'inversion doit être quelque chose comme:

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

J'ai pris ce que JMBR fait et fondamentalement renversée. Il prit le rayon de la déformation de l'image (Rd, qui est, la distance en pixels depuis le centre de l'image) et a trouvé une formule pour le Ru, le rayon de la distortion de l'image.

Vous voulez aller dans l'autre sens. Pour chaque pixel de la non faussée (image traitée), vous voulez savoir ce que le pixel correspondant est dans la déformation de l'image.
En d'autres termes, compte tenu de (xu, yu) --> (xd, yd). Vous replacez ensuite chaque pixel dans une image avec son pixel correspondant de l'image déformée.

À partir du moment où JMBR fait, je fais l'inverse, de trouver des Rd en fonction de Ru. J'obtiens:

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

où f est la longueur focale en pixels (j'expliquerai plus tard), et r = Ru/f.

La longueur focale de mon appareil photo était de 2,5 mm. La taille de chaque pixel sur mon CCD de 6 um carré. f est donc 2500/6 = 417 pixels. Ce qui peut être trouvé par essai et erreur.

Trouver Rd vous permet de trouver le pixel correspondant dans l'image déformée à l'aide des coordonnées polaires.

L'angle de chaque pixel à partir du point central est le même:

theta = arctan( (yu-yc)/(xu-xc) ) où xc, yc sont les points centraux de.

Puis,

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

Assurez-vous que vous savez quel quadrant vous êtes dans.

Voici le code C#, j'ai utilisé

 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;
      }
}