Changeset 2750 for inundation/fit_interpolate/least_squares.py

Timestamp:

Apr 21, 2006, 5:20:45 PM (19 years ago)

Author:

duncan

Message:

Do interpolation using fit_interpolate/interpolate.py, instead of least_squares.py.

File:

• inundation/fit_interpolate/least_squares.py (modified) (2 diffs)

inundation/fit_interpolate/least_squares.py

@@ -534 +534 @@
         """
 
-
-
-        #FIXME (Ole): Check that this function is memeory efficient.
-        #6 million datapoints and 300000 basis functions
-        #causes out-of-memory situation
-        #First thing to check is whether there is room for self.A and self.AtA
-        #
-        #Maybe we need some sort of blocking
-
         from pyvolution.quad import build_quadtree
         from utilities.polygon import inside_polygon
         
-
-        if data_origin is None:
-            data_origin = self.data_origin #Use the one from
-                                           #interpolation instance
-
         #Convert input to Numeric arrays just in case.
         point_coordinates = ensure_numeric(point_coordinates, Float)
-
-        #Keep track of discarded points (if any).
-        #This is only registered if precrop is True
-        self.cropped_points = False
-
-        #Shift data points to same origin as mesh (if specified)
-
-        #FIXME this will shift if there was no geo_ref.
-        #But all this should be removed anyhow.
-        #change coords before this point
-        mesh_origin = self.mesh.geo_reference.get_origin()
-        if point_coordinates is not None:
-            if data_origin is not None:
-                if mesh_origin is not None:
-
-                    #Transformation:
-                    #
-                    #Let x_0 be the reference point of the point coordinates
-                    #and xi_0 the reference point of the mesh.
-                    #
-                    #A point coordinate (x + x_0) is then made relative
-                    #to xi_0 by
-                    #
-                    # x_new = x + x_0 - xi_0
-                    #
-                    #and similarly for eta
-
-                    x_offset = data_origin[1] - mesh_origin[1]
-                    y_offset = data_origin[2] - mesh_origin[2]
-                else: #Shift back to a zero origin
-                    x_offset = data_origin[1]
-                    y_offset = data_origin[2]
-
-                point_coordinates[:,0] += x_offset
-                point_coordinates[:,1] += y_offset
-            else:
-                if mesh_origin is not None:
-                    #Use mesh origin for data points
-                    point_coordinates[:,0] -= mesh_origin[1]
-                    point_coordinates[:,1] -= mesh_origin[2]
-
-
-
-        #Remove points falling outside mesh boundary
-        #This reduced one example from 1356 seconds to 825 seconds
-
-        
-        if precrop is True:
-            from Numeric import take
-
-            if verbose: print 'Getting boundary polygon'
-            P = self.mesh.get_boundary_polygon()
-
-            if verbose: print 'Getting indices inside mesh boundary'
-            indices = inside_polygon(point_coordinates, P, verbose = verbose)
-
-
-            if len(indices) != point_coordinates.shape[0]:
-                self.cropped_points = True
-                if verbose:
-                    print 'Done - %d points outside mesh have been cropped.'\
-                          %(point_coordinates.shape[0] - len(indices))
-
-            point_coordinates = take(point_coordinates, indices)
-            self.point_indices = indices
-
-
-
 
         #Build n x m interpolation matrix
@@ -629 +547 @@
         if verbose: print 'Number of basis functions: %d' %m
 
-        #FIXME (Ole): We should use CSR here since mat-mat mult is now OK.
-        #However, Sparse_CSR does not have the same methods as Sparse yet
-        #The tests will reveal what needs to be done
-
-        #
-        #self.A = Sparse_CSR(Sparse(n,m))
-        #self.AtA = Sparse_CSR(Sparse(m,m))
         self.A = Sparse(n,m)
         self.AtA = Sparse(m,m)
 
-        #Build quad tree of vertices (FIXME: Is this the right spot for that?)
+        #Build quad tree of vertices
         root = build_quadtree(self.mesh,
                               max_points_per_cell = max_points_per_cell)