.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/example_curvature.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_auto_examples_example_curvature.py: .. _example_curvature: =================================== example of curvature estimation in slam =================================== .. GENERATED FROM PYTHON SOURCE LINES 8-15 .. code-block:: Python # Authors: Guillaume Auzias # Julien Barrès # License: MIT # sphinx_gallery_thumbnail_number = 2 .. GENERATED FROM PYTHON SOURCE LINES 16-20 NOTE: there is no visualization tool in slam, but we provide at the end of this script exemplare code to do the visualization with an external solution ############################################################################## .. GENERATED FROM PYTHON SOURCE LINES 20-28 .. code-block:: Python # importation of slam modules import slam.utils as ut import numpy as np import slam.generate_parametric_surfaces as sgps import slam.io as sio import slam.curvature as scurv .. GENERATED FROM PYTHON SOURCE LINES 29-30 loading an examplar mesh .. GENERATED FROM PYTHON SOURCE LINES 30-33 .. code-block:: Python mesh_file = "../examples/data/example_mesh.gii" mesh = sio.load_mesh(mesh_file) .. GENERATED FROM PYTHON SOURCE LINES 34-35 Comptue estimations of principal curvatures .. GENERATED FROM PYTHON SOURCE LINES 35-38 .. code-block:: Python PrincipalCurvatures, PrincipalDir1, PrincipalDir2 \ = scurv.curvatures_and_derivatives(mesh) .. rst-class:: sphx-glr-script-out .. code-block:: none Calculating vertex normals .... Please wait Finished calculating vertex normals Calculating curvature tensors ... Please wait Finished Calculating curvature tensors Calculating Principal Components ... Please wait Finished Calculating principal components .. GENERATED FROM PYTHON SOURCE LINES 39-40 Comptue Gauss curvature from principal curvatures .. GENERATED FROM PYTHON SOURCE LINES 40-42 .. code-block:: Python gaussian_curv = PrincipalCurvatures[0, :] * PrincipalCurvatures[1, :] .. GENERATED FROM PYTHON SOURCE LINES 43-44 Comptue mean curvature from principal curvatures .. GENERATED FROM PYTHON SOURCE LINES 44-46 .. code-block:: Python mean_curv = 0.5 * (PrincipalCurvatures[0, :] + PrincipalCurvatures[1, :]) .. GENERATED FROM PYTHON SOURCE LINES 47-50 Decomposition of the curvatures into ShapeIndex and Curvedness Based on 'Surface shape and curvature scales Jan JKoenderink & Andrea Jvan Doorn' .. GENERATED FROM PYTHON SOURCE LINES 50-52 .. code-block:: Python shapeIndex, curvedness = scurv.decompose_curvature(PrincipalCurvatures) .. GENERATED FROM PYTHON SOURCE LINES 53-54 Estimation error on the principal curvature length .. GENERATED FROM PYTHON SOURCE LINES 54-67 .. code-block:: Python K = [1, 0] quadric = sgps.generate_quadric( K, nstep=[20, 20], ax=3, ay=3, random_sampling=False, ratio=0.3, random_distribution_type="gamma", equilateral=True, ) .. GENERATED FROM PYTHON SOURCE LINES 68-69 Estimated computation of the Principal curvature, K_gauss, K_mean .. GENERATED FROM PYTHON SOURCE LINES 69-77 .. code-block:: Python p_curv, d1_estim, d2_estim = scurv.curvatures_and_derivatives(quadric) k1_estim, k2_estim = p_curv[0, :], p_curv[1, :] k_gauss_estim = k1_estim * k2_estim k_mean_estim = 0.5 * (k1_estim + k2_estim) .. rst-class:: sphx-glr-script-out .. code-block:: none Calculating vertex normals .... Please wait Finished calculating vertex normals Calculating curvature tensors ... Please wait Finished Calculating curvature tensors Calculating Principal Components ... Please wait Finished Calculating principal components .. GENERATED FROM PYTHON SOURCE LINES 78-79 Analytical computation of the curvatures .. GENERATED FROM PYTHON SOURCE LINES 79-97 .. code-block:: Python k_mean_analytic = sgps.quadric_curv_mean(K)( np.array(quadric.vertices[:, 0]), np.array(quadric.vertices[:, 1]) ) k_gauss_analytic = sgps.quadric_curv_gauss(K)( np.array(quadric.vertices[:, 0]), np.array(quadric.vertices[:, 1]) ) k1_analytic = np.zeros((len(k_mean_analytic))) k2_analytic = np.zeros((len(k_mean_analytic))) for i in range(len(k_mean_analytic)): a, b = np.roots((1, -2 * k_mean_analytic[i], k_gauss_analytic[i])) k1_analytic[i] = min(a, b) k2_analytic[i] = max(a, b) .. GENERATED FROM PYTHON SOURCE LINES 98-99 Error computation .. GENERATED FROM PYTHON SOURCE LINES 99-107 .. code-block:: Python k_mean_relative_change = abs( (k_mean_analytic - k_mean_estim) / k_mean_analytic) k_mean_absolute_change = abs((k_mean_analytic - k_mean_estim)) k1_relative_change = abs((k1_analytic - k1_estim) / k1_analytic) k1_absolute_change = abs((k1_analytic - k1_estim)) .. GENERATED FROM PYTHON SOURCE LINES 108-112 Estimation error on the curvature directions commented because there is a bug: ValueError: shapes (3,2) and (3,2) not aligned: 2 (dim 1) != 3 (dim 0) actually, vec1.shape=(3,) while vec2.shape=(3,2) .. GENERATED FROM PYTHON SOURCE LINES 112-126 .. code-block:: Python K = [1, 0] quadric = sgps.generate_quadric( K, nstep=[20, 20], ax=3, ay=3, random_sampling=False, ratio=0.3, random_distribution_type="gamma", equilateral=True, ) .. GENERATED FROM PYTHON SOURCE LINES 127-128 Estimated computation of the Principal curvature, Direction1, Direction2 .. GENERATED FROM PYTHON SOURCE LINES 128-130 .. code-block:: Python p_curv_estim, d1_estim, d2_estim = scurv.curvatures_and_derivatives(quadric) .. rst-class:: sphx-glr-script-out .. code-block:: none Calculating vertex normals .... Please wait Finished calculating vertex normals Calculating curvature tensors ... Please wait Finished Calculating curvature tensors Calculating Principal Components ... Please wait Finished Calculating principal components .. GENERATED FROM PYTHON SOURCE LINES 131-132 Analytical computation of the directions .. GENERATED FROM PYTHON SOURCE LINES 132-149 .. code-block:: Python analytical_directions = sgps.compute_all_principal_directions_3D( K, quadric.vertices) estimated_directions = np.zeros(analytical_directions.shape) estimated_directions[:, :, 0] = d1_estim estimated_directions[:, :, 1] = d2_estim angular_error_0, dotprods = ut.compare_analytic_estimated_directions( analytical_directions[:, :, 0], estimated_directions[:, :, 0] ) angular_error_0 = 180 * angular_error_0 / np.pi angular_error_1, dotprods = ut.compare_analytic_estimated_directions( analytical_directions[:, :, 1], estimated_directions[:, :, 1] ) angular_error_1 = 180 * angular_error_1 / np.pi .. GENERATED FROM PYTHON SOURCE LINES 150-216 VISUALIZATION USING EXTERNAL TOOLS ############################################################################ import visbrain # visu using visbrain Plot mean curvature visb_sc = splt.visbrain_plot( mesh=mesh, tex=mean_curv, caption="mean curvature", cblabel="mean curvature" ) visb_sc.preview() ############################################################################ # Plot Gauss curvature visb_sc = splt.visbrain_plot( mesh=mesh, tex=gaussian_curv, caption="Gaussian curvature", cblabel="Gaussian curvature", cmap="hot", ) visb_sc.preview() ############################################################################## Plot of ShapeIndex and Curvedness visb_sc = splt.visbrain_plot( mesh=mesh, tex=shapeIndex, caption="ShapeIndex", cblabel="ShapeIndex", cmap="coolwarm", ) visb_sc.preview() visb_sc = splt.visbrain_plot( mesh=mesh, tex=curvedness, caption="Curvedness", cblabel="Curvedness", cmap="hot" ) visb_sc.preview() ############################################################################## Error plot visb_sc = splt.visbrain_plot( mesh=quadric, tex=k_mean_absolute_change, caption="K_mean absolute error", cblabel="K_mean absolute error", ) visb_sc.preview() ############################################################################### # Error plot visb_sc = splt.visbrain_plot( mesh=quadric, tex=angular_error_0, caption="Angular error 0", cblabel="Angular error 0", ) visb_sc.preview() visb_sc = splt.visbrain_plot( mesh=quadric, tex=angular_error_1, caption="Angular error 1", cblabel="Angular error 1", ) visb_sc.preview() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 4.652 seconds) .. _sphx_glr_download_auto_examples_example_curvature.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_curvature.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_curvature.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_curvature.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_