Skip to content

Commit

Permalink
Synchronize gwcs region with changes in jwst/romancal
Browse files Browse the repository at this point in the history
  • Loading branch information
WilliamJamieson committed Oct 18, 2024
1 parent 9006b03 commit 23200c5
Showing 1 changed file with 67 additions and 47 deletions.
114 changes: 67 additions & 47 deletions gwcs/region.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,7 @@
__all__ = ['Region', 'Edge', 'Polygon']


class Region(metaclass=abc.ABCMeta):

class Region:
"""
Base class for regions.
Expand All @@ -33,14 +32,14 @@ def __init__(self, rid, coordinate_frame):
self._rid = rid

@abc.abstractmethod
def __contains__(self, x, y):
def __contains__(self, px):
"""
Determines if a pixel is within a region.
Parameters
----------
x, y : float
x , y values of a pixel
px : tuple[float, float]
A pixel coordinate (x, y)
Returns
-------
Expand All @@ -49,6 +48,7 @@ def __contains__(self, x, y):
Subclasses must define this method.
"""

@abc.abstractmethod
def scan(self, mask):
"""
Sets mask values to region id for all pixels within the region.
Expand All @@ -68,7 +68,6 @@ def scan(self, mask):


class Polygon(Region):

"""
Represents a 2D polygon region with multiple vertices.
Expand All @@ -86,11 +85,11 @@ class Polygon(Region):
"""

def __init__(self, rid, vertices, coord_frame="detector"):
def __init__(self, rid, vertices, coord_system="Cartesian"):
if len(vertices) < 4:
raise ValueError("Expected vertices to be "
"a list of minimum 4 tuples (x,y)")
super(Polygon, self).__init__(rid, coord_frame)
raise ValueError("Expected vertices to be a list of minimum 4 tuples (x,y)")

Check warning on line 90 in gwcs/region.py

View check run for this annotation

Codecov / codecov/patch

gwcs/region.py#L90

Added line #L90 was not covered by tests

super().__init__(rid, coord_system)

# self._shiftx & self._shifty are introduced to shift the bottom-left
# corner of the polygon's bounding box to (0,0) as a (hopefully
Expand All @@ -114,8 +113,9 @@ def __init__(self, rid, vertices, coord_frame="detector"):
self._shifty = int(round(self._shifty))

self._bbox = self._get_bounding_box()
self._scan_line_range = \
list(range(self._bbox[1], self._bbox[3] + self._bbox[1] + 1))
self._scan_line_range = list(
range(self._bbox[1], self._bbox[3] + self._bbox[1] + 1)
)
# constructs a Global Edge Table (GET) in bbox coordinates
self._GET = self._construct_ordered_GET()

Expand All @@ -130,7 +130,7 @@ def _construct_ordered_GET(self):
"""
Construct a Global Edge Table (GET)
The GET is an OrderedDict. Keys are scan line numbers,
The GET is an OrderedDict. Keys are scan line numbers,
ordered from bbox.ymin to bbox.ymax, where bbox is the
bounding box of the polygon.
Values are lists of edges for which edge.ymin==scan_line_number.
Expand All @@ -149,7 +149,7 @@ def _construct_ordered_GET(self):
ymin_ind = (ymin == i).nonzero()[0]
# a hack for incomplete filling .any() fails if 0 is in ymin_ind
# if ymin_ind.any():
yminindlen, = ymin_ind.shape
(yminindlen,) = ymin_ind.shape
if yminindlen:
GET[i] = [edges[ymin_ind[0]]]
for j in ymin_ind[1:]:
Expand All @@ -160,8 +160,10 @@ def get_edges(self):
"""
Create a list of Edge objects from vertices
"""
return [Edge(name=f'E{i - 1}', start=self._vertices[i - 1], stop=self._vertices[i])
for i in range(1, len(self._vertices))]
return [
Edge(name=f"E{i - 1}", start=self._vertices[i - 1], stop=self._vertices[i])
for i in range(1, len(self._vertices))
]

def scan(self, data):
"""
Expand All @@ -186,10 +188,11 @@ def scan(self, data):
5. Set elements between pairs of X in the AET to the Edge's ID
"""
# TODO: 1.This algorithm does not mark pixels in the top row and left most column.
# Pad the initial pixel description on top and left with 1 px to prevent this.
# 2. Currently it uses intersection of the scan line with edges. If this is
# too slow it should use the 1/m increment (replace 3 above) (or the increment
# TODO: 1.This algorithm does not mark pixels in the top
# row and left most column. Pad the initial pixel description on
# top and left with 1 px to prevent this. 2. Currently it uses
# intersection of the scan line with edges. If this is too slow it
# should use the 1/m increment (replace 3 above) (or the increment
# should be removed from the GET entry).

# see comments in the __init__ function for the reason of introducing
Expand All @@ -204,18 +207,23 @@ def scan(self, data):
scline = self._scan_line_range[-1]

while y <= scline:

if y < scline:
AET = self.update_AET(y, AET)

if self._bbox[2] <= 0:
y += 1
continue

scan_line = Edge('scan_line', start=[self._bbox[0], y],
stop=[self._bbox[0] + self._bbox[2], y])
x = [int(np.ceil(e.compute_AET_entry(scan_line)[1]))
for e in AET if e is not None]
scan_line = Edge(
"scan_line",
start=[self._bbox[0], y],
stop=[self._bbox[0] + self._bbox[2], y],
)
x = [
int(np.ceil(e.compute_AET_entry(scan_line)[1]))
for e in AET
if e is not None
]
xnew = np.sort(x)
ysh = y + self._shifty

Expand All @@ -226,7 +234,7 @@ def scan(self, data):
for i, j in zip(xnew[::2], xnew[1::2]):
xstart = max(0, i + self._shiftx)
xend = min(j + self._shiftx, nx - 1)
data[ysh][xstart:xend + 1] = self._rid
data[ysh][xstart : xend + 1] = self._rid

y += 1

Expand Down Expand Up @@ -254,13 +262,16 @@ def update_AET(self, y, AET):
return AET

def __contains__(self, px):
"""even-odd algorithm or smth else better sould be used"""
return px[0] >= self._bbox[0] and px[0] <= self._bbox[0] + self._bbox[2] and \
px[1] >= self._bbox[1] and px[1] <= self._bbox[1] + self._bbox[3]
"""even-odd algorithm or something else better should be used"""
return (

Check warning on line 266 in gwcs/region.py

View check run for this annotation

Codecov / codecov/patch

gwcs/region.py#L266

Added line #L266 was not covered by tests
px[0] >= self._bbox[0]
and px[0] <= self._bbox[0] + self._bbox[2]
and px[1] >= self._bbox[1]
and px[1] <= self._bbox[1] + self._bbox[3]
)


class Edge:

"""
Edge representation.
Expand All @@ -271,7 +282,7 @@ class Edge:
"""

def __init__(self, name=None, start=None, stop=None, next=None):
def __init__(self, name=None, start=None, stop=None, next=None): # noqa: A002
self._start = None
if start is not None:
self._start = np.asarray(start)
Expand Down Expand Up @@ -329,8 +340,12 @@ def compute_GET_entry(self):
if np.diff(earr[:, 1]).item() == 0:
return None
else:
entry = [self._ymax, self._yminx,
(np.diff(earr[:, 0]) / np.diff(earr[:, 1])).item(), None]
entry = [
self._ymax,
self._yminx,
(np.diff(earr[:, 0]) / np.diff(earr[:, 1])).item(),
None,
]
return entry

def compute_AET_entry(self, edge):
Expand All @@ -347,19 +362,20 @@ def __repr__(self):
fmt = ""
if self._name is not None:
fmt += self._name
next = self.next
while next is not None:
next_edge = self.next
while next_edge is not None:

Check warning on line 366 in gwcs/region.py

View check run for this annotation

Codecov / codecov/patch

gwcs/region.py#L365-L366

Added lines #L365 - L366 were not covered by tests
fmt += "-->"
fmt += next._name
next = next.next
fmt += next_edge._name
next_edge = next_edge.next

Check warning on line 369 in gwcs/region.py

View check run for this annotation

Codecov / codecov/patch

gwcs/region.py#L368-L369

Added lines #L368 - L369 were not covered by tests

return fmt

@property
def next(self):
def next(self): # noqa: A003
return self._next

@next.setter
def next(self, edge):
def next(self, edge): # noqa: A003
if self._name is None:
self._name = edge._name
self._stop = edge._stop
Expand All @@ -372,21 +388,25 @@ def intersection(self, edge):
u = self._stop - self._start
v = edge._stop - edge._start
w = self._start - edge._start
D = _cross(u, v)

if abs(D) <= 1e2 * np.finfo(float).eps:
# Find the determinant of the matrix formed by the vectors u and v
# Note: Originally this was computed using a numpy "2D" cross product,
# however, this functionality has been deprecated and slated for
# removal.
D = np.linalg.det([u, v])

if np.allclose(D, 0, rtol=0, atol=1e2 * np.finfo(float).eps):
return np.array(self._start)

return _cross(v, w) / D * u + self._start
# See note above
return np.linalg.det([v, w]) / D * u + self._start

def is_parallel(self, edge):
u = self._stop - self._start
v = edge._stop - edge._start
return abs(_cross(u, v)) <= 1e2 * np.finfo(float).eps


def _cross(u, v):
return u[0] * v[1] - u[1] * v[0]
return np.allclose(

Check warning on line 407 in gwcs/region.py

View check run for this annotation

Codecov / codecov/patch

gwcs/region.py#L407

Added line #L407 was not covered by tests
np.linalg.det([u, v]), 0, rtol=0, atol=1e2 * np.finfo(float).eps
)


def _round_vertex(v):
Expand Down

0 comments on commit 23200c5

Please sign in to comment.