Lens Correction in PanoTools

Lens Correction in PanoTools

The PanoTools library implements an effective, but rather idiosyncratic method for correcting lens projections, that causes a good deal of puzzlement. Lens parameters optimized for one image format generally do not work for a different format; even rotating a set of images 90 degrees before aligning them produces different and incompatible lens parameters. One would expect that there must be a way to convert either of those parameter sets to a common form, that would apply equally well to both formats, or indeed to any image taken with the same lens. To see how that might be done, I have made a detailed analysis of PanoTools lens correction computations, based on the code in historic as well as current versions of libpano and helpful discussions with Helmut Dersch.

Why Lens Correction?

To make a panoramic image from photographs, it is essential to be able to calculate the direction in space corresponding to any given position in a given image. specifically, we need to know the angles between the view directions of the photos (the alignment of the images), and a radial projection function that relates the distance of a point from image center to the true angle of view, measured from the optical axis of the lens. Given a set of control points linking the images, PanoTools estimates both the alignment and the lens projection by a nonlinear least squares fitting procedure -- optimization. Using the fitted lens parameters, the stitcher can correct each image to match the ideal geometry of the scene, according to whatever projection is chosen for the panorama. Done right, that makes all the images fit together perfectly; moreover, it yields a panoramic image that seems to have been made with a perfect lens.

Mapping View Angle <=> Radius

The radial projection of a real lens may approximate some known mathematical function, but in practice it must be determined experimentally, a process known as calibrating the lens. A calibration is a parametrized mathematical model, fitted to experimental data. The typical model consists of an ideal angle-to-radius function, and the coefficients of a polynomial that converts the ideal radius to the actual radius measured on the image.

Like many lens calibration programs, libpano uses just two model projection functions: rectilinear, for 'normal' lenses, and 'fisheye', for all others. The rectilinear projection has radius proportional to the tangent of the view angle. PT's 'fisheye', better known as the equal-angle spherical projection, has radius proportional to angle. The constant of proportionality is the lens focal length, F. With angle A in radians, and R the ideal radius, the formulas are

Rectilinear: R / F = tan( A )

Equal-angle: R / F = A

Of course R and F have to be measured in the same units. If we have F in mm, then R is in mm also. If we want to measure R in pixels, then we need F in pixels; for example

F (pixels) = (FL in mm) * (image width in pixels) / (sensor width in mm).

In any case, F is the constant of proportionality between the actual radius and the value of a trigonometric function that defines the basic shape of the projection.

In physical optics, focal length is defined as the first derivative of R by A, at A = 0. That is easy to see if we write R = F A or R = F tan(A), because the slopes of A and tan(A) are both 1 at A = 0. This is also true of other trigonometric functions commonly used as ideal lens projections:

Equal-Area: R / F = 2 sin( A / 2 )

Stereographic: R / F = 2 tan( A / 2 ).

The dimensionless quantity N = R / F is the normalized ideal radius. Multiplying N by the focal length, in any units, gives the ideal image radius in the same units.

Generic Correction Scheme

The difference between the real lens projection and the ideal one is modeled by an adjustable correction function that gives the observed radius as a function of the ideal radius. The adjustable part is almost always a polynomial, because it it easy to fit polynomials to experimental data. The argument to the polynomial is usually the normalized ideal radius,

N = R / F,

because that makes the polynomial coefficients independent of how image size is measured. The constant term is 0 because both radii are zero at the same point. If the coefficient of the linear term is 1, so that the first derivative at 0 is 1, then the value of the polynomial will be the normalized observed radius, n = r / F.

Many calibration packages use a polynomial with only even order terms beyond the first:

n = N + a N ^2 + b N ^4 + c N ^6.

Equivalently

n = N ( 1 + a N + b N^3 + c N^5 )

The expression in parentheses is the ratio of observed to ideal radius, which is expected to be close to 1 everywhere if the ideal model function is well chosen.

Multiplying n by the focal length, in any units, gives the observed image radius in the same units:

r = F n.

PanoTools Correction Scheme

Lens correction in PanoTools is unusual in several respects. First, it ignores the physical parameters of the lens (focal length) and camera (pixel size). Instead, it computes angle-to-radius scale factors from image dimensions and fields of view, as described below. All correction computations are in terms of image radii, measured in pixels, rather than the normalized radii described above. However, normalized radii are evaluated implicitly.

Second, the correction is computed in ideal equal-angle spherical coordinates, rather than camera coordinates. Observed image points are found by remapping their ideal coordinates to camera coordinates, according to the ideal lens projection, and rescaling them according to the ratio of pixel sizes in the source and ideal images.

Third, the adjustable polynomial (which gives the ratio of observed to ideal radius) is normalized to hold a certain radius constant:

r = R poly( R/r0 ),

where r0 is the constant radius. The polynomial is normalized so that its value is 1 when its argument is 1 (that is, when R = r0). With

X = R / r0

The correction polynomial is

x = (1 - a - b - c) X + a X^2 + b X^3 + c X^4,

and the observed radius is given by

r = R x.

The observed radius is thus a 5th order polynomial in R with first order coefficient zero:

r = s R^2 + t R^3 + u R^4 + v R^5,

where s = (1-a-b-c) / r0, t = a / r0^2, u = b / r0^3, v = c / r0^4.

The normalization makes this relation stable under optimization of a, b, and c. It is also essential to the correctness of the result, which can be seen as follows. The ideal radius is

R = F N

where F is the ideal focal length in pixels. We can thus write the adjusted radius as

r = F A poly( F A / r0 ),

The normalized observed radius n = r / F, so

n = A poly( F A / r0 ).

Since n is a dimensionless number, this can only be true if the argument to poly() is also a dimensionless number. F has the dimension of pixels, therefore it must be divided by another parameter dimensioned in pixels.

The overall computation proceeds as follows. PanoTools computes the ideal radius R by mapping a point in the panorama (which plays the role of the ideal image) to equal angle spherical projection. Then R = sqrt( h^2 + v^2 ), where h and v are the pixel coordinates relative to the center of the equal-angle projection. Then PT's radius() function computes x as described, and returns scaled coordinates ( h x, v x ). If the lens is rectilinear, PT next remaps those coordinates to rectilinear; if it is a fisheye, no remapping is needed. In either case the coordinates are finally rescaled to account for any difference in resolution between the panorama and the source image. The scale factor is computed from the dimensions and angular fields of view of the panorama and the source image, as follows.

d = (half width of pano) / A2Npano( half hfov of pano ),

e = (half width of source) / A2Nsource( half hfov of source ),

where A2Npano and A2Nsource are the ideal functions for panorama and lens. Then

Rsource = Rpano e / d.

The scale factors d and e are focal lengths in pixels, because A2N() yields the normalized radius, equal to R / F. For the panorama, which follows an ideal projection, d is identical to F. In fact d, under the name “distance factor”, is used in many of libpano's coordinate transformation functions to convert radius in pixels to the ideal normalized radius in trigonometric units.

For the source image, whose true projection is only approximately known, e is an estimate of F according to the fitted correction parameters. Since hfov is one of those parameters, the fitted value of e will be proportional to the true F; the constant of proportionality will approach 1 as the fitted polynomial coefficients approach 0.

Portable Lens Parameters

The focal length in pixels, or something proportional to it, must be known at a certain stage of lens calibration. This quantity clearly depends on both lens and camera properties. In most cases today, equipment manufacturers' specifications can provide the needed data:

F = (FL in mm) * (image width in pixels) / (sensor width in mm).

This is normally close to the true value, the main uncertainty being how accurately the nominal focal length reflects the true one. Computing F this way makes it possible for the fitted correction coefficients to be independent of camera format, as explained above.

In any practical calibration scheme F is also an adjustable parameter. However the data used to determine F are often different from the data used to fit the lens curve. PanoTools is somewhat unique in fitting all lens parameters to one set of experimental values.

PanoTools does not use physical lens or camera specifications. Instead, it takes image dimensions and angular field of view (fov) as the primary data, and determines the effective value of F implicitly during the parameter fitting process. For an ideal equal-angle spherical projection, it is true that

F = (image width in pixels) / (horizontal fov in radians).

This formula is not correct for any other projection, but for purposes of lens calibration it can be considered a useful estimate for initializing an optimization procedure. PT actually adjusts hfov (along with a, b, and c) during optimization, rather than F.

The normalization of the radial polynomial makes it possible to equate the correction computed by PT to the generic one. As discussed above, that is true for any choice of r0 that is proportional to image size. The best choice would be

r0 = Fpano = “distance parameter” d, defined above.

That would make the argument of the polynomial

Rpano / Fpano

which is equal to the ideal normalized radius, N. Then the polynomial coefficients fitted by PanoTools would be independent of image dimensions, and could be used directly in a generic function to compute the radius correction in camera coordinates. The value actually used is merely proportional to Fpano; the proportionality factor is a function of image dimensions, so the correction coefficients depend on the image format and are not portable.

To convert existing PT lens parameters to portable form requires two things: computing portable radial coefficients, and defining F in a way that is not tied to an image format. The latter is easy, but does require at least one physical parameter, not in PT's parameter set. Ideally, at least when dealing with a digital camera, that would be the physical size of a pixel in the images used for calibration. With w that pixel width in mm,

F(mm) = w F(pixels) = w (scale factor e, defined above)

This could be calculated by front-end software using the scale factor e from libpano.

The coefficients can be converted using data available inside libpano. With

k = d / r0,

a' = a k^2,

b' = b k^3,

c' = c k^4

These are the coefficients of a polynomial in the normalized ideal radius N, with first order term

w' = 1 – a' – b' – c'

that computes the same radius correction factor as the PT polynomial.

Along with the ideal mapping function A2Nsource(), these coefficients constitute a portable radial projection function that would be very useful even without a portable version of PT's focal length estimate. Focal length is just a linear scale factor, that can be fitted or estimated in many different ways according to the problem at hand, while the radial function represents the essential 'lens curve' independently of image scale.

-- 18 Jan 2010 TKSharpless