diff --git a/kwave/utils/mapgen.py b/kwave/utils/mapgen.py index 3c957be96..4b8a76a1d 100644 --- a/kwave/utils/mapgen.py +++ b/kwave/utils/mapgen.py @@ -954,371 +954,279 @@ def create_pixel_dim(Nx: int, origin_size: float, shift: float) -> Tuple[np.ndar @typechecker -def make_line( +def make_line_straight( grid_size: Vector, startpoint: Union[Tuple[Int[kt.ScalarLike, ""], Int[kt.ScalarLike, ""]], Int[np.ndarray, "2"]], endpoint: Optional[Union[Tuple[Int[kt.ScalarLike, ""], Int[kt.ScalarLike, ""]], Int[np.ndarray, "2"]]] = None, - angle: Optional[Float[kt.ScalarLike, ""]] = None, - length: Optional[Int[kt.ScalarLike, ""]] = None, ) -> kt.NP_ARRAY_BOOL_2D: - """ - Generate a line shape with a given start and end point, angle, or length. - - Args: - grid_size: The size of the grid in pixels. - startpoint: The start point of the line, given as a tuple of x and y coordinates. - endpoint: The end point of the line, given as a tuple of x and y coordinates. - If not specified, the line is drawn from the start point at a given angle and length. - angle: The angle of the line in radians, measured counterclockwise from the x-axis. - If not specified, the line is drawn from the start point to the end point. - length: The length of the line in pixels. - If not specified, the line is drawn from the start point to the end point. - - Returns: - line: A 2D array of the same size as the input parameters, - with a value of 1 for pixels that are part of the line and 0 for pixels that are not. - """ - assert len(grid_size) == 2, "Grid size must be a 2-element vector." - - startpoint = np.array(startpoint, dtype=int) - if endpoint is not None: - endpoint = np.array(endpoint, dtype=int) + a, b = np.array(startpoint) - 1, np.array(endpoint) - 1 # Adjust for zero-based indexing + + if np.all(a == b): + raise ValueError("The first and last points cannot be the same.") + if not (0 <= a).all() or not (b < grid_size).all(): + raise ValueError("Both the start and end points must lie within the grid.") + + line = np.zeros(grid_size, dtype=bool) # Initialize the grid to hold the line + m = (b[1] - a[1]) / (b[0] - a[0]) if b[0] != a[0] else np.inf # Calculate gradient, handle vertical line + c = a[1] - m * a[0] if np.isfinite(m) else 0 # Calculate y-intercept if m is finite + + # Generate line points based on slope + if abs(m) < 1 or np.isinf(m): # Handle more horizontal or perfectly vertical lines + x_start, x_end = sorted([a[0], b[0]]) + for x in range(x_start, x_end + 1): + y = round(m * x + c) if np.isfinite(m) else a[1] + line[x, y] = True + else: # Handle more vertical lines + y_start, y_end = sorted([a[1], b[1]]) + for y in range(y_start, y_end + 1): + x = round((y - c) / m) + line[x, y] = True - if len(startpoint) != 2: - raise ValueError("startpoint should be a two-element vector.") + return line - if np.any(startpoint < 1) or startpoint[0] > grid_size.x or startpoint[1] > grid_size.y: - ValueError("The starting point must lie within the grid, between [1 1] and [grid_size.x grid_size.y].") - # ========================================================================= - # LINE BETWEEN TWO POINTS OR ANGLED LINE? - # ========================================================================= +@typechecker +def make_line_angled( + grid_size: Vector, + startpoint: Union[Tuple[Int[kt.ScalarLike, ""], Int[kt.ScalarLike, ""]], Int[np.ndarray, "2"]], + angle: Optional[Float[kt.ScalarLike, ""]] = None, + length: Optional[Int[kt.ScalarLike, ""]] = None, +) -> kt.NP_ARRAY_BOOL_2D: + angle, linelength = angle, length - if endpoint is not None: - linetype = "AtoB" - a, b = startpoint, endpoint + # angle must lie between -np.pi and np.pi + angle = angle % (2 * np.pi) + if angle > np.pi: + angle = angle - (2 * np.pi) + elif angle < -np.pi: + angle = angle + (2 * np.pi) - # Addition => Fix Matlab2Python indexing - a -= 1 - b -= 1 - else: - linetype = "angled" - angle, linelength = angle, length + # define an empty grid to hold the line + line = np.zeros(grid_size, dtype=bool) - # ========================================================================= - # MORE INPUT CHECKING - # ========================================================================= + # start at the atart + x, y = startpoint - if linetype == "AtoB": - # a and b must be different points - if np.all(a == b): - raise ValueError("The first and last points cannot be the same.") - - # end point must be a two-element row vector - if len(b) != 2: - raise ValueError("endpoint should be a two-element vector.") - - # a and b must be within the grid - xx = np.array([a[0], b[0]], dtype=int) - yy = np.array([a[1], b[1]], dtype=int) - if np.any(a < 0) or np.any(b < 0) or np.any(xx > grid_size.x - 1) or np.any(yy > grid_size.y - 1): - raise ValueError("Both the start and end points must lie within the grid.") - - if linetype == "angled": - # angle must lie between -np.pi and np.pi - angle = angle % (2 * np.pi) - if angle > np.pi: - angle = angle - (2 * np.pi) - elif angle < -np.pi: - angle = angle + (2 * np.pi) + # fill in the first point + line[x - 1, y - 1] = 1 - # ========================================================================= - # CALCULATE A LINE FROM A TO B - # ========================================================================= + # initialise the current length of the line + line_length = 0 - if linetype == "AtoB": - # define an empty grid to hold the line - line = np.zeros(grid_size, dtype=bool) + if abs(angle) == np.pi: + # Calculate and draw the line to the maximum allowed y position within the grid and line length constraints + max_y = min(y + linelength, grid_size.y) + line[x - 1, y : max_y - 1] = True - # find the equation of the line - m = (b[1] - a[1]) / (b[0] - a[0]) # gradient of the line - c = a[1] - m * a[0] # where the line crosses the y axis + # Update y to the last position and calculate the actual length of the line + y = max_y - 1 + line_length = np.hypot(x - startpoint[0], y - startpoint[1]) - if abs(m) < 1: - # start at the end with the smallest value of x - if a[0] < b[0]: - x, y = a - x_end = b[0] - else: - x, y = b - x_end = a[0] - - # fill in the first point - line[x, y] = 1 - - while x < x_end: - # next points to try are - poss_x = [x, x, x + 1, x + 1, x + 1] - poss_y = [y - 1, y + 1, y - 1, y, y + 1] - - # find the point closest to the line - true_y = m * poss_x + c - diff = (poss_y - true_y) ** 2 - index = matlab_find(diff == min(diff))[0] - - # the next point - x = poss_x[index[0] - 1] - y = poss_y[index[0] - 1] - - # add the point to the line - line[x - 1, y - 1] = 1 - - elif not np.isinf(abs(m)): - # start at the end with the smallest value of y - if a[1] < b[1]: - x = a[0] - y = a[1] - y_end = b[1] - else: - x = b[0] - y = b[1] - y_end = a[1] - - # fill in the first point - line[x, y] = 1 - - while y < y_end: - # next points to try are - poss_y = [y, y, y + 1, y + 1, y + 1] - poss_x = [x - 1, x + 1, x - 1, x, x + 1] - - # find the point closest to the line - true_x = (poss_y - c) / m - diff = (poss_x - true_x) ** 2 - index = matlab_find(diff == min(diff))[0] - - # the next point - x = poss_x[index[0] - 1] - y = poss_y[index[0] - 1] - - # add the point to the line - line[x, y] = 1 - - else: # m = +-Inf - # start at the end with the smallest value of y - if a[1] < b[1]: - x = a[0] - y = a[1] - y_end = b[1] - else: - x = b[0] - y = b[1] - y_end = a[1] + elif (angle < np.pi) and (angle > np.pi / 2): + # define the equation of the line + m = -np.tan(angle - np.pi / 2) # gradient of the line + c = y - m * x # where the line crosses the y axis - # fill in the first point - line[x, y] = 1 + while line_length < linelength: + # next points to try are + poss_x = np.array([x - 1, x - 1, x]) + poss_y = np.array([y, y + 1, y + 1]) - while y < y_end: - # next point - y = y + 1 + # find the point closest to the line + true_y = m * poss_x + c + diff = (poss_y - true_y) ** 2 + index = matlab_find(diff == min(diff))[0] - # add the point to the line - line[x, y] = 1 + # the next point + x = poss_x[index[0] - 1] + y = poss_y[index[0] - 1] - # ========================================================================= - # CALCULATE AN ANGLED LINE - # ========================================================================= + # stop the points incrementing at the edges + if (x < 0) or (y > grid_size.y - 1): + break - elif linetype == "angled": - # define an empty grid to hold the line - line = np.zeros(grid_size, dtype=bool) + # add the point to the line + line[x - 1, y - 1] = 1 - # start at the atart - x, y = startpoint + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - # fill in the first point - line[x - 1, y - 1] = 1 + elif angle == np.pi / 2: + while line_length < linelength: + # next point + x = x - 1 - # initialise the current length of the line - line_length = 0 + # stop the points incrementing at the edges + if x < 1: + break - if abs(angle) == np.pi: - while line_length < linelength: - # next point - y = y + 1 + # add the point to the line + line[x - 1, y - 1] = 1 - # stop the points incrementing at the edges - if y > grid_size.y: - break + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - # add the point to the line - line[x - 1, y - 1] = 1 + elif (angle < np.pi / 2) and (angle > 0): + # define the equation of the line + m = np.tan(np.pi / 2 - angle) # gradient of the line + c = y - m * x # where the line crosses the y axis - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + while line_length < linelength: + # next points to try are + poss_x = np.array([x - 1, x - 1, x]) + poss_y = np.array([y, y - 1, y - 1]) - elif (angle < np.pi) and (angle > np.pi / 2): - # define the equation of the line - m = -np.tan(angle - np.pi / 2) # gradient of the line - c = y - m * x # where the line crosses the y axis + # find the point closest to the line + true_y = m * poss_x + c + diff = (poss_y - true_y) ** 2 + index = matlab_find(diff == min(diff))[0] - while line_length < linelength: - # next points to try are - poss_x = np.array([x - 1, x - 1, x]) - poss_y = np.array([y, y + 1, y + 1]) + # the next point + x = poss_x[index[0] - 1] + y = poss_y[index[0] - 1] - # find the point closest to the line - true_y = m * poss_x + c - diff = (poss_y - true_y) ** 2 - index = matlab_find(diff == min(diff))[0] + # stop the points incrementing at the edges + if (x < 1) or (y < 1): + break - # the next point - x = poss_x[index[0] - 1] - y = poss_y[index[0] - 1] + # add the point to the line + line[x - 1, y - 1] = 1 - # stop the points incrementing at the edges - if (x < 0) or (y > grid_size.y - 1): - break - - # add the point to the line - line[x - 1, y - 1] = 1 - - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - - elif angle == np.pi / 2: - while line_length < linelength: - # next point - x = x - 1 - - # stop the points incrementing at the edges - if x < 1: - break - - # add the point to the line - line[x - 1, y - 1] = 1 - - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + + elif angle == 0: + while line_length < linelength: + # next point + y = y - 1 - elif (angle < np.pi / 2) and (angle > 0): - # define the equation of the line - m = np.tan(np.pi / 2 - angle) # gradient of the line - c = y - m * x # where the line crosses the y axis + # stop the points incrementing at the edges + if y < 1: + break - while line_length < linelength: - # next points to try are - poss_x = np.array([x - 1, x - 1, x]) - poss_y = np.array([y, y - 1, y - 1]) + # add the point to the line + line[x - 1, y - 1] = 1 - # find the point closest to the line - true_y = m * poss_x + c - diff = (poss_y - true_y) ** 2 - index = matlab_find(diff == min(diff))[0] + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - # the next point - x = poss_x[index[0] - 1] - y = poss_y[index[0] - 1] + elif (angle < 0) and (angle > -np.pi / 2): + # define the equation of the line + m = -np.tan(np.pi / 2 + angle) # gradient of the line + c = y - m * x # where the line crosses the y axis - # stop the points incrementing at the edges - if (x < 1) or (y < 1): - break + while line_length < linelength: + # next points to try are + poss_x = np.array([x + 1, x + 1, x]) + poss_y = np.array([y, y - 1, y - 1]) - # add the point to the line - line[x - 1, y - 1] = 1 + # find the point closest to the line + true_y = m * poss_x + c + diff = (poss_y - true_y) ** 2 + index = matlab_find(diff == min(diff))[0] - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + # the next point + x = poss_x[index[0] - 1] + y = poss_y[index[0] - 1] - elif angle == 0: - while line_length < linelength: - # next point - y = y - 1 + # stop the points incrementing at the edges + if (x > grid_size.x) or (y < 1): + break - # stop the points incrementing at the edges - if y < 1: - break + # add the point to the line + line[x - 1, y - 1] = 1 - # add the point to the line - line[x - 1, y - 1] = 1 + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + elif angle == -np.pi / 2: + while line_length < linelength: + # next point + x = x + 1 - elif (angle < 0) and (angle > -np.pi / 2): - # define the equation of the line - m = -np.tan(np.pi / 2 + angle) # gradient of the line - c = y - m * x # where the line crosses the y axis + # stop the points incrementing at the edges + if x > grid_size.x: + break - while line_length < linelength: - # next points to try are - poss_x = np.array([x + 1, x + 1, x]) - poss_y = np.array([y, y - 1, y - 1]) + # add the point to the line + line[x - 1, y - 1] = 1 - # find the point closest to the line - true_y = m * poss_x + c - diff = (poss_y - true_y) ** 2 - index = matlab_find(diff == min(diff))[0] + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - # the next point - x = poss_x[index[0] - 1] - y = poss_y[index[0] - 1] + elif (angle < -np.pi / 2) and (angle > -np.pi): + # define the equation of the line + m = np.tan(-angle - np.pi / 2) # gradient of the line + c = y - m * x # where the line crosses the y axis - # stop the points incrementing at the edges - if (x > grid_size.x) or (y < 1): - break + while line_length < linelength: + # next points to try are + poss_x = np.array([x + 1, x + 1, x]) + poss_y = np.array([y, y + 1, y + 1]) - # add the point to the line - line[x - 1, y - 1] = 1 + # find the point closest to the line + true_y = m * poss_x + c + diff = (poss_y - true_y) ** 2 + index = matlab_find(diff == min(diff))[0] - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + # the next point + x = poss_x[index[0] - 1] + y = poss_y[index[0] - 1] - elif angle == -np.pi / 2: - while line_length < linelength: - # next point - x = x + 1 + # stop the points incrementing at the edges + if (x > grid_size.x) or (y > grid_size.y): + break - # stop the points incrementing at the edges - if x > grid_size.x: - break + # add the point to the line + line[x - 1, y - 1] = 1 - # add the point to the line - line[x - 1, y - 1] = 1 + # calculate the current length of the line + line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + return line - elif (angle < -np.pi / 2) and (angle > -np.pi): - # define the equation of the line - m = np.tan(-angle - np.pi / 2) # gradient of the line - c = y - m * x # where the line crosses the y axis - while line_length < linelength: - # next points to try are - poss_x = np.array([x + 1, x + 1, x]) - poss_y = np.array([y, y + 1, y + 1]) +@typechecker +def make_line( + grid_size: Vector, + startpoint: Union[Tuple[Int[kt.ScalarLike, ""], Int[kt.ScalarLike, ""]], Int[np.ndarray, "2"]], + endpoint: Optional[Union[Tuple[Int[kt.ScalarLike, ""], Int[kt.ScalarLike, ""]], Int[np.ndarray, "2"]]] = None, + angle: Optional[Float[kt.ScalarLike, ""]] = None, + length: Optional[Int[kt.ScalarLike, ""]] = None, +) -> kt.NP_ARRAY_BOOL_2D: + """ + Generate a line shape with a given start and end point, angle, or length. - # find the point closest to the line - true_y = m * poss_x + c - diff = (poss_y - true_y) ** 2 - index = matlab_find(diff == min(diff))[0] + Args: + grid_size: The size of the grid in pixels. + startpoint: The start point of the line, given as a tuple of x and y coordinates. + endpoint: The end point of the line, given as a tuple of x and y coordinates. + If not specified, the line is drawn from the start point at a given angle and length. + angle: The angle of the line in radians, measured counterclockwise from the x-axis. + If not specified, the line is drawn from the start point to the end point. + length: The length of the line in pixels. + If not specified, the line is drawn from the start point to the end point. - # the next point - x = poss_x[index[0] - 1] - y = poss_y[index[0] - 1] + Returns: + line: A 2D array of the same size as the input parameters, + with a value of 1 for pixels that are part of the line and 0 for pixels that are not. + """ + assert len(grid_size) == 2, "Grid size must be a 2-element vector." - # stop the points incrementing at the edges - if (x > grid_size.x) or (y > grid_size.y): - break + startpoint = np.array(startpoint, dtype=int) + if endpoint is not None: + endpoint = np.array(endpoint, dtype=int) - # add the point to the line - line[x - 1, y - 1] = 1 + if len(startpoint) != 2: + raise ValueError("startpoint should be a two-element vector.") - # calculate the current length of the line - line_length = np.sqrt((x - startpoint[0]) ** 2 + (y - startpoint[1]) ** 2) + if np.any(startpoint < 1) or startpoint[0] > grid_size.x or startpoint[1] > grid_size.y: + ValueError("The starting point must lie within the grid, between [1 1] and [grid_size.x grid_size.y].") - return line + if endpoint is not None: + return make_line_straight(grid_size=grid_size, startpoint=startpoint, endpoint=endpoint) + else: + return make_line_angled(grid_size=grid_size, startpoint=startpoint, angle=angle, length=length) @typechecker diff --git a/tests/matlab_test_data_collectors/matlab_collectors/collect_makeLine.m b/tests/matlab_test_data_collectors/matlab_collectors/collect_makeLine.m index 6936243ed..1bb2a2037 100644 --- a/tests/matlab_test_data_collectors/matlab_collectors/collect_makeLine.m +++ b/tests/matlab_test_data_collectors/matlab_collectors/collect_makeLine.m @@ -2,6 +2,8 @@ all_params = { ... {20, 20, [1, 1], [19, 19]}, ... {20, 20, [1, 1], [5, 5]}, ... + {20, 20, [1, 1], [5, 1]}, ... + {20, 20, [1, 1], [1, 5]}, ... {10, 20, [4, 15], [10, 1]}, ... {20, 20, [15, 15], 0.75 * pi, 5}, ... {20, 20, [15, 15], 0.75 * pi, -5}, ... diff --git a/tests/matlab_test_data_collectors/python_testers/makeLine_test.py b/tests/matlab_test_data_collectors/python_testers/makeLine_test.py index 7197dccc3..3b98fe99b 100644 --- a/tests/matlab_test_data_collectors/python_testers/makeLine_test.py +++ b/tests/matlab_test_data_collectors/python_testers/makeLine_test.py @@ -3,6 +3,7 @@ from pathlib import Path import numpy as np +import pytest from scipy.io import loadmat from kwave.data import Vector @@ -43,3 +44,17 @@ def test_makeLine(): assert np.allclose(expected_line, line) logging.log(logging.INFO, "make_line(..) works as expected!") + + +def test_start_greater_grid_size(): + with np.testing.assert_raises(ValueError): + make_line(Vector([1, 1]), (-10, 10), (10, 10)) + + +def test_a_b_same(): + with np.testing.assert_raises(ValueError): + make_line(Vector([1, 1]), (10, 10), (10, 10)) + + +if __name__ == "__main__": + pytest.main([__file__])