Skip to content

Commit

Permalink
ENH: updated anchor point finding
Browse files Browse the repository at this point in the history
  • Loading branch information
@marklescroart
marklescroart committed May 16, 2024
1 parent a424150 commit 8c28079
Showing 1 changed file with 93 additions and 44 deletions.
137 changes: 93 additions & 44 deletions cortex/dartboards.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -160,6 +160,7 @@ def vertical_flip_path(path, overlay):


def determine_path_hemi(overlay, path, percent_cutoff=.4):
"""Returns 0 (False) for left, 1 (True) for right hemisphere"""
if isinstance(path, matplotlib.path.Path):
path = path.vertices
middle_x = overlay.svgshape[0]/2
Expand All @@ -178,7 +179,9 @@ def center_of_mass(X):
-------
center_of_mass
array of shape (2,) that contains (x, y) center of mass
"""
"""
if not np.all(np.isclose(X[0], X[-1])):
X = np.vstack([X, X[0]])
# calculate center of mass of a closed polygon
x = X[:, 0]
y = X[:, 1]
Expand Down Expand Up @@ -378,6 +381,7 @@ def _angles_from_vertex(overlay, start_idx, target_idxs=None, period=2*np.pi, th


def get_roi_paths(overlay, roi, cleaned=True, filter_zeros=True, vertical_flip=True, overlay_file=None):
"""returns left, right hemisphere paths"""
if roi in list(overlay.rois.shapes.keys()):
paths = overlay.rois[roi].splines
elif roi in list(overlay.sulci.shapes.keys()):
Expand Down Expand Up @@ -670,7 +674,7 @@ def show_dartboard(data,
"""
if max_radius is None:
max_radius = 1
data = np.array(data).astype(np.float)
data = np.array(data).astype(float)
if isinstance(data2, np.ndarray):
data = np.stack([data, data2], axis=0)
else:
Expand Down Expand Up @@ -748,6 +752,27 @@ def show_dartboard(data,
return axis


def clockwise_test(pts):
"""Test for clockwise ordering of points
Parameters
----------
pts : array
2D array (n_pts, xy) defining polygon points to be evaluated for clockwise rotation
Returns
-------
is_clockwise : bool
True for clockwise, False for counter-clockwise
"""
if not np.all(pts[0] == pts[-1]):
pts = np.vstack([pts, pts[0]])
sum = 0
for c0, c1 in zip(pts[:-1], pts[1:]):
sum += (c1[0] - c0[0]) * (c1[1] + c0[1])
return sum > 0


def interpolate_outlines(phi, rho, resolution=50):
new_phi = np.linspace(0,2*np.pi,resolution)
new_rho = np.interp(new_phi, phi, rho, period=2*np.pi)
Expand Down Expand Up @@ -873,40 +898,54 @@ def sort_roi_border(xy, angles, start_angle=0, start_index=None, max_deg_from_ze
if start_index is None:
start_index = np.argmin(
np.abs(circ_dist(angles, start_angle, degrees=False)))
vertex_order = [start_index]
ct = 0
for i_vert in range(len(xy)):
# Special case for first one: go clockwise
if i_vert == 0:
new_verts = np.argsort(dsts[start_index])[:5]
new_verts = np.array([x for x in new_verts if x not in vertex_order])
angles_from_origin = circ_dist(angles[new_verts], angles[start_index], degrees=False)
is_clockwise = angles_from_origin > 0
#if ~any(is_clockwise):
# plt.plot(*xy.T)
# plt.plot()
min_clockwise_angle = np.min(angles_from_origin[is_clockwise])
j = new_verts[angles_from_origin == min_clockwise_angle][0]
vertex_order.append(j)
else:
success = False
n = 3
while (not success) and (n<15):
new_verts = np.argsort(dsts[vertex_order[-1]])[:n]
new_verts = np.array([x for x in new_verts if x not in vertex_order])
success = len(new_verts==1)
n += 1
#if n > 4:
# print(n)
if len(new_verts) > 0:
vertex_order.append(new_verts[-1])
else:
ct += 1
if verbose:
print(f"Skipped a vertex ({ct} skipped)")
#1/0
vertex_order.append(vertex_order[0])
return np.array(vertex_order)
is_clockwise = clockwise_test(xy)
n = len(xy)
if is_clockwise:
print('Detected CLOCKWISE order')
vertex_order = np.hstack([np.arange(start_index, n),
np.arange(0, start_index),
start_index])
else:
print('Detected COUNTER-CLOCKWISE order')
vertex_order = np.hstack([np.arange(start_index, 0, -1),
np.arange(n-1, start_index, -1),
start_index])

# vertex_order = [start_index]
# ct = 0
# for i_vert in range(len(xy)):
# # Special case for first one: go clockwise
# if i_vert == 0:
# new_verts = np.argsort(dsts[start_index])[:5]
# new_verts = np.array([x for x in new_verts if x not in vertex_order])
# angles_from_origin = circ_dist(angles[new_verts], angles[start_index], degrees=False)
# is_clockwise = angles_from_origin > 0
# #if ~any(is_clockwise):
# # plt.plot(*xy.T)
# # plt.plot()
# min_clockwise_angle = np.min(angles_from_origin[is_clockwise])
# j = new_verts[angles_from_origin == min_clockwise_angle][0]
# vertex_order.append(j)
# else:
# success = False
# n = 3
# while (not success) and (n<15):
# new_verts = np.argsort(dsts[vertex_order[-1]])[:n]
# new_verts = np.array([x for x in new_verts if x not in vertex_order])
# success = len(new_verts==1)
# n += 1
# #if n > 4:
# # print(n)
# if len(new_verts) > 0:
# vertex_order.append(new_verts[-1])
# else:
# ct += 1
# if verbose:
# print(f"Skipped a vertex ({ct} skipped)")
# #1/0
# vertex_order.append(vertex_order[0])
# return np.array(vertex_order)
return vertex_order


def get_interpolated_outlines(overlay,
Expand All @@ -926,15 +965,18 @@ def get_interpolated_outlines(overlay,
svgoverlay for a given subject
outline_roi : str
name of ROI to plot (must sexist in overlays.svg)
center_roi : str
name of center ROI for dartboard
anchor_angles : array
array of angles ??
geodesic_distances : bool, optional
Flag for whether to compute geodesic distances for eccentricity bins,
by default True
resolution : int, optional
number of points estimated for outline, by default 100
every_n : int, optional
smooth path by fitting a cubic spline to `every_n` vertices on border
recache : bool, optional
whether to force re-computation (recache=True) or allow loading from cached
file (recache=False)
verbose : bool, optional
verbose output or not
Returns
-------
Expand Down Expand Up @@ -1163,8 +1205,8 @@ def show_dartboard_pair(dartboard_data,
space, a little smoothing usually helps aesthetically. 1 is no smoothing
, by default 5
even_sampling_over : str, optional
How to resample ROI paths, by angle or along path length. 'angle' is perhaps slightly
more principled for convex ROIs, but 'path_length' gives better results for non-convex
How to resample ROI paths, by angle or along path length. 'polar angle' is slightly
more principled for convex ROIs, but 'path length' gives better results for non-convex
ROIs, by default 'path length'
roi_linewidth : int, optional
width of lines for drawn ROIs, by default 1
Expand Down Expand Up @@ -1566,10 +1608,17 @@ def _get_anchor_points(svg, center_roi, anchors, return_indices=True):
anchor_name, anchor_type = anchor, 'centroid'
if j == 0:
center_name, center_type = anchor_name, anchor_type
if anchor_type == 'nearest':
centroids[anchor_name] = _get_closest_vertex_to_roi(svg, anchor_name, center_name, return_indices=return_indices)
if 'nearest' in anchor_type:
if '_to_' in anchor_type:
_, _, roi_to = anchor_type.split('_')
centroids[f'{anchor_name}-{roi_to}'] = _get_closest_vertex_to_roi(svg, anchor_name, roi_to, return_indices=return_indices)
else:
centroids[anchor_name] = _get_closest_vertex_to_roi(svg, anchor_name, center_name, return_indices=return_indices)
elif anchor_type == 'centroid':
centroids[anchor_name] = get_roi_centroids(svg, anchor_name, return_indices=return_indices)
elif anchor_type == 'superior':
pass
#centroids[anchor_name] = _get_closest_vertex_to_roi(svg, anchor_name, center_name, return_indices=return_indices)
else:
raise ValueError("unknown anchor type specified: %s\n(Must be 'nearest' or 'centroid')"%(anchor_type))
return centroids
Expand Down

0 comments on commit 8c28079

Please sign in to comment.