feature(measure): Adds the height_profiles sub-module

[ns-rse]

Further work on height profiles this PR adds the topostats.measure.height_profiles sub-module which for a given pair of feret coordinates...

1. Determines the orientation of of the feret co-ordinates relative to horizontal.
2. Rotates the image so that the feret is horizontal.
3. Recalculates the co-ordinates of the feret after rotation.
4. Extracts the height profile along the feret.

Includes a function in topostats.plotting.plot_height_profiles() which produces a line-plot of multiple height profiles.

Test shapes defined tests/measure/test_feret.py moved to fixtures in tests/measure/conftest.py to avoid duplication and tests for all functions are included. Because these fixtures are then used in @pytest.mark.parameterize() it has necessitated using the request fixture and its .getfixturevalue() method.

More on this can be read in this blog post.

ToDo...

Still a fair few steps to integrate this into the processing.

• Add configuration options to topostats/default_config.yaml of whether to calculate height_profiles.
• Update GrainStats() to calculate the height profile for the image being processed if required based on configuration.
• Return the height_profile (1-D Numpy array) as part of the tuple GrainStats() returns.
• Collect height_profile across grains into a dictionary (may require current code as written to be adapted to work with dictionaries, currently works with lists in plot_height_profiles()).
• Add functionality to write the profiles to JSON for subsequent use/plotting (e.g. customising style/axis labels/etc.
• Currently only configured to give height profiles across maximum feret, would need additional work and consideration for profiles perpendicular to the mid-point of the maximum feret (but code is structured such that the extraction of this would be relatively easy, but would require additional work on plotting and returning data).

Related : https://github.com/AFM-SPM/TopoStats/issues/748 https://github.com/AFM-SPM/TopoStats/pull/755

[ns-rse]

@llwiggins has identified a scenario where behaviour is not as expected.

import numpy as np
from topostats.measure import height_profiles

circle = np.array([
    [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0],
    [0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0],
    [0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
    [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
    [0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0],
    [0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0],
    [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0]
])

binary_mask = circle

height_profiles.extract_feret_profiles(circle, binary_mask)

This fails the following test...

import pytest
import numpy.typing as npt

@pytest.mark.parametrize(
    ("img", "target"),
    [
        pytest.param(
            {
                "img": np.asarray(
                    [
                        [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0],
                        [0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0],
                        [0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0],
                        [0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0],
                        [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0],
                    ]
                ),
                "skeleton": np.asarray(
                    [
                        [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0],
                        [0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0],
                        [0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0],
                        [0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0],
                        [0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0],
                        [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0],
                    ]
                ),
            },
            np.asarray([1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]),
        ),
    ],
)
def test_extract_feret_profiles2(img: npt.NDArray, target: npt.NDArray, request) -> None:
    """Test extract_feret_profiles()."""
    # _img = request.getfixturevalue(img)
    # Convert boolean to 0/1 and extract profiles
    feret_profiles = height_profiles.extract_feret_profiles(img=img["img"], skeleton=np.where(img["skeleton"], 1, 0))
    # This fails...
    # feret_profiles = height_profiles.extract_feret_profiles(img=_img["img"], skeleton=_img["skeleton"])
    np.testing.assert_array_almost_equal(feret_profiles, target, decimal=21)

E           AssertionError: 
E           Arrays are not almost equal to 21 decimals
E           
E           Mismatched elements: 10 / 11 (90.9%)
E           Max absolute difference: 1
E           Max relative difference: 1.
E            x: array([0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0])
E            y: array([1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1])

Why does this fail?

Whilst the image doesn't need rotating this fails because the maximum feret co-ordinates after rotating (which wasn't in effect performed) are a vertical line rather than the expected horizontal line...

rotated_feret_stats={'max_feret': 10.0, 
                     'min_feret': 8.48528137423857, 
                     'max_feret_coords': array([[10,  5],
                                                [ 0,  5]]),
                     'min_feret_coords': array([[1., 3.],
                                                [7., 9.]])}

And the profile() function is very basic, taking an image and the row from which the profile is to be returned...

def profile(img: npt.NDArray, row: int) -> npt.NDArray:
    """
    ...
    """
    try:
        return img[row]
    except IndexError as e:
        raise IndexError("The slice (row) is outside the image boundary.") from e

And so because the row being used is rotated_feret_stats["max_feret_coords"][1, 0] a slice is taken across the first row which is why we see the returned value.

[ns-rse]

@llwiggins I've revised how the profiles are obtained and ditched the rotation method instead using scipy.interpolate.RegularGridInterpolator.

As a consequence the topostats.measure.heights_profile sub-module is considerably smaller and consists of a single function interpolate_height_profiles() which performs interpolation of heights directly between two points by first getting evenly spaced points on the x-axis between the maximum feret co-ordinates and then interpolating the corresponding y-coordinates. These x and y co-ordinates then have the height interpolated using RegularGridInterpolator.

Arguments can be passed via **kwargs to RegularGridInterpolator() and the effect of modifying the method is tested. I've included your example of a basic circle which caused problems previously and those tests pass now (whenever a scenario that causes failure is encountered its a good candidate for a test so thank you).

I've found for some reason that interpolation (using the method = linear) on M$-Win gives different results so have marked those tests to be conditionally skipped if running on a Windows system.

[SylviaWhittle]

Been reading the PR and looks good to me. Much prefer the interpolation method, so thank you. Not given it a thorough enough read & tinker with to approve yet, but I like the method and implementation! I'll test it out soon (once my report is done at the latest). Apologies for not being active much on TopoStats recently.