"""Alpha shape """ import exceptions from Numeric import array, Float from load_mesh.loadASCII import load_xya_file, export_boundary_file class PointError(exceptions.Exception): pass OUTPUT_FILE_TITLE = "# The alpha shape boundary defined by point index pairs of edges" def alpha_shape_via_files(point_file, boundary_file, alpha= None): from load_mesh.loadASCII import load_xya_file point_dict = load_xya_file(point_file) points = point_dict['pointlist'] #title_string = point_dict['title'] alpha = Alpha_Shape(points) alpha.write_boundary(boundary_file) class Alpha_Shape: def __init__(self, points, alpha = None): """ Build interpolation matrix mapping from function values at vertices to function values at data points Inputs: points: List of coordinate pairs [[x1, y1],[x2, y2]..] alpha: alpha shape parameter """ self._set_points(points) self._alpha_shape_algorithm() def _set_points(self, points): """ """ if len (points) <= 2: raise PointError, "Too few points to find an alpha shape" #Convert input to Numeric arrays self.points = array(points).astype(Float) if len (points) <= 2: raise PointError, "Too few points to find an alpha shape" def write_boundary(self,file_name): """ Write the boundary to a file """ #print " this info will be in the file",boundary export_boundary_file(file_name, self.get_boundary(), OUTPUT_FILE_TITLE, delimiter = ',') def get_boundary(self): """ """ return self.boundary def _alpha_shape_algorithm(self): # A python style guide suggests using _[function name] to # specify interal functions #- this is the first time I've used it though - DSG # this produces a baaad boundary boundary = [] for point_index in range(len(self.points)-1): boundary.append([point_index,point_index +1]) boundary.append([len(self.points)-1,0]) self.boundary = boundary #------------------------------------------------------------- if __name__ == "__main__": """ Load in a data point file. Determine the alpha shape boundary Save the boundary to a file. usage: alpha_shape.py point_file.xya boundary_file.bnd [alpha] The alpha value is optional. """ import os, sys usage = "usage: %s point_file.xya boundary_file.bnd [alpha]"\ %os.path.basename(sys.argv[0]) # I made up the .bnd affix. Other ideas welcome. -DSG if len(sys.argv) < 3: print usage else: point_file = sys.argv[1] boundary_file = sys.argv[2] if len(sys.argv) > 4: alpha = sys.argv[3] else: alpha = None alpha_shape_via_files(point_file, boundary_file, alpha)