A REVIEW OF THE ONE-PARAMETER DIVISION UNDISTORTION MODEL

The one-parameter division undistortion model by (Lenz, 1987) and (Fitzgibbon, 2001) is a simple radial distortion model with beneficial algebraic properties that allows to reason about some problems analytically that can only be handled numerically in other distortion models. One property of this distortion model is that straight lines in the undistorted image correspond to circles in the distorted image. These circles are fully described by their center point, as the radius can be calculated from the position of the center and the distortion parameter only. This publication collects the properties of this distortion model from several sources and reviews them. Moreover, we show in this publication that the space of this center is projectively isomorphic to the dual space of the undistorted image plane, i.e. its line space. Therefore, projective invariant measurements on the undistorted lines are possible by the according measurements on the centers of the distorted circles. As an example of application, we use this to find the metric distance of two parallel straight rails with known track gauge in a single uncalibrated camera image with significant radial distortion.


INTRODUCTION
A main objective of photogrammetry is to obtain metric information about physical objects from images thereof. The more accurate the image generation process can be modeled, the more accurate results can be expected. To achieve this, sophisticated camera models were devised. E.g. OpenCV (Bradski, 2000) uses a camera model with up to 18 parameters for camera calibration and Metashape (Agisoft LLC, 2021) uses still up to 12 parameters for the frame camera model. These models can lead to very accurate results for laboratory calibrations of measuring cameras or for auto-calibration during bundle adjustment in structure from motion applications, if sufficiently many well distributed feature correspondences in enough frames are available. However, sometimes, a laboratory calibration is not possible, e.g. because the camera that took an image does not exist anymore or is not at reach and there are not enough images or not enough precise feature correspondences to estimate that many parameters reliably. In this case, simpler models are advantageous as they allow for simpler and more robust calibration.
As an example, in the Building Rome in a Day project, (Agarwal et al., 2009) reconstructed popular sights of Rome from tourist photos found publicly available on the internet. Since each photo was basically taken by an individual camera, a full calibration of each camera would not have been feasible. Instead, the photos were prealigned using only the values in the EXIF metadata as a rough estimate of the focal length, which was then refined during bundle adjustment together with only two more parameters for the radial distortion. This simple three parameter camera model was still sufficiently accurate to generate an impressive large scale sparse 3D reconstruction.
In this article, we discuss an even simpler camera model requiring only two parameters, one for the adjusted focal length and a single one for radial distortion. The distortion model was * bastian.erdnuess@kit.edu introduced by (Lenz, 1987) and (Lenz and Fritsch, 1990) and further spread by (Fitzgibbon, 2001) who called it the division model. Here, we'll call it more specifically the one-parameter division undistortion model in order to emphasize that the division formula models the undistortion process instead of the distortion process. This model offers some unique algebraic properties, like being analytically invertible and distorting lines to perfect circles, that can be (and have been) beneficially used in algorithms. This article reviews some known properties of the (un)distortion model, devises some new properties, and shows relations to other theories.
As two main results, this article shows, (1) that for cameras with tissue distortion, the relation between distorted and undistorted coordinates is deeply related to the transformation between the Poincaré and Beltrami-Klein disks in hyperbolic geometry, cf. Sec. 4.6, and (2) that the line space of the undistorted image is projective isomorph to the space of the center points of the circles that correspond to the lines in the distorted image, cf. Sec. 4.9. The second result is illustrated in an application where the distance of two parallel straight railroad tracks is extracted from a single distorted image by only using the center points of the circles that contain the distorted rails and the known track gauge of the rails, cf. Sec. 4.10. (Lenz, 1987) motivates the one-parameter division undistortion model by it's algebraic properties that allow to solve some algorithms faster and direct that would usually have to be solved iteratively. However, he also noted, that adding parameters for higher order distortion didn't improve the results significantly for the investigated CCD sensors but made the calibration procedure more unstable. Later, (Fitzgibbon, 2001) compares the model to the more common one-parameter polynomial distortion model and finds them to perform almost equal, with the division model slightly ahead. About at the same time, (Geyer and Daniilidis, 2001a, Geyer and Daniilidis, 2001b, Geyer and Daniilidis, 2002 studied catadioptric cameras with parabolic mirrors and telecentric cameras which lead to an algebraically equivalent distortion of the images. In their papers, they present deep insights into the algebraic properties of the distortion model for the case of barrel distortion. Later, at the same institute, (Barreto andDaniilidis, 2005, Barreto, 2006) used the algebraic properties of the model to introduce a lifting method and used it to define a (so called) radial variant of the fundamental matrix, that also includes the distortion parameter. Furthermore, they adopted the algorithm to solve for the fundamental matrix accordingly, which led to a new direct (i.e. non-iterative) algorithm for solving for the epipolar geometry and the (un)distortion parameter simultaneously, cf. Sec. 4.5. (Brito et al., 2012, Brito et al., 2013b, Brito et al., 2013a, Brito, 2017 applied this lifting technique to greater extent to further improve direct algorithms for simultaneous reconstruction and camera calibration. They could significantly improve the performance of algorithms that usually assume a pinhole camera model without distortion.

CAMERA MODEL
We assume a simple camera model that focuses on the lens distortion, see Fig. 1. We assume a digital frame sensor with W ×H pixels and introduce normalized coordinates x with origin in the center of the sensor such that the sensor fits exactly in the unit disk. I.e. pixel P = (R, C) at row R and column C (counting from zero) gets mapped to x = (x, y), where Here, we assume a rectangular grid of square pixels. If this is not the case (e.g. pixels are not square, like e.g. on DVCAM), these formulas have to be adjusted accordingly in order to make the x-y-axis orthogonal and of equal scale.
We assume that the sensor sits (rotated upside-down) at distance c behind a lens. Opposite of the lens at distance c, we introduce a (virtual) scene plane. We assume the whole scene to be projected through the center of the lens onto the scene plane and consider the mapping between the scene plane and the sensor plane to describe the lens distortion.
A scene point X = (X, Y, Z) gets projected onto the point ξ = (ξ, η) on the scene plane, where That point sits at (ξ, η, c) = (ξ, c) = c Z X in the global coordinate system. It gets refracted at the lens in the origin and mapped to the point x = (x, y) in the sensor plane which is located at (−x, −y, −c) = −(x, c) in the global coordinate system. Note that the sensor plane's x-y-coordinate axes are anti-parallel to the world's X-Y-axis, or in other words, they are rotated by 180 • around the Z-axis. The distortion is the map ξ → x that describes, how a point (ξ, c) of the scene plane is mapped to a point −(x, c) on the sensor plane and the choice of the rotated x-y-coordinate system for the sensor ensures that the identity map x = ξ corresponds to the case of no distortion.
x x x x x y y y y y y y y y y y y y y y y y Figure 1. Camera model showing (1) the sensor with integer R-C-pixel grid and centered x-y-coordinate system, (2) the lens and world X-Y-Z-coordinate system, (3) the scene plane (red) and sensor plane (blue, rotated upside down) with their ξ-ηand x-y-coordinate systems, resp., and (4) a 3D scene point X, it's projection ξ onto the scene plane and x onto the sensor plane with distortion ϕ to ϕ . Notes: The world, sensor, and scene coordinate systems are scaled such that the sensor has a diagonal of 2 units length, i.e. the whole sensor fits exactly in the unit disk.
x is located on the plane containing X and the world's Z-axis. The sensor, as shown on the left hand side, is illuminated from the backside. The world's X-axis and scene plane's ξ-axis point towards the reader and the sensor plane's x-axis points opposite; these axes are omitted in the image. The distance c between the lens and the sensor plane is the adjusted focal length (relative to the sensor size) and determines the field of view of the camera, the distance c between the lens and the scene plane is chosen for symmetry.

Distortion Model
Due to optical effects the image gets distorted while projecting from the scene plane to the sensor plane. As suggested by (Lenz, 1987) and (Fitzgibbon, 2001) we assume the simple oneparameter division undistortion model given by the formula where x 2 = x 2 denotes the squared norm of a vector and α is the distortion parameter. In opposite to more complex distortion models like e.g. the distortion model of (Brown, 1971) this model only handles radial distortion and uses only one parameter for it. Therefore, it is not as expressive and can't model the imaging process as accurately, but on the other hand the parameter is much easier and often more stable to estimate.
Note that since Eq. (3) gives the undistorted coordinates ξ in terms of the distorted coordinates x, it would actually better be referred to as an undistortion model, i.e. a map x → ξ that retrieves the undistorted point (ξ, c) on the scene plane given the distorted point −(x, c) on the sensor plane.
Taking the norm on both sides of Eq.
(3) and rearranging leads to This fraction on the right hand side describes the relative amount a point ξ in the scene plane gets distorted by the lens. By construction of Eq. (1) the very corners of the image have coordinates x with x 2 = 1 and substituting x 2 on the left hand side shows that α corresponds to the distortion in the image corners, i.e. the maximal distortion for this distortion model. For images with tissue distortion, α is positive and for images with barrel distortion, α is negative. For undistorted images, α is zero. In the following, we assume α = 0 without further mention. α is sometimes denoted as percentage, i.e. a barrel distortion of 5% corresponds to the distortion parameter α = −5% = −0.05.

The Adjusted Focal Length c
We choose to scale the sensor to fit precisely in the unit disk in order to interpret the parameter α as maximal distortion. The second parameter of the camera model is the adjusted focal length c which approximately corresponds to the focal length (Förstner and Wrobel, 2016, p. 461) of the lens measured in multiples of the half sensor diagonal.
If the adjusted focal length is known in advance, another choice would have been to scale Fig. 1 such that c = 1. This would affect the distortion parameter which would bê in the scaled model, where α and c are from the unscaled model. If c is unknown, α can still be determined first from e.g. some known lines in the image as discussed below. However, in the undistorted image the adjusted focal length c has then still to be determined if necessary for the application. Sometimes, it is more natural to determine c first and then α orα, but sometimes it is beneficial to determine α first, e.g. when c is not needed because one only needs to extract some projective invariant measures from the image in which c would cancel out anyways.
The scene plane is chosen to be also at distance c of the lens in order to ensure that the distorted and undistorted image almost agree around the image center. Another choice could have been to place it at distance 2c. This would make the undistorted image twice as large around the image center compared to the distorted image. It would simplify some formulas at the expense of some minor complications in other formulas by changing some factors. It would also increase the similarity to the transformation between the Poincaré and Beltrami-Klein models of hyperbolic geometry discussed in Sec. 4.6, as well as to the discussion of the parabolic catadioptric camera in (Geyer and Daniilidis, 2001b).

The Inverse
This model is algebraically invertible and it's inverse is given by cf. (Lenz, 1987). This is actually the corresponding distortion model ξ → x.

Lines Distort to Circles
Under this distortion model, a line in the scene plane gets projected to a circle or a line on the sensor plane, cf. (Geyer and Daniilidis, 2001b, Fig. 1). More specific, lines through the origin are unaffected by the distortion (as a whole, which holds in general for radial distortion) and lines not going to the origin are projected onto circles.
A line on the scene plane not passing through the origin can be parameterized by a vector p such that for all ξ in the scene plane it holds p, ξ = 1 , where ·, · denotes the dot product. Substituting Eq.
(3) and rearranging the terms for x and x 2 leads after completing the square to With the substitutions this is the equation of a circle with center m and radius r. Such circle will be denoted as C(m, r 2 ) in the following. When the arc of a circle C(m, r 2 ) has been extracted from an image, it belongs to a straight line in the scene plane precisely if m 2 − r 2 = 1 α . In this case, the circle corresponds to the line in the scene plane consisting of all ξ for which This can be seen by dividing Eq. (7) by 2α and substituting Eq. (9) in the result.
If the distortion parameter α is unknown but several arcs from circles C(m k , r 2 k ) have been extracted from a photo of a scene with linear structures, each arc votes for the distortion para- If the α k cluster around a value, this value would be a good guess for the distortion parameter α.

Transformation Under Scaling
By scaling x and ξ tō x andξ are related bȳ where the sign in front of thex 2 is that of α and opposite for the sign in front ofξ 2 . Therefore, up to scale, it is often sufficient to consider the two cases α = ±1.

Interpretation in Homogeneous Coordinates
Both, the undistortion formula in Eq.
(3) as well as Eq. (6) can be understood as dehomogenization of an homogeneous vector, see e.g. (Förstner and Wrobel, 2016, p. 199) for homogeneous coordinates. I.e. the following vectors are multiples of each other: In the first case, the point on the right hand side is on a paraboloid. In the second case its on an ellipsoid for positive α or on a hyperboloid sheet for negative α. This allows for a different interpretation of the mapping process, see Fig. 2 and Fig. 3 for a discussion of the cases α = ±1.
x ξ z −1 −0.5 0.5 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 If α is positive but not 1, the paraboloid in Fig. 2 is more broad or narrow and the sphere turns into a spheroid, still of height 1. When α approaches 0, the paraboloid and ellipsoid become more and more flat until finally, for α = 0 the upper half ellipsoid and the paraboloid both agree with the plane at z = 1. When α turns negative, the paraboloid flips over and bends down and the ellipsoid breaks up to a hyperbola with the upper half bending up and the lower half going through the origin and bending down. Fig. 4 illustrates this for several values of α.
In a similar way, one can also define a radial trifocal tensor or a radial homography matrix (Brito et al., 2013b, Kukelova et al., 2015.  Fig. 2 shows that for α = 1 the distortion direction ξ → x can be decomposed in a parallel projection from the plane to the sphere followed by a central projection back to the plane. As the central projection has its origin in the south pole of the sphere, this happens to be a stereographic projection. The composition of parallel and stereographic projection is well known from hyperbolic geometry and transforms the Beltrami-Klein disk into the Poincaré disk. With the pointx in the Poincaré disk and the according pointξ in the Beltrami-Klein disk, the transformation formula between them is (Richter-Gebert, 2011, p. 480

Relation to Hyperbolic Geometry
This formula agrees with Eq. (12) for α = 1 up to the substitutionξ = 2ξ. Therefore, up to scale, undistorting a tissue distorted image with the one-parameter division undistortion model is the same as the transformation from the Poincaré disk to the Beltrami-Klein disk, see also Fig. 5.
A C B Figure 5. Difference between the Beltrami-Klein disk and the undistorted camera image for α = 1. The line A gets distorted to circle B which is close to A near the origin. This ensures that the distorted and undistorted images almost agree at the image center. Circle B in the Poincaré disk corresponds to line C in the Beltrami-Klein disk which both intersect at the disk border. Line C is twice as long and twice as far away from the origin as line A. Instead of being close at near the origin, C shares other properties with B, e.g. the polar of C with respect to the disk border is the center of the circle B, see Fig. 6.
For α > 0, a circle C(m, r 2 ) is a distorted line, precisely when it intersects Cα orthogonally, cf. (Richter-Gebert, 2011, p. 480). This is a corollary from the observation that the distorted image can be interpreted as Poincaré disk. See also Fig. 6.
See also Fig. 6. In the case α < 0, Cα is the distorted image of the plane at infinity Z = 0 in the global coordinate system (cf. Fig. 1). Cα

Poles and Polars
Given a circle K = C(m, r 2 ) and a point p = m, the line L consisting of all points l with is the polar of p (with respect to K) and p is the pole of L, cf. (Koecher and Krieg, 2000, IV.1.4).
If p is outside K, L intersects K in two points and the line through p and each of them is tangent to K, cf. (Koecher and Krieg, 2000).
In the special case K = Cα, Eq. (18) simplifies to Comparing this with Eq. (10) shows that is the pole of L with respect to Cα, where L is the undistorted line corresponding to the distorted circle with center m. The sign in Eq. (20) is the sign of α.
In Fig. 5, let m be the center of the circle B. Since B intersects the outer circle orthogonally, the tangents at the intersections intersect in m and therefore, m is the pole of the line C. Since A is half way from C to the origin, the pole of A is 2m, twice that of C. This is in accordance with Eq. (20).

Cross-Ratio and Duality in Projective Geometry
The cross-ratio of four points a, b, c, and d on a line L is the number where o is an arbitrary point not on L and | · | denotes the determinant of a square matrix. The result of this formula is independent of the choice of o, cf. (Richter-Gebert, 2011, Corollary 4.1 and Lemma 4.6) and (Erdnüß, 2017, Sec. 5.3 and Eq. (57)). The cross-ratio is invariant under projective transformations, cf. (Richter-Gebert, 2011, Lemma 4.5). Perspective transformations include especially scaling, rotation, and translation.
In projective geometry, points and lines are dual to each other. Therefore, each theorem about points exists in a slightly modified version for lines as well, where mainly the line through two points is replaced by the intersection of two lines and vice versa, cf. also (Richter-Gebert, 2011, Sec. 3.5).
Therefore, there exists a projective invariant cross-ratio for four lines L k with k = 1, 2, 3, 4, that intersect in a common point η. By (Thomas, 2020, Prop. 4), the cross-ratio of the four lines agrees with the cross-ratio of their four poles p k with respect to any fixed circle K, i.e.
Let L k be lines in the scene plane and p k their polars with respect to Cα. Let C k = C(m k , m 2 k − 1 α ) be the circles in the sensor plane that are the distorted images of the lines L k . By Eq. (20), p k = ±2m k with the sign depending on the sign of α. Therefore, m k and p k agree up to scale and since the crossratio is invariant under scale, it follows that (L1, L2; L3, L4) = (m1, m2; m3, m4) .
In other words, the space of center points m of the circles (in the sensor plane) that correspond to distorted lines (of the scene plane), is projectively isomorphic to the space of the corresponding lines L (in the scene plane).
On the left hand side is a fact about a set of lines in the real world, that can be calculated by the term on the right hand side, which contains only the coordinates of the centers of circles as observed in the distorted sensor image.

Application to the Distance of Parallel Railroads
Given two parallel railroads with the four rails S k for k = 1, 2, 3, 4 in order and four points s k on them on a line perpendicular to the lines S k on the ground plane, cf. Fig. 7. Let g be  s1 s1 s1 s1 s1  s1  s1  s1  s1  s1  s1  s1  s1 s1 s1 s1 s1  s2 s2 s2 s2 s2  s2  s2  s2  s2  s2  s2  s2  s2 s2 s2 s2 s2  s3 s3 s3 s3 s3  s3  s3  s3  s3  s3  s3  s3  s3 s3 s3 s3 s3  s4 s4 s4 s4 s4  s4  s4  s4  s4  s4  s4  s4 s4 s4 s4 s4 s4 Figure 7. Railroads Distance. Let s k = (s k , 1) for k = 1, 2, 3, 4 such that in a local coordinate system with origin o = (0, 0) the distances s2 − s1 = s4 − s3 = g are the known distance between two rails on a railroad in world units (e. g. g = 1.5 m) and s3 − s1 = s4 − s2 = d are the unknown distance between the railroads. the known distance between the two rails of each railroad track and d the unknown distance between the two railroad tracks (left rail to left rail or right rail to right rail). The ratio d g can be measured (in some sense) from the left as the fraction s 3 −s 1 s 2 −s 1 but also from the right as s 4 −s 2 s 4 −s 2 . The ratio d g is therefore also the geometric mean of these two measurements However, under the square root is now the projectively invariant cross-ratio, cf. Eq. (21) q = (s1, s4; s3, s2) = s2 s4 0 1 1 0 1 1 1 Note that each of the determinants is just the difference between the upper two numbers.
By duality, this is also the cross-ratio q = (S1, S4; S3, S2) of the entire rails S k , cf. (Richter-Gebert, 2011, p. 77) and by Eq. (25), this is also the cross-ratio q = (m1, m4; m3, m2) of the centers m k of the circles C k that are the distorted images of the rails S k .
Therefore, the distance d between the two railroads is given by where m k can be directly measured in the uncalibrated distorted sensor image and g is known. Most railroads use the Stephenson gauge of 1.435 m, which specifies the inner distance between the two rails. Most rails have a width of around 0.07 m at the head, which makes a total distance g of almost exactly 1.5 m between the centers of the rails, that are visually best to detect.

Circles Through a Point
Let x be a point in the distorted image and C(m, r 2 ) a circle with r 2 = m 2 − 1 α through x, then m sits on the line with normal equation This is a corollary from Eq. (10) by substituting ξ with Eq. (3).
The direction in which the arc goes through x depends on the position where m is located on the line. The other way around, m can be determined by the direction of the arc passing through x.
Let w be the unit vector with w, x > 0 that is perpendicular to the arc through x, then where κ is the curvature of the arc. This is by the definition of curvature, cf. e.g. (Kühnel, 2010).
Substituting Eq. (30) into Eq. (29) and solving for κ leads to and back substituting this result in Eq. (30) gives Solving Eq. (31) for α −1 gives The equivalent of a line detector in the distorted image would be a circle arc detector. Most circle detectors look for small circles completely contained within the image. In this case, the arcs of interest have usually only a small curvature. Therefore, it is better to start with a line detector and aggregate line segments with a compatible curvature. w can be estimated from the image gradients and κ from the change of w in a local neighborhood of x. If α is unknown, each triplet x, w, κ votes for one α by Eq. (33).

Extension to Higher Order Radial Distortion
If analytical invertibility is required one can construct an invertible higher order undistortion model by iterating the oneparameter division undistortion model, i.e.
for k = 1, . . . , n with parameters α1, . . . , αn and inverse For a second order model with n = 2, this yields Neglecting the x 6 term in the denominator, the model is roughly equivalent to a classical second order radial distortion model with parameters k1 = α1 + α2 and k2 = −α1α2.

CONCLUSION
This article reviews the one-parameter division undistortion model suggested by (Lenz, 1987) and (Fitzgibbon, 2001). It cites some of its main sources, states some of its main properties, and shows a relations to hyperbolic geometry. The article suggests a new method to extract projective invariants about scene lines from uncalibrated distorted images by means of a circular arc detector and shows an application to measure the distance of two parallel railroads. It also points out, that a higher order radial distortion model can be build by iterating the one-parameter division undistortion model, without loosing analytic invertibility.