class BoundaryRegion:
    def __init__(self, points, xRange = lambda x: x, yRange = lambda y: y):
        self.rawPoints = np.array(points) # must be of shape mxn where m is the number of points and n is the number of dimensions
        self.affinePoints = np.array([[xRange(x), yRange(y)] for x, y in zip(points[:, 0], points[:, 1])]) 
        self.sort_points()
        self.polygon = Polygon(self.points)
        self.tck = None
        self.u = None
        
    def sort_points(self):
        """ Sort points to form a valid polygon using nearest neighbor heuristic. """
        if len(self.rawPoints) < 3:
            return  # Not enough points to form a polygon

        sorted_points = [self.rawPoints[0]]
        
        # Iterate to find the nearest unvisited point
        current_point = sorted_points[0]
        visitedIndexs = [0]
        while True:
            dx = (self.rawPoints[:, 0] - current_point[0])
            dy = (self.rawPoints[:, 1] - current_point[1])
            distances = np.sqrt(dx**2 + dy**2)
            distances[visitedIndexs] = np.inf
            nearest_index = np.argmin(distances)
            if nearest_index not in visitedIndexs:
                visitedIndexs.append(nearest_index)
                current_point = self.rawPoints[nearest_index]
                sorted_points.append(current_point)
            else:
                break
        
        self.rawPoints = np.array(sorted_points)
        self.affinePoints = self.affinePoints[visitedIndexs]

    @property
    def points(self):
        return self.affinePoints

    def is_point_inside(self, test_point):
        """ Check if a point is inside the polygon defined by the boundary points using shapely. """
        point = Point(test_point)
        return self.polygon.contains(point)

    def parameterize_boundary(self):
        """ Parameterize the boundary using spline interpolation. """
        t = np.linspace(0, 1, len(self.points))
        x, y = self.points[:, 0], self.points[:, 1]
        self.tck, self.u = interp1d(t, x, kind='cubic'), interp1d(t, y, kind='cubic')
        return self.tck, self.u

    def get_point_on_boundary(self, t):
        """ Get a point on the parameterized boundary. """
        if self.tck is None or self.u is None:
            self.parameterize_boundary()
        x = self.tck(t)
        y = self.u(t)
        return x, y

    def plot_region(self, ax):
        """ Plot the region on a given matplotlib axis. """
        from matplotlib.patches import Polygon as MplPolygon
        polygon = MplPolygon(self.points, closed=True, edgecolor='blue', facecolor='none')
        ax.add_patch(polygon)