Isis 3 Programmer Reference
TaylorCameraDistortionMap.cpp
Go to the documentation of this file.
1 
21 #include <cmath>
22 #include "IString.h"
24 
25 using namespace std;
26 namespace Isis {
27 
41  TaylorCameraDistortionMap::TaylorCameraDistortionMap(Camera *parent, double zDirection):
42  CameraDistortionMap(parent, zDirection) {
43  }
44 
75  void TaylorCameraDistortionMap::SetDistortion(const int naifIkCode) {
76  QString odtxkey = "INS" + toString(naifIkCode) + "_OD_T_X";
77  QString odtykey = "INS" + toString(naifIkCode) + "_OD_T_Y";
78  for(int i = 0; i < 10; ++i) {
79  p_odtx.push_back(p_camera->getDouble(odtxkey, i));
80  p_odty.push_back(p_camera->getDouble(odtykey, i));
81  }
82  }
83 
104  bool TaylorCameraDistortionMap::SetFocalPlane(const double dx, const double dy) {
105  p_focalPlaneX = dx;
106  p_focalPlaneY = dy;
107 
108  // No coefficients == no distortion
109  if(p_odtx.size() <= 0 && p_odty.size() <= 0) {
112  return true;
113  }
114 
115  // Solve the distortion equation using the Newton-Raphson method.
116  // Set the error tolerance to about one millionth of a NAC pixel.
117  const double tol = 1.4E-5;
118 
119  // The maximum number of iterations of the Newton-Raphson method.
120  const int maxTries = 20;
121 
122  double x;
123  double y;
124  double fx;
125  double fy;
126  double Jxx;
127  double Jxy;
128  double Jyx;
129  double Jyy;
130 
131  // Initial guess at the root
132  x = dx;
133  y = dy;
134 
135  this->DistortionFunction(x, y, &fx, &fy);
136 
137  for(int count = 1; ((fabs(fx) + fabs(fy)) > tol) && (count < maxTries); count++) {
138 
139  this->DistortionFunction(x, y, &fx, &fy);
140 
141  fx = dx - fx;
142  fy = dy - fy;
143 
144  this->DistortionFunctionJacobian(x, y, &Jxx, &Jxy, &Jyx, &Jyy);
145 
146  double determinant = Jxx * Jyy - Jxy * Jyx;
147  if(determinant < 1E-6) {
148  //
149  // Near-zero determinant. Add error handling here.
150  //
151  //-- Just break out and return with no convergence
152  break;
153  }
154 
155  x = x + (Jyy * fx - Jxy * fy) / determinant;
156  y = y + (Jxx * fy - Jyx * fx) / determinant;
157  }
158 
159  if((fabs(fx) + fabs(fy)) <= tol) {
160  // The method converged to a root.
163  }
164  else {
165  // The method did not converge to a root within the maximum
166  // number of iterations. Return with no distortion.
169  }
170 
171  return true;
172  }
173 
190  const double uy) {
193 
194  // No coefficients == nodistortion
195  if(p_odtx.size() <= 0 && p_odty.size() <= 0) {
196  p_focalPlaneX = ux;
197  p_focalPlaneY = uy;
198  return true;
199  }
200 
202 
203  return true;
204  }
205 
215  void TaylorCameraDistortionMap::DistortionFunction(double ux, double uy, double *dx, double *dy) {
216 
217  double f[10];
218  f[0] = 1;
219  f[1] = ux;
220  f[2] = uy;
221  f[3] = ux * ux;
222  f[4] = ux * uy;
223  f[5] = uy * uy;
224  f[6] = ux * ux * ux;
225  f[7] = ux * ux * uy;
226  f[8] = ux * uy * uy;
227  f[9] = uy * uy * uy;
228 
229  *dx = 0.0;
230  *dy = 0.0;
231 
232  for(int i = 0; i < 10; i++) {
233  *dx = *dx + f[i] * p_odtx[i];
234  *dy = *dy + f[i] * p_odty[i];
235  }
236 
237  }
238 
251  void TaylorCameraDistortionMap::DistortionFunctionJacobian(double x, double y, double *Jxx, double *Jxy, double *Jyx, double *Jyy) {
252 
253  double d_dx[10];
254  d_dx[0] = 0;
255  d_dx[1] = 1;
256  d_dx[2] = 0;
257  d_dx[3] = 2 * x;
258  d_dx[4] = y;
259  d_dx[5] = 0;
260  d_dx[6] = 3 * x * x;
261  d_dx[7] = 2 * x * y;
262  d_dx[8] = y * y;
263  d_dx[9] = 0;
264  double d_dy[10];
265  d_dy[0] = 0;
266  d_dy[1] = 0;
267  d_dy[2] = 1;
268  d_dy[3] = 0;
269  d_dy[4] = x;
270  d_dy[5] = 2 * y;
271  d_dy[6] = 0;
272  d_dy[7] = x * x;
273  d_dy[8] = 2 * x * y;
274  d_dy[9] = 3 * y * y;
275 
276  *Jxx = 0.0;
277  *Jxy = 0.0;
278  *Jyx = 0.0;
279  *Jyy = 0.0;
280 
281  for(int i = 0; i < 10; i++) {
282  *Jxx = *Jxx + d_dx[i] * p_odtx[i];
283  *Jxy = *Jxy + d_dy[i] * p_odtx[i];
284  *Jyx = *Jyx + d_dx[i] * p_odty[i];
285  *Jyy = *Jyy + d_dy[i] * p_odty[i];
286  }
287  }
288 
289 }
double p_focalPlaneX
Distorted focal plane x.
bool SetFocalPlane(const double dx, const double dy)
Compute undistorted focal plane x/y.
std::vector< double > p_odty
distortion y coefficients
Namespace for the standard library.
bool SetUndistortedFocalPlane(const double ux, const double uy)
Compute distorted focal plane x/y.
QString toString(bool boolToConvert)
Global function to convert a boolean to a string.
Definition: IString.cpp:226
void DistortionFunctionJacobian(double x, double y, double *Jxx, double *Jxy, double *Jyx, double *Jyy)
Jacobian of the distortion function.
double p_undistortedFocalPlaneX
Undistorted focal plane x.
Camera * p_camera
The camera to distort/undistort.
std::vector< double > p_odtx
distortion x coefficients
Unless noted otherwise, the portions of Isis written by the USGS are public domain.
Distort/undistort focal plane coordinates.
const double E
Sets some basic constants for use in ISIS programming.
Definition: Constants.h:55
double p_focalPlaneY
Distorted focal plane y.
double p_undistortedFocalPlaneY
Undistorted focal plane y.
Namespace for ISIS/Bullet specific routines.
Definition: Apollo.h:31
void SetDistortion(const int naifIkCode)
Load distortion coefficients.
void DistortionFunction(double ux, double uy, double *dx, double *dy)
Compute distorted focal plane dx,dy given an undistorted focal plane ux,uy.
SpiceDouble getDouble(const QString &key, int index=0)
This returns a value from the NAIF text pool.
Definition: Spice.cpp:963