.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "../../docs/source/_gallery/optimize_parameters.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_.._.._docs_source__gallery_optimize_parameters.py: .. _sphx_glr__gallery_optimize_parameters.py: Optimizing distortion parameters with least squares ====================================================================================================== This example illustrate how to use the ``optimize_parameters_trf`` function to optimize the distortion parameters of a camera model, which includes the intrinsic and distortion transformations. .. seealso:: - :func:`pycvcam.optimize_parameters_trf`: Optimize the parameters of a transformation using Trust Region Reflective optimization. - :func:`pycvcam.optimize_camera_trf`: Optimize the parameters of a complete camera model using Trust Region Reflective optimization. .. GENERATED FROM PYTHON SOURCE LINES 18-31 Optimizing Parameters to fit a Distortion model (Image Alignment) ------------------------------------------------------------------- In this example, we consider a random image, an intrinsic matrix and a distortion model. We distort the image by the distortion model and then add noise. The objective is to optimize the distortion parameters of the distortion model to fit the distorted image with noise. To do this, we will compute the flow between the both images and then optimize the distortion parameters to minimize the flow in the normalized image space. First create the images and compute the real flow between the original and distorted images. .. GENERATED FROM PYTHON SOURCE LINES 31-110 .. code-block:: Python import numpy import pycvcam import cv2 import matplotlib.pyplot as plt image = pycvcam.get_lena_image() image_height, image_width = image.shape[:2] print("Image shape:", image.shape) pixel_points = numpy.indices((image_height, image_width), dtype=numpy.float64) pixel_points = pixel_points.reshape(2, -1).T # shape (H*W, 2) WARNING: [H, W -> Y, X] image_points = pixel_points[:, [1, 0]] # Swap to [X, Y] format # Create an intrisic transformation intrinsic = pycvcam.Cv2Intrinsic.from_matrix( [[1000.0, 0.0, image_width / 2], [0.0, 1000.0, image_height / 2], [0.0, 0.0, 1.0]] ) # Create a distortion transformation in the image space and convert it to the normalized image space distortion = pycvcam.ZernikeDistortion( parameters=[ 0.8541972545746392, -5.468596289790535, -5.974287819021697, 14.292956075116104, 2.1403205479372627, 4.544169430137205, -0.10099732464199339, 0.4363509204067417, -0.5106374355681896, -5.770087687650705, -0.39147505788710696, 11.699411273002498, ] # In pixels units, example values for a small distortion ) distortion.center = ((image_width - 1) / 2, (image_height - 1) / 2) distortion.radius = numpy.sqrt( (distortion.center[0]) ** 2 + (distortion.center[1]) ** 2 ) distortion.parameters_x = distortion.parameters_x / intrinsic.fx distortion.parameters_y = distortion.parameters_y / intrinsic.fy distortion.radius_x = distortion.radius_x / intrinsic.fx distortion.radius_y = distortion.radius_y / intrinsic.fy distortion.center_x = (distortion.center_x - intrinsic.cx) / intrinsic.fx distortion.center_y = (distortion.center_y - intrinsic.cy) / intrinsic.fy true_distorted_points = pycvcam.distort_points( image_points, intrinsic, distortion, P=intrinsic ) true_flow = true_distorted_points - image_points # shape (H*W, 2) true_flow_magnitude = numpy.linalg.norm(true_flow, axis=1) print("True flow magnitude statistics:") print(" Min:", numpy.min(true_flow_magnitude)) print(" Max:", numpy.max(true_flow_magnitude)) print(" Mean:", numpy.mean(true_flow_magnitude)) print(" Median:", numpy.median(true_flow_magnitude)) print(" Std:", numpy.std(true_flow_magnitude)) # Create the distorted image with noise (5% of noise in the GL units) distorted_image = pycvcam.distort_image( image, intrinsic, distortion, interpolation="spline3" ) noise_GLpc = 5.0 noise = numpy.random.normal(0, noise_GLpc / 100, distorted_image.shape) distorted_image = distorted_image * (1 + noise) # Add multiplicative noise distorted_image = numpy.clip(distorted_image, 0, 255).astype(numpy.uint8) fig = plt.figure(figsize=(10, 5)) ax1 = fig.add_subplot(1, 2, 1) ax1.imshow(image, cmap="gray") ax1.set_title("Original Image") ax1.axis("off") ax2 = fig.add_subplot(1, 2, 2) ax2.imshow(distorted_image, cmap="gray") ax2.set_title("Distorted Image with Noise") ax2.axis("off") .. image-sg:: /_gallery/images/sphx_glr_optimize_parameters_001.png :alt: Original Image, Distorted Image with Noise :srcset: /_gallery/images/sphx_glr_optimize_parameters_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Image shape: (474, 474) True flow magnitude statistics: Min: 0.012265189243667858 Max: 24.30234213301494 Mean: 6.324732682604791 Median: 3.94728868310939 Std: 5.5042902818601105 (np.float64(-0.5), np.float64(473.5), np.float64(473.5), np.float64(-0.5)) .. GENERATED FROM PYTHON SOURCE LINES 111-113 Then compute the optical flow between both images using the ``compute_optical_flow`` function, which computes the dense optical flow between two images using the DIS algorithm. .. GENERATED FROM PYTHON SOURCE LINES 113-168 .. code-block:: Python flow_x, flow_y = pycvcam.compute_optical_flow( image, distorted_image, disflow_params={ "PatchSize": 15, "PatchStride": 1, }, ) flow_x = flow_x.reshape(-1, 1) # shape (H*W, 1) flow_y = flow_y.reshape(-1, 1) # shape (H*W, 1) flow = numpy.hstack((flow_x, flow_y)) # shape (H*W, 2) # Compare the computed flow with the true flow by visualizing the flow magnitude and # components for both the computed flow and the true flow to ensure consistency. fx = flow[:, 0].reshape(image_height, image_width) fy = flow[:, 1].reshape(image_height, image_width) F = numpy.sqrt(fx**2 + fy**2) tfx = true_flow[:, 0].reshape(image_height, image_width) tfy = true_flow[:, 1].reshape(image_height, image_width) tF = numpy.sqrt(tfx**2 + tfy**2) tF_max = numpy.max(tF) tfx_min, tfx_max = numpy.min(tfx), numpy.max(tfx) tfy_min, tfy_max = numpy.min(tfy), numpy.max(tfy) fig = plt.figure(figsize=(15, 10)) ax1 = fig.add_subplot(2, 3, 1) ax1.imshow(F, cmap="inferno", vmin=0, vmax=tF_max) ax1.set_title("Flow Magnitude") ax1.axis("off") ax2 = fig.add_subplot(2, 3, 2) ax2.imshow(fx, cmap="inferno", vmin=tfx_min, vmax=tfx_max) ax2.set_title("Flow X Component") ax2.axis("off") ax3 = fig.add_subplot(2, 3, 3) ax3.imshow(fy, cmap="inferno", vmin=tfy_min, vmax=tfy_max) ax3.set_title("Flow Y Component") ax3.axis("off") ax4 = fig.add_subplot(2, 3, 4) ax4.imshow(tF, cmap="inferno", vmin=0, vmax=tF_max) ax4.set_title("True Flow Magnitude") ax4.axis("off") ax5 = fig.add_subplot(2, 3, 5) ax5.imshow(tfx, cmap="inferno", vmin=tfx_min, vmax=tfx_max) ax5.set_title("True Flow X Component") ax5.axis("off") ax6 = fig.add_subplot(2, 3, 6) ax6.imshow(tfy, cmap="inferno", vmin=tfy_min, vmax=tfy_max) ax6.set_title("True Flow Y Component") ax6.axis("off") # plt.show() .. image-sg:: /_gallery/images/sphx_glr_optimize_parameters_002.png :alt: Flow Magnitude, Flow X Component, Flow Y Component, True Flow Magnitude, True Flow X Component, True Flow Y Component :srcset: /_gallery/images/sphx_glr_optimize_parameters_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none (np.float64(-0.5), np.float64(473.5), np.float64(473.5), np.float64(-0.5)) .. GENERATED FROM PYTHON SOURCE LINES 169-174 To finish, we optimize the distortion parameters to minimize the flow in the normalized image space. To avoid the border effects of the distortion and undistortion, we will consider only the inner part of the image for the optimization by applying a mask to the points. .. GENERATED FROM PYTHON SOURCE LINES 174-221 .. code-block:: Python images_mask = numpy.ones_like(image, dtype=bool) border_size = 50 images_mask[0:border_size, :] = False images_mask[-border_size:, :] = False images_mask[:, 0:border_size] = False images_mask[:, -border_size:] = False distorted_points = image_points + flow # shape (H*W, 2) normalized_points = intrinsic.inverse_transform(image_points).transformed_points distorted_points = intrinsic.inverse_transform(distorted_points).transformed_points normalized_points = normalized_points[images_mask.flatten()] distorted_points = distorted_points[images_mask.flatten()] initial_distortion = distortion.copy() initial_distortion.parameters = numpy.zeros_like(distortion.parameters) print("\n") parameters, result = pycvcam.optimize_parameters_lm( initial_distortion, normalized_points, distorted_points, auto=True, # Set ftol, xtol and gtol to 1e-8 return_result=True, verbose_level=3, ) print("\n") optimize_distortion = initial_distortion.copy() optimize_distortion.parameters = parameters optimized_flow = ( pycvcam.distort_points(image_points, intrinsic, optimize_distortion, P=intrinsic) - image_points ) flow_error = numpy.linalg.norm(optimized_flow - true_flow, axis=1) # shape (H*W,) rmse_flow = numpy.sqrt(numpy.mean(flow_error[images_mask.flatten()] ** 2)) params_error = numpy.linalg.norm(parameters - distortion.parameters) params_rel_error = params_error / numpy.linalg.norm(distortion.parameters) print("Optimization success:", result.success) print("Optimization message:", result.message) print("Optimization cost:", result.cost) print("Flow RMSE:", rmse_flow) print("Parameters Error:", params_error, f"({params_rel_error:.2%})") .. rst-class:: sphx-glr-script-out .. code-block:: none ================================================== -------------------------------------------------- Initial Jacobian (Iter 0) analysis of the least squares problem -------------------------------------------------- Jacobian shape: 279752 x 12 (equations x parameters) Overdetermined system (more residuals than parameters), the optimization is likely to converge to a unique solution. Density: 50.00% Singular values (max/min): 4.359e+02 / 4.929e+01 Condition number: 8.843e+00 Singular values and their contribution to the variance: | Index | Singular Value λ | Var = 1/λ^2 | | 0 | 4.359e+02 | 5.263e-06 | | 1 | 4.359e+02 | 5.263e-06 | | 2 | 1.207e+02 | 6.861e-05 | | 3 | 1.207e+02 | 6.861e-05 | | 4 | 1.207e+02 | 6.861e-05 | | 5 | 1.207e+02 | 6.861e-05 | | 6 | 8.459e+01 | 1.398e-04 | | 7 | 8.459e+01 | 1.398e-04 | | 8 | 7.794e+01 | 1.646e-04 | | 9 | 7.794e+01 | 1.646e-04 | | 10 | 4.929e+01 | 4.115e-04 | | 11 | 4.929e+01 | 4.115e-04 | Higher and Smallers singular values directions vectors (V^T rows): | Parameter | Vt[0] | Vt[9] | Vt[10] | Vt[11] | | 0 | -0.000e+00 | 0.000e+00 | -0.000e+00 | 0.000e+00 | | 1 | 8.519e-01 | -9.145e-17 | 5.737e-16 | 3.896e-23 | | 2 | -3.123e-19 | 3.271e-17 | 1.068e-18 | -1.110e-16 | | 3 | 4.984e-17 | 4.192e-16 | -1.350e-16 | 0.000e+00 | | 4 | 2.094e-19 | 4.552e-17 | 1.487e-18 | -5.551e-17 | | 5 | -2.220e-16 | 1.223e-14 | 3.287e-17 | -2.303e-28 | | 6 | 1.943e-18 | 2.273e-31 | -2.246e-30 | 2.359e-16 | | 7 | -7.345e-19 | -1.000e+00 | -7.190e-16 | -1.543e-33 | | 8 | 0.000e+00 | 0.000e+00 | 0.000e+00 | -3.864e-17 | | 9 | -5.236e-01 | -1.474e-16 | 9.699e-16 | -2.395e-23 | | 10 | -5.740e-20 | -4.976e-32 | -2.136e-32 | -1.000e+00 | | 11 | 1.909e-17 | -8.749e-16 | 1.000e+00 | -3.749e-32 | Estimated variances of the parameters: | Parameter | Value P | Var = σ^2 (J.T J)^-1 | Ratio √V/|P| | | 0 | 0.000e+00 | 8.028e-10 | > 1000 % | | 1 | 0.000e+00 | 8.028e-10 | > 1000 % | | 2 | 0.000e+00 | 1.307e-09 | > 1000 % | | 3 | 0.000e+00 | 1.307e-09 | > 1000 % | | 4 | 0.000e+00 | 1.307e-09 | > 1000 % | | 5 | 0.000e+00 | 1.307e-09 | > 1000 % | | 6 | 0.000e+00 | 3.136e-09 | > 1000 % | | 7 | 0.000e+00 | 3.136e-09 | > 1000 % | | 8 | 0.000e+00 | 1.960e-09 | > 1000 % | | 9 | 0.000e+00 | 1.960e-09 | > 1000 % | | 10 | 0.000e+00 | 7.840e-09 | > 1000 % | | 11 | 0.000e+00 | 7.840e-09 | > 1000 % | -------------------------------------------------- Optimization in progress... -------------------------------------------------- `gtol` termination condition is satisfied. Function evaluations 2, initial cost 2.6644e+00, final cost 1.8526e-03, first-order optimality 9.31e-12. -------------------------------------------------- Jacobian analysis of the least squares problem (End of optimization) -------------------------------------------------- Singular values (max/min): 4.359e+02 / 4.929e+01 Condition number: 8.843e+00 Singular values and their contribution to the variance: | Index | Singular Value λ | Var = 1/λ^2 | | 0 | 4.359e+02 | 5.263e-06 | | 1 | 4.359e+02 | 5.263e-06 | | 2 | 1.207e+02 | 6.861e-05 | | 3 | 1.207e+02 | 6.861e-05 | | 4 | 1.207e+02 | 6.861e-05 | | 5 | 1.207e+02 | 6.861e-05 | | 6 | 8.459e+01 | 1.398e-04 | | 7 | 8.459e+01 | 1.398e-04 | | 8 | 7.794e+01 | 1.646e-04 | | 9 | 7.794e+01 | 1.646e-04 | | 10 | 4.929e+01 | 4.115e-04 | | 11 | 4.929e+01 | 4.115e-04 | Higher and Smallers singular values directions vectors (V^T rows): | Parameter | Vt[0] | Vt[9] | Vt[10] | Vt[11] | | 0 | -0.000e+00 | 0.000e+00 | -0.000e+00 | 0.000e+00 | | 1 | 8.519e-01 | -9.145e-17 | 5.737e-16 | 3.896e-23 | | 2 | -3.123e-19 | 3.271e-17 | 1.068e-18 | -1.110e-16 | | 3 | 4.984e-17 | 4.192e-16 | -1.350e-16 | 0.000e+00 | | 4 | 2.094e-19 | 4.552e-17 | 1.487e-18 | -5.551e-17 | | 5 | -2.220e-16 | 1.223e-14 | 3.287e-17 | -2.303e-28 | | 6 | 1.943e-18 | 2.273e-31 | -2.246e-30 | 2.359e-16 | | 7 | -7.345e-19 | -1.000e+00 | -7.190e-16 | -1.543e-33 | | 8 | 0.000e+00 | 0.000e+00 | 0.000e+00 | -3.864e-17 | | 9 | -5.236e-01 | -1.474e-16 | 9.699e-16 | -2.395e-23 | | 10 | -5.740e-20 | -4.976e-32 | -2.136e-32 | -1.000e+00 | | 11 | 1.909e-17 | -8.749e-16 | 1.000e+00 | -3.749e-32 | Estimated variances of the parameters: | Parameter | Value P | Var = σ^2 (J.T J)^-1 | Ratio √V/|P| | | 0 | 8.530e-04 | 5.582e-13 | 0.088 % | | 1 | -5.470e-03 | 5.582e-13 | 0.014 % | | 2 | -5.963e-03 | 9.088e-13 | 0.016 % | | 3 | 1.428e-02 | 9.088e-13 | 0.007 % | | 4 | 2.142e-03 | 9.088e-13 | 0.045 % | | 5 | 4.524e-03 | 9.088e-13 | 0.021 % | | 6 | -1.123e-04 | 2.180e-12 | 1.315 % | | 7 | 4.616e-04 | 2.180e-12 | 0.320 % | | 8 | -5.153e-04 | 1.363e-12 | 0.227 % | | 9 | -5.768e-03 | 1.363e-12 | 0.020 % | | 10 | -4.051e-04 | 5.451e-12 | 0.576 % | | 11 | 1.173e-02 | 5.451e-12 | 0.020 % | ================================================== Optimization success: True Optimization message: `gtol` termination condition is satisfied. Optimization cost: 0.0018525991108891717 Flow RMSE: 0.011271562303385487 Parameters Error: 4.856649693275793e-05 (0.22%) .. GENERATED FROM PYTHON SOURCE LINES 222-231 Optimizing Parameters of a complete Camera Model (PnP problem) ----------------------------------------------------------------- In this example, we optimize the parameters of a complete camera model, which includes the intrinsic, distortion and extrinsic transformations. To do this, consider a PnP problem with a set of 3D points and their corresponding 2D projections in the image. For example, create a point cloud in a rectangle with bounds :math:`[-1, 1]` meters in the :math:`x` and :math:`y` directions and with a :math:`z` component in the range :math:`[4.5, 5.5]` meters. Then place a camera near the origin with a small rotation and translation and project the 3D points to 2D image points using a pinhole camera model with Brown-Conrady distortion with a focal length of 1000 pixels and a principal point at (320, 240) pixels. .. GENERATED FROM PYTHON SOURCE LINES 231-330 .. code-block:: Python print("\n\n\n") # Define the 3D points in the world coordinate system x = numpy.random.uniform(-1.0, 1.0, (100, 1)) # shape (100, 1) y = numpy.random.uniform(-1.0, 1.0, (100, 1)) # shape (100, 1) z = numpy.random.uniform(4.5, 5.5, (100, 1)) # shape (100, 1) world_points = numpy.hstack((x, y, z)) # shape (100, 3) # Define the Extrinsic transformation, example rotation vector and translation vector rvec = numpy.array([0.01, 0.02, 0.03]) # small rotation tvec = numpy.array([0.01, -0.05, 0.04]) # small translation extrinsic = pycvcam.Cv2Extrinsic.from_rt(rvec, tvec) # Define the Intrinsic transformation, example intrinsic camera matrix K = numpy.array([[1000.0, 0.0, 320.0], [0.0, 1000.0, 240.0], [0.0, 0.0, 1.0]]) intrinsic = pycvcam.Cv2Intrinsic.from_matrix(K) # Define the Distortion transformation, example Brown-Conrady 5 parameters distortion = pycvcam.Cv2Distortion(parameters=[0.1, 0.05, 0.01, 0.01, 0.01]) # Project the 3D points to 2D image points (without Jacobians) result = pycvcam.project_points( world_points, intrinsic=intrinsic, distortion=distortion, extrinsic=extrinsic, ) # pycvcam.TransformResult data class image_points = result.image_points # shape (100, 2) # Optimize the camera parameters to fit the noisy image points initial_distortion = distortion.copy() initial_distortion.parameters = numpy.zeros_like(distortion.parameters) initial_extrinsic = extrinsic.copy() initial_extrinsic.parameters = numpy.zeros_like(extrinsic.parameters) distortion_bounds = ( numpy.array([-0.1, -0.1, -0.1, -0.1, -0.1]), # k1, k2, p1, p2, k3 lower bounds numpy.array([0.1, 0.1, 0.1, 0.1, 0.1]), # k1, k2, p1, p2, k3 upper bounds ) extrinsic_bounds = ( # rvec and tvec bounds numpy.array([-0.1, -0.1, -0.1, -0.1, -0.1, -0.1]), # rvec and tvec lower bounds numpy.array([0.1, 0.1, 0.1, 0.1, 0.1, 0.1]), # rvec and tvec upper bounds ) print("\n") params, result = pycvcam.optimize_camera_trf( intrinsic, initial_distortion, initial_extrinsic, world_points, image_points, mask_intrinsic=[False for _ in range(4)], # Do not optimize intrinsic bounds_distortion=distortion_bounds, bounds_extrinsic=extrinsic_bounds, auto=True, # Set ftol, xtol and gtol to 1e-8 return_result=True, verbose_level=3, ) print("\n") optimized_intrinsic_params, optimized_distortion_params, optimized_extrinsic_params = ( params ) print("Optimized Intrinsic Parameters:", optimized_intrinsic_params) print("Optimized Distortion Parameters:", optimized_distortion_params) print("Optimized Extrinsic Parameters:", optimized_extrinsic_params) optimized_distortion = initial_distortion.copy() optimized_distortion.parameters = optimized_distortion_params optimized_extrinsic = initial_extrinsic.copy() optimized_extrinsic.parameters = optimized_extrinsic_params optimized_image_points = pycvcam.project_points( world_points, intrinsic=intrinsic, distortion=optimized_distortion, extrinsic=optimized_extrinsic, ).image_points error_points = numpy.linalg.norm(optimized_image_points - image_points, axis=1) rmse_points = numpy.sqrt(numpy.mean(error_points**2)) error_distortion = numpy.linalg.norm( optimized_distortion.parameters - distortion.parameters ) error_extrinsic = numpy.linalg.norm( optimized_extrinsic.parameters - extrinsic.parameters ) rel_error_distortion = error_distortion / numpy.linalg.norm(distortion.parameters) rel_error_extrinsic = error_extrinsic / numpy.linalg.norm(extrinsic.parameters) print("Optimization success:", result.success) print("Optimization message:", result.message) print("Optimization cost:", result.cost) print("RMSE Real Points:", rmse_points) print("Error Distortion:", error_distortion, f"({rel_error_distortion:.2%})") print("Error Extrinsic:", error_extrinsic, f"({rel_error_extrinsic:.2%})") .. rst-class:: sphx-glr-script-out .. code-block:: none ================================================== -------------------------------------------------- Initial Jacobian (Iter 0) analysis of the least squares problem -------------------------------------------------- 5 Distortion parameters to optimize - Parameters 0 to 4 6 Extrinsic parameters to optimize - Parameters 5 to 10 Jacobian shape: 200 x 11 (equations x parameters) Overdetermined system (more residuals than parameters), the optimization is likely to converge to a unique solution. Density: 90.91% Singular values (max/min): 1.036e+04 / 9.951e-03 Condition number: 1.041e+06 Singular values and their contribution to the variance: | Index | Singular Value λ | Var = 1/λ^2 | | 0 | 1.036e+04 | 9.316e-09 | | 1 | 1.035e+04 | 9.337e-09 | | 2 | 1.658e+03 | 3.637e-07 | | 3 | 4.893e+02 | 4.177e-06 | | 4 | 4.596e+02 | 4.735e-06 | | 5 | 3.294e+02 | 9.217e-06 | | 6 | 1.085e+02 | 8.501e-05 | | 7 | 1.060e+02 | 8.906e-05 | | 8 | 2.499e+01 | 1.601e-03 | | 9 | 5.173e-01 | 3.737e+00 | | 10 | 9.951e-03 | 1.010e+04 | Higher and Smallers singular values directions vectors (V^T rows): | Parameter | Vt[0] | Vt[8] | Vt[9] | Vt[10] | | 0 | -7.378e-04 | -9.782e-01 | 7.658e-02 | 3.763e-03 | | 1 | -3.284e-05 | -7.618e-02 | -9.905e-01 | -1.143e-01 | | 2 | 5.576e-02 | -3.776e-03 | 1.579e-04 | -3.909e-06 | | 3 | -3.889e-04 | -6.711e-03 | -3.757e-04 | 8.491e-06 | | 4 | -1.473e-06 | -5.026e-03 | -1.142e-01 | 9.934e-01 | | 5 | -9.793e-01 | -4.963e-03 | 3.184e-05 | -3.023e-06 | | 6 | -5.323e-03 | 1.598e-03 | 9.892e-05 | -2.128e-06 | | 7 | 9.237e-03 | 2.033e-04 | -2.560e-06 | 9.475e-08 | | 8 | -1.057e-03 | -6.410e-03 | -4.250e-04 | 9.377e-06 | | 9 | 1.941e-01 | -2.453e-02 | 1.398e-04 | -1.429e-05 | | 10 | 3.227e-03 | -1.912e-01 | 5.889e-03 | 1.608e-04 | Estimated variances of the parameters: | Parameter | Value P | Var = σ^2 (J.T J)^-1 | Ratio √V/|P| | | 0 | 0.000e+00 | 8.204e+01 | > 1000 % | | 1 | 0.000e+00 | 6.678e+04 | > 1000 % | | 2 | 0.000e+00 | 2.578e-03 | > 1000 % | | 3 | 0.000e+00 | 2.960e-03 | > 1000 % | | 4 | 0.000e+00 | 4.912e+06 | > 1000 % | | 5 | 0.000e+00 | 1.751e-03 | > 1000 % | | 6 | 0.000e+00 | 1.689e-03 | > 1000 % | | 7 | 0.000e+00 | 1.838e-04 | > 1000 % | | 8 | 0.000e+00 | 4.117e-02 | > 1000 % | | 9 | 0.000e+00 | 4.329e-02 | > 1000 % | | 10 | 0.000e+00 | 2.257e-01 | > 1000 % | -------------------------------------------------- Optimization in progress... -------------------------------------------------- Iteration Total nfev Cost Cost reduction Step norm Optimality 0 1 4.6567e+04 2.37e+05 1 2 1.0851e+03 4.55e+04 4.16e-02 2.89e+04 2 3 2.0405e+01 1.06e+03 2.94e-02 9.54e+02 3 4 4.1792e+00 1.62e+01 1.23e-01 3.21e+01 4 5 2.2662e-01 3.95e+00 2.80e-02 3.51e+00 5 6 1.1181e-02 2.15e-01 1.09e-02 2.88e-01 6 7 1.6713e-03 9.51e-03 4.96e-03 8.22e-03 7 8 4.7294e-04 1.20e-03 1.68e-03 1.45e-03 8 9 1.7819e-04 2.95e-04 1.92e-01 3.36e-01 9 10 4.7250e-05 1.31e-04 1.54e-02 3.26e-04 10 11 6.4814e-06 4.08e-05 1.20e-02 2.06e-04 11 12 1.2486e-06 5.23e-06 4.24e-03 2.62e-05 12 13 5.5937e-07 6.89e-07 1.73e-03 4.67e-06 13 14 1.2333e-07 4.36e-07 5.93e-02 2.86e-06 14 15 2.7711e-08 9.56e-08 2.73e-02 6.02e-07 15 16 6.5798e-09 2.11e-08 1.24e-02 1.24e-07 16 17 1.6247e-09 4.96e-09 5.87e-03 2.76e-08 17 18 4.0558e-10 1.22e-09 2.89e-03 6.93e-09 `gtol` termination condition is satisfied. Function evaluations 18, initial cost 4.6567e+04, final cost 4.0558e-10, first-order optimality 6.93e-09. -------------------------------------------------- Jacobian analysis of the least squares problem (End of optimization) -------------------------------------------------- Singular values (max/min): 1.036e+04 / 9.186e-03 Condition number: 1.128e+06 Singular values and their contribution to the variance: | Index | Singular Value λ | Var = 1/λ^2 | | 0 | 1.036e+04 | 9.316e-09 | | 1 | 1.033e+04 | 9.364e-09 | | 2 | 1.653e+03 | 3.658e-07 | | 3 | 5.725e+02 | 3.052e-06 | | 4 | 4.852e+02 | 4.247e-06 | | 5 | 2.725e+02 | 1.347e-05 | | 6 | 1.086e+02 | 8.484e-05 | | 7 | 1.058e+02 | 8.933e-05 | | 8 | 2.433e+01 | 1.689e-03 | | 9 | 5.171e-01 | 3.740e+00 | | 10 | 9.186e-03 | 1.185e+04 | Higher and Smallers singular values directions vectors (V^T rows): | Parameter | Vt[0] | Vt[8] | Vt[9] | Vt[10] | | 0 | -2.385e-03 | -9.796e-01 | 7.712e-02 | 3.599e-03 | | 1 | -1.285e-04 | -7.676e-02 | -9.908e-01 | -1.115e-01 | | 2 | 2.824e-02 | -2.149e-02 | -1.686e-04 | 6.976e-07 | | 3 | -5.226e-02 | 1.390e-02 | 1.302e-04 | 2.720e-06 | | 4 | -7.162e-06 | -5.039e-03 | -1.115e-01 | 9.938e-01 | | 5 | -4.446e-01 | -4.617e-03 | -9.360e-05 | -1.482e-06 | | 6 | -8.723e-01 | 1.087e-03 | -9.606e-05 | 1.643e-07 | | 7 | -9.951e-03 | 2.512e-04 | 1.986e-06 | -7.632e-08 | | 8 | -1.742e-01 | -7.572e-03 | 4.633e-04 | -9.670e-07 | | 9 | 8.538e-02 | -1.971e-02 | -3.999e-04 | -7.164e-06 | | 10 | 9.588e-03 | -1.828e-01 | 5.933e-03 | 1.486e-04 | Estimated variances of the parameters: | Parameter | Value P | Var = σ^2 (J.T J)^-1 | Ratio √V/|P| | | 0 | 9.999e-02 | 7.612e-13 | 0.001 % | | 1 | 5.034e-02 | 6.482e-10 | 0.051 % | | 2 | 1.000e-02 | 2.445e-17 | 0.000 % | | 3 | 1.000e-02 | 2.171e-17 | 0.000 % | | 4 | 7.136e-03 | 5.023e-08 | 3.141 % | | 5 | 1.000e-02 | 1.513e-17 | 0.000 % | | 6 | 2.000e-02 | 1.441e-17 | 0.000 % | | 7 | 3.000e-02 | 1.622e-18 | 0.000 % | | 8 | 1.000e-02 | 3.543e-16 | 0.000 % | | 9 | -5.000e-02 | 3.730e-16 | 0.000 % | | 10 | 4.000e-02 | 1.980e-15 | 0.000 % | ================================================== Optimized Intrinsic Parameters: [1000. 1000. 320. 240.] Optimized Distortion Parameters: [0.09998801 0.05034095 0.01 0.00999999 0.00713579] Optimized Extrinsic Parameters: [ 0.01000001 0.02 0.03 0.00999999 -0.04999997 0.03999944] Optimization success: True Optimization message: `gtol` termination condition is satisfied. Optimization cost: 4.055768652509038e-10 RMSE Real Points: 2.8480760707920136e-06 Error Distortion: 0.002884455814632005 (2.55%) Error Extrinsic: 5.652400610347772e-07 (0.00%) .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 2.322 seconds) .. _sphx_glr_download_.._.._docs_source__gallery_optimize_parameters.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: optimize_parameters.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: optimize_parameters.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: optimize_parameters.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_