Potential issue with radial distortion computation

[oleg-alexandrov]

I think I found an issue with how radial distortion is computed, in the code: https://github.com/DOI-USGS/usgscsm/blob/main/src/Distortion.cpp#L91

The code is:

 ux = dx;
 uy = dy;

  switch (distortionType) {
    // Compute undistorted focal plane coordinate given a distorted
    // coordinate set and the distortion coefficients
    case RADIAL: {
      double rr = dx * dx + dy * dy;

      if (rr > tolerance) {
        double dr = opticalDistCoeffs[0] +
                    (rr * (opticalDistCoeffs[1] + rr * opticalDistCoeffs[2]));

        ux = dx * (1.0 - dr);
        uy = dy * (1.0 - dr);
      }
    } break;

The intention is to avoid dealing with small numbers, that is why the check has rr > tolerance (usually tolerance is 1e-6).

The problem is that the formula

double dr = opticalDistCoeffs[0] + (rr * (opticalDistCoeffs[1] + rr * opticalDistCoeffs[2]));

does not get applied when rr is very small, even though it has the term opticalDistCoeffs[0], which need not be small even if rr is small. This will create a discontinuity.

This kind of check for tolerance is not needed anyway, since nothing can go wrong when multiplying numbers close to 0. It is needed only when reversing this process, when one has to do a division.

That undistortion code, at https://github.com/DOI-USGS/usgscsm/blob/main/src/Distortion.cpp#L329, seems to suffer from the same problem. One has to first handle the opticalDistCoeffs[0] term (then have to divide by dx and dy) before deciding if it is worth dividing by rr.

In either case the tolerance of 1e-6 for rr is kind of high (https://github.com/DOI-USGS/usgscsm/blob/main/include/usgscsm/Distortion.h#L23). Since this is dx dx + dy dy, that means tolerance in dx and dy is 1e-3.

I am not sure if dx and dy are measured in meters or millimeters or pixels (lens quantities are usually in millimeters), but then 1e-3 times that is either a millimeter or a micron, and pixels on a CCD sensor can be on the order of 1-10 microns. So then this is about 1 pixel which may be too close for comfort.

Not sure if this will make a difference in practice.

[acpaquette]

@oleg-alexandrov Can this be closed as of merging #464?

[oleg-alexandrov]

Yes, indeed. Thanks.