.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/example_dpf_Boucher.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_dpf_Boucher.py: .. _example_dpf_Boucher: =================================== Example of depth potential function in slam =================================== .. GENERATED FROM PYTHON SOURCE LINES 8-15 .. code-block:: Python # Authors: Julien Lefevre # 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 22-23 Import of modules .. GENERATED FROM PYTHON SOURCE LINES 23-29 .. code-block:: Python import slam.curvature as sc import slam.differential_geometry as sdg import trimesh import numpy as np from scipy.spatial import Delaunay .. GENERATED FROM PYTHON SOURCE LINES 30-34 Define the example provided in Figure 5 of Depth potential function for folding pattern representation, registration and analysis Maxime Boucher a,b, * , Sue Whitesides a , Alan Evans .. GENERATED FROM PYTHON SOURCE LINES 34-81 .. code-block:: Python def boucher_surface(params, ax, ay, nstep): # Parameters xmin, xmax = [-ax, ax] ymin, ymax = [-ay, ay] # Define the sampling stepx = (xmax - xmin) / nstep stepy = stepx * np.sqrt(3) / 2 # to ensure equilateral faces # Coordinates x = np.arange(xmin, xmax, stepx) y = np.arange(ymin, ymax, stepy) X, Y = np.meshgrid(x, y) X[::2] += stepx / 2 X = X.flatten() Y = Y.flatten() # Delaunay faces_tri = Delaunay(np.vstack((X, Y)).T, qhull_options="QJ Qt Qbb") # Equation for Z M = params[0] sigma = params[1] # called sigma_y in the paper Z = ( M / sigma * np.exp(-(X**2) - Y**2 / (2 * sigma**2)) * (Y**2 - sigma**2) ) # Mesh coords = np.array([X, Y, Z]).transpose() mesh = trimesh.Trimesh( faces=faces_tri.simplices, vertices=coords, process=False) return mesh params = [4, 0.25] ax = 2 ay = 1 nstep = 50 mesh = boucher_surface(params, ax, ay, nstep) .. GENERATED FROM PYTHON SOURCE LINES 82-83 Compute dpf for various alpha .. GENERATED FROM PYTHON SOURCE LINES 83-97 .. code-block:: Python res = sc.curvatures_and_derivatives(mesh) mean_curvature = res[0].sum(axis=0) alphas = [0.001, 0.01, 0.1, 1, 10, 100] dpfs = sdg.depth_potential_function( mesh, curvature=mean_curvature, alphas=alphas) amplitude_center = [] amplitude_peak = [] index_peak_pos = np.argmax(mesh.vertices[:, 2]) index_peak_neg = np.argmin(mesh.vertices[:, 2]) for i in range(len(dpfs)): amplitude_center.append(dpfs[i][index_peak_neg]) amplitude_peak.append(dpfs[i][index_peak_pos]) .. 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 324 = 3.7465309898242367 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 .. GENERATED FROM PYTHON SOURCE LINES 98-99 Fix alpha and vary M = params[0] .. GENERATED FROM PYTHON SOURCE LINES 99-110 .. code-block:: Python all_M = [0.5, 0.75, 1, 1.25, 1.5, 1.75, 2] all_amplitudes = [] for M in all_M: mesh = boucher_surface([M, 0.25], ax, ay, nstep) res = sc.curvatures_and_derivatives(mesh) mean_curvature = res[0].sum(axis=0) dpfs = sdg.depth_potential_function( mesh, curvature=mean_curvature, alphas=[0.0015]) all_amplitudes.append(dpfs[0][len(mesh.vertices) // 2]) .. 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 110 = 1.2719703977798336 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 112 = 1.2950971322849214 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 112 = 1.2950971322849214 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 112 = 1.2950971322849214 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 110 = 1.2719703977798336 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 110 = 1.2719703977798336 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 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 Computing Laplacian Computing mesh weights of type fem -edge length threshold needed for 0 values = 0.0 % -number of Nan in weights: 0 = 0.0 % -number of Negative values in weights: 110 = 1.2719703977798336 % -nb Nan in Laplacian : 0 -nb Inf in Laplacian : 0 .. GENERATED FROM PYTHON SOURCE LINES 111-155 VISUALIZATION USING EXTERNAL TOOLS ############################################################################ import visbrain # visu using visbrain import slam.plot as splt import matplotlib.pyplot as plt ######################################################################### # Visualization of the mesh visb_sc = splt.visbrain_plot( mesh=mesh, caption="Boucher mesh", bgcolor=[ 0.3, 0.5, 0.7]) visb_sc visb_sc.preview() ############################################################################## plt.figure() plt.semilogx(alphas, amplitude_center) plt.semilogx(alphas, amplitude_peak) plt.semilogx(alphas, len(alphas) * [params[0] * (1 + 2 * np.exp(-3 / 2))], "--") plt.xlabel("alpha") plt.ylabel("amplitude") plt.legend(["DPF at center", "DPF (secondary peaks)", "True amplitude"]) plt.show() ##################################### # Display dpfs on the surfaces visb_sc = splt.visbrain_plot( mesh=mesh, tex=dpfs[0], caption="Boucher mesh", bgcolor="white" ) visb_sc = splt.visbrain_plot( mesh=mesh, tex=dpfs[5], caption="Boucher mesh", visb_sc=visb_sc ) visb_sc.preview() ############################################################################# plt.figure() plt.plot(all_M, all_amplitudes, "+-") plt.xlabel("M") plt.ylabel("Amplitude of DPF") plt.show() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 41.308 seconds) .. _sphx_glr_download_auto_examples_example_dpf_Boucher.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_dpf_Boucher.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_dpf_Boucher.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_dpf_Boucher.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_