.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plot_dendrite.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_plot_dendrite.py: Dendrite-Centric Coordinate System ================================== This example demonstrates how to fit a spline to a dendrite and align the image coordinate system with the spline. Data Source: `DeepD3 `_ Data Used: Crop of DeepD3_Benchmark.tif .. GENERATED FROM PYTHON SOURCE LINES 10-23 .. code-block:: Python import matplotlib.animation import matplotlib.pyplot as plt import numpy as np import pyvista as pv import scipy import skan import skimage import splinebox splinebox_color = plt.rcParams["axes.prop_cycle"].by_key()["color"][0] .. GENERATED FROM PYTHON SOURCE LINES 25-28 1. Load and Inspect the Data ---------------------------- We begin by loading the TIFF data, then visualize the image stack through the z-axis. .. GENERATED FROM PYTHON SOURCE LINES 28-43 .. code-block:: Python img = skimage.io.imread("dendrite.tif") def _update0(i): mpl_img.set_array(img[i]) return (mpl_img,) fig, ax = plt.subplots(figsize=(7, 3)) mpl_img = ax.imshow(img[0], cmap="Greys_r", vmin=160, vmax=2000) ax.set(xlim=(0, img.shape[2]), ylim=(0, img.shape[1])) animation = matplotlib.animation.FuncAnimation(fig, _update0, len(img), interval=100, blit=True) plt.show() .. container:: sphx-glr-animation .. raw:: html

.. GENERATED FROM PYTHON SOURCE LINES 44-47 **Visualize in 3D** Since we are fitting a 3D spline, let's visualize the image stack in 3D. .. GENERATED FROM PYTHON SOURCE LINES 47-57 .. code-block:: Python grid = pv.ImageData(dimensions=np.array(img.shape) + 1) grid.origin = (0, 0, 0) grid.spacing = (1, 1, 1) grid.cell_data["values"] = img.flatten(order="F") plotter = pv.Plotter() plotter.add_volume(grid, cmap="bone", clim=(160, 2000)) plotter.camera_position = "yz" plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_plot_dendrite_002.png :alt: plot dendrite :srcset: /auto_examples/images/sphx_glr_plot_dendrite_002.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/splinebox/checkouts/v0.4.1/docs/auto_examples/images/sphx_glr_plot_dendrite_002.vtksz .. GENERATED FROM PYTHON SOURCE LINES 58-63 NOTE: if you don't see the image after switching to 'Interactive Scene' you might have to click and drag on the white space once. 2. Segmentation and Skeletonization ----------------------------------- We segment the dendrite using Otsu's method and skeletonize it to obtain the pixel coordinates for spline fitting. .. GENERATED FROM PYTHON SOURCE LINES 63-85 .. code-block:: Python thresh = skimage.filters.threshold_otsu(img) mask = img > thresh # Keep only the largest connected component label_img = skimage.measure.label(mask) label_biggest = np.argmax(np.bincount(label_img.flatten())[1:]) + 1 mask = label_img == label_biggest # Skeletonize skeleton = skimage.morphology.skeletonize(mask) # Get skeleton coordinates skeleton_points = np.stack(np.where(skeleton), axis=-1) # Visualize skeleton skeleton_point_cloud = pv.PolyData(skeleton_points.astype(float)) plotter = pv.Plotter() plotter.add_volume(grid, cmap="bone", clim=(160, 2000)) plotter.add_mesh(skeleton_point_cloud, color=splinebox_color, point_size=10, render_points_as_spheres=True) plotter.camera_position = "yz" plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_plot_dendrite_003.png :alt: plot dendrite :srcset: /auto_examples/images/sphx_glr_plot_dendrite_003.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/splinebox/checkouts/v0.4.1/docs/auto_examples/images/sphx_glr_plot_dendrite_003.vtksz .. GENERATED FROM PYTHON SOURCE LINES 86-89 **Extract the Longest Path of the Skeleton** We convert the skeleton into a graph and extract the longest path, which corresponds to the main dendrite. .. GENERATED FROM PYTHON SOURCE LINES 89-119 .. code-block:: Python # Returns a sparse connectivity matrix (graph) and the corresponding pixel coordinates # for each knot in the graph. graph, coords = skan.csr.skeleton_to_csgraph(skeleton) coords = np.stack(coords, axis=-1) # Use the shortest path algorithm to find the distances between all knots(skeleton points). dist_matrix, predecessors = scipy.sparse.csgraph.shortest_path(graph, return_predecessors=True) # Extract the index of the start and end knots of the longest path start_index, stop_index = np.unravel_index(np.argmax(dist_matrix), dist_matrix.shape) # Reconstruct the longest path using the predecessor matrix i = start_index skeleton_points = [] while i != stop_index: skeleton_points.append(coords[i]) i = predecessors[stop_index, i] skeleton_points = np.array(skeleton_points) # Visualize longest path skeleton_point_cloud = pv.PolyData(skeleton_points.astype(float)) plotter = pv.Plotter() plotter.add_volume(grid, cmap="bone", clim=(160, 2000)) plotter.add_mesh(skeleton_point_cloud, color=splinebox_color, point_size=10, render_points_as_spheres=True) plotter.camera_position = "yz" plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_plot_dendrite_004.png :alt: plot dendrite :srcset: /auto_examples/images/sphx_glr_plot_dendrite_004.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/splinebox/checkouts/v0.4.1/docs/auto_examples/images/sphx_glr_plot_dendrite_004.vtksz .. GENERATED FROM PYTHON SOURCE LINES 120-123 3. Fit a Spline --------------- Now that we have the main points of the dendrite, we fit a spline. .. GENERATED FROM PYTHON SOURCE LINES 123-129 .. code-block:: Python M = 20 basis_function = splinebox.basis_functions.B3() spline = splinebox.spline_curves.Spline(M=M, basis_function=basis_function, closed=False) spline.fit(skeleton_points) .. GENERATED FROM PYTHON SOURCE LINES 130-133 4. Plot the Fitted Spline ------------------------- Let's visualize the spline and its knots along with the segmented dendrite. .. GENERATED FROM PYTHON SOURCE LINES 133-154 .. code-block:: Python # Creat meshes for the spline and the knots of the spline t = np.linspace(0, M - 1, M * 15) spline_mesh = pv.MultipleLines(points=spline(t)) knots_point_cloud = pv.PolyData(spline.knots) # Prepare segmentation mesh grid = pv.ImageData(dimensions=mask.shape) mesh = grid.contour([0.5], mask.flatten(order="F"), method="marching_cubes") mesh = mesh.clean() mesh = mesh.decimate(0.98) mesh = mesh.smooth(100) plotter = pv.Plotter() plotter.add_mesh(mesh, style="wireframe", color="black") plotter.add_mesh(spline_mesh, color=splinebox_color, line_width=10) plotter.add_mesh(knots_point_cloud, color="red", point_size=10, render_points_as_spheres=True) plotter.camera_position = "yz" plotter.zoom_camera(2) plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_plot_dendrite_005.png :alt: plot dendrite :srcset: /auto_examples/images/sphx_glr_plot_dendrite_005.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/splinebox/checkouts/v0.4.1/docs/auto_examples/images/sphx_glr_plot_dendrite_005.vtksz .. GENERATED FROM PYTHON SOURCE LINES 155-160 5. Compute Normal Planes ------------------------ To align the image coordinate system with the spline, we compute normal planes along the spline. The normal planes are spanned by two vectors, which are normal to the local derivative vector of the spline (i.e. the vector pointing in the local direction of the spline). .. GENERATED FROM PYTHON SOURCE LINES 160-164 .. code-block:: Python # Compute derivative vectors along the spline deriv = spline(t, derivative=1) .. GENERATED FROM PYTHON SOURCE LINES 165-169 **Compute the first normal vector** We select the normal vector that lies in the x-y plane (i.e. we set the z component to zero). .. GENERATED FROM PYTHON SOURCE LINES 169-173 .. code-block:: Python normal1 = np.zeros((len(t), 3)) normal1[:, 1] = deriv[:, 2] normal1[:, 2] = -deriv[:, 1] .. GENERATED FROM PYTHON SOURCE LINES 174-179 **Compute the second normal vector** The second normal vector can be obtained using the cross product. This yields a vector that is perpendicular to the two input vectors ``deriv`` and ``normal1``. .. GENERATED FROM PYTHON SOURCE LINES 179-186 .. code-block:: Python normal2 = np.zeros((len(t), 3)) normal2 = np.cross(deriv, normal1) # Normalize vectors normal1 /= np.linalg.norm(normal1, axis=1)[:, np.newaxis] normal2 /= np.linalg.norm(normal2, axis=1)[:, np.newaxis] .. GENERATED FROM PYTHON SOURCE LINES 187-190 **Visualize Normal Planes** We scale the vectors for better visibility and plot them to verify they are perpendicular to the spline. .. GENERATED FROM PYTHON SOURCE LINES 190-203 .. code-block:: Python spline_mesh["normal1"] = normal1 * 7 spline_mesh["normal2"] = normal2 * 7 plotter = pv.Plotter() spline_mesh.set_active_vectors("normal1") plotter.add_mesh(spline_mesh.arrows, lighting=False, color="black") spline_mesh.set_active_vectors("normal2") plotter.add_mesh(spline_mesh.arrows, lighting=False, color="red") plotter.add_mesh(spline_mesh, color=splinebox_color, line_width=10) plotter.camera_position = "yz" plotter.zoom_camera(2) plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_plot_dendrite_006.png :alt: plot dendrite :srcset: /auto_examples/images/sphx_glr_plot_dendrite_006.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/splinebox/checkouts/v0.4.1/docs/auto_examples/images/sphx_glr_plot_dendrite_006.vtksz .. GENERATED FROM PYTHON SOURCE LINES 204-207 **Extract Pixel Values in Normal Planes** Finally, we interpolate pixel values from the original image along the computed normal planes. .. GENERATED FROM PYTHON SOURCE LINES 207-243 .. code-block:: Python # Centers of the normal planes spline_coordinates = spline(t) # Coefficients for scaling the normal vectors half_window_size = 25 window_range = np.arange(-half_window_size, half_window_size) ii, jj = np.meshgrid(window_range, window_range) # Compute pixel coordinates using scaled normal vectors normal_planes = np.multiply.outer(ii, normal1) + np.multiply.outer(jj, normal2) # Fix the order of the axes (spline position first, before the normal directions) normal_planes = np.rollaxis(normal_planes, 2, 0) # Position normal planes on spline normal_planes += spline_coordinates[:, np.newaxis, np.newaxis] # Interpolate pixel values shape = normal_planes.shape vals = scipy.ndimage.map_coordinates( img, normal_planes.reshape(-1, 3).T, order=1, ) vals = vals.reshape(shape[:-1]).astype(np.float64) # Mask out pixels outside the volume mask = ( (np.min(normal_planes, axis=3) < 0) | (normal_planes[:, :, :, 0] > img.shape[0] - 1) | (normal_planes[:, :, :, 1] > img.shape[1] - 1) | (normal_planes[:, :, :, 2] > img.shape[2] - 1) ) vals[mask] = np.nan .. GENERATED FROM PYTHON SOURCE LINES 244-247 6. Animate the Dendrite-Centric Image ------------------------------------- We create an animation showing the dendrite as seen along the fitted spline. .. GENERATED FROM PYTHON SOURCE LINES 247-268 .. code-block:: Python def _update1(i): mpl_point.set_offsets( spline_coordinates[ i, 2:0:-1, ].T ) mpl_img.set_array(vals[i]) return (mpl_point, mpl_img) fig, axes = plt.subplots(1, 2, figsize=(7, 3)) axes[0].imshow(np.max(img, axis=0), cmap="Greys_r", vmin=160, vmax=2000) mpl_point = axes[0].scatter((spline_coordinates[0, 2],), (spline_coordinates[0, 1],)) mpl_img = axes[1].imshow(vals[0], cmap="Greys_r", vmin=160, vmax=2000) axes[1].set(xlim=(0, vals.shape[2]), ylim=(0, vals.shape[1])) axes[1].scatter((half_window_size,), (half_window_size,)) animation = matplotlib.animation.FuncAnimation(fig, _update1, len(vals), interval=100, blit=True) plt.show() .. container:: sphx-glr-animation .. raw:: html

.. rst-class:: sphx-glr-timing **Total running time of the script:** (1 minutes 12.697 seconds) .. _sphx_glr_download_auto_examples_plot_dendrite.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_dendrite.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_dendrite.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_dendrite.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_