Changeset 1670

Timestamp:

Aug 2, 2005, 2:27:25 PM (20 years ago)

Author:

ole

Message:

Refactored Interpolation_function and tested

Location:

inundation/ga/storm_surge/pyvolution

Files:

• least_squares.py (modified) (3 diffs)
• mesh.py (modified) (2 diffs)
• test_least_squares.py (modified) (1 diff)
• test_mesh.py (modified) (1 diff)
• test_util.py (modified) (3 diffs)
• util.py (modified) (2 diffs)
• wiki/issues.txt (modified) (1 diff)

inundation/ga/storm_surge/pyvolution/least_squares.py

@@ -412 +412 @@
         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)
@@ -463 +466 @@
             indices = inside_polygon(point_coordinates, P, verbose = verbose)
 
-            if verbose:
-                print 'Done'
-                if len(indices) != point_coordinates.shape[0]:
-                    print '%d points outside mesh have been cropped.'\
+
+            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
@@ -870 +875 @@
 
 
+
+
+class Interpolation_function:
+    """Interpolation_function - creates callable object f(t, id) or f(t,x,y)
+    which is interpolated from time series defined at vertices of
+    triangular mesh (such as those stored in sww files)
+
+    Let m be the number of vertices, n the number of triangles
+    and p the number of timesteps. 
+
+    Mandatory input
+        time:               px1 array of monotonously increasing times (Float)
+        quantities:         Dictionary of pxm arrays or 1 pxm array (Float)
+        
+    Optional input:
+        quantity_names:     List of keys into quantities dictionary 
+        vertex_coordinates: mx2 array of coordinates (Float)
+        triangles:          nx3 array of indices into vertex_coordinates (Int)
+        interpolation_points: array of coordinates to be interpolated to 
+        verbose:            Level of reporting
+    
+    
+    The quantities returned by the callable object are specified by
+    the list quantities which must contain the names of the
+    quantities to be returned and also reflect the order, e.g. for
+    the shallow water wave equation, on would have
+    quantities = ['stage', 'xmomentum', 'ymomentum']
+
+    The parameter interpolation_points decides at which points interpolated
+    quantities are to be computed whenever object is called.
+    If None, return average value
+    """
+
+    
+    
+    def __init__(self,
+                 time,
+                 quantities,
+                 quantity_names = None,  
+                 vertex_coordinates = None,
+                 triangles = None,
+                 interpolation_points = None,
+                 verbose = False):
+        """Initialise object and build spatial interpolation if required
+        """
+
+        from Numeric import array, zeros, Float, alltrue, concatenate,\
+             reshape, ArrayType
+
+        from util import mean, ensure_numeric
+        from config import time_format
+        import types
+
+
+
+        #Check temporal info
+        time = ensure_numeric(time)        
+        msg = 'Time must be a monotonuosly '
+        msg += 'increasing sequence %s' %time
+        assert alltrue(time[1:] - time[:-1] > 0 ), msg
+
+
+        #Check if quantities is a single array only
+        if type(quantities) != types.DictType:
+            quantities = ensure_numeric(quantities)
+            quantity_names = ['Attribute']
+
+            #Make it a dictionary
+            quantities = {quantity_names[0]: quantities}
+
+
+        #Use keys if no names are specified
+        if quantity_names is not None:
+            self.quantity_names = quantity_names
+        else:
+            self.quantity_names = quantities.keys()
+
+
+        #Check spatial info
+        if vertex_coordinates is None:
+            self.spatial = False
+        else:    
+            vertex_coordinates = ensure_numeric(vertex_coordinates)
+
+            assert triangles is not None, 'Triangles array must be specified'
+            triangles = ensure_numeric(triangles)
+            self.spatial = True            
+            
+  
+        #      
+        self.interpolation_points = interpolation_points #FIXWME Needed?
+        self.T = time[:]  #Time assumed to be relative to starttime
+        self.index = 0    #Initial time index
+        self.precomputed_values = {}
+            
+
+
+        #Precomputed spatial interpolation if requested
+        if interpolation_points is not None:
+            if self.spatial is False:
+                raise 'Triangles and vertex_coordinates must be specified'
+            
+
+            try:
+                interpolation_points = ensure_numeric(interpolation_points)
+            except:
+                msg = 'Interpolation points must be an N x 2 Numeric array '+\
+                      'or a list of points\n'
+                msg += 'I got: %s.' %( str(interpolation_points)[:60] + '...')
+                raise msg
+
+
+            for name in quantity_names:
+                self.precomputed_values[name] =\
+                                              zeros((len(self.T),
+                                                     len(interpolation_points)),
+                                                    Float)
+
+            #Build interpolator
+            interpol = Interpolation(vertex_coordinates,
+                                     triangles,
+                                     point_coordinates = interpolation_points,
+                                     alpha = 0,
+                                     precrop = False, 
+                                     verbose = verbose)
+
+            #if interpol.cropped_points is True:
+            #    raise 'Some interpolation points were outside mesh'
+            #FIXME: This will be raised if triangles are listed as
+            #discontinuous even though there is no need to stop
+            #(precrop = True above)
+
+            if verbose: print 'Interpolate'
+            for i, t in enumerate(self.T):
+                #Interpolate quantities at this timestep
+                if verbose: print ' time step %d of %d' %(i, len(self.T))
+                for name in quantity_names:
+                    self.precomputed_values[name][i, :] =\
+                    interpol.interpolate(quantities[name][i,:])
+
+            #Report
+            if verbose:
+                x = vertex_coordinates[:,0]
+                y = vertex_coordinates[:,1]                
+            
+                print '------------------------------------------------'
+                print 'Interpolation_function statistics:'
+                print '  Extent:'
+                print '    x in [%f, %f], len(x) == %d'\
+                      %(min(x), max(x), len(x))
+                print '    y in [%f, %f], len(y) == %d'\
+                      %(min(y), max(y), len(y))
+                print '    t in [%f, %f], len(t) == %d'\
+                      %(min(self.T), max(self.T), len(self.T))
+                print '  Quantities:'
+                for name in quantity_names:
+                    q = quantities[name][:].flat
+                    print '    %s in [%f, %f]' %(name, min(q), max(q))
+                print '  Interpolation points (xi, eta):'\
+                      ' number of points == %d ' %interpolation_points.shape[0]
+                print '    xi in [%f, %f]' %(min(interpolation_points[:,0]),
+                                             max(interpolation_points[:,0]))
+                print '    eta in [%f, %f]' %(min(interpolation_points[:,1]),
+                                              max(interpolation_points[:,1]))
+                print '  Interpolated quantities (over all timesteps):'
+                
+                for name in quantity_names:
+                    q = self.precomputed_values[name][:].flat
+                    print '    %s at interpolation points in [%f, %f]'\
+                          %(name, min(q), max(q))
+                print '------------------------------------------------'
+            
+        else:
+            #Store quantitites as is
+            for name in quantity_names:
+                self.precomputed_values[name] = quantities[name]
+
+
+        #else:
+        #    #Return an average, making this a time series
+        #    for name in quantity_names:
+        #        self.values[name] = zeros(len(self.T), Float)
+        #
+        #    if verbose: print 'Compute mean values'
+        #    for i, t in enumerate(self.T):
+        #        if verbose: print ' time step %d of %d' %(i, len(self.T))
+        #        for name in quantity_names:
+        #           self.values[name][i] = mean(quantities[name][i,:])
+
+
+
+
+    def __repr__(self):
+        return 'Interpolation function (spation-temporal)'
+
+    def __call__(self, t, point_id = None, x = None, y = None):
+        """Evaluate f(t), f(t, point_id) or f(t, x, y)
+
+        Inputs:
+          t: time - Model time. Must lie within existing timesteps
+          point_id: index of one of the preprocessed points.
+          x, y:     Overrides location, point_id ignored
+          
+          If spatial info is present and all of x,y,point_id
+          are None an exception is raised 
+                    
+          If no spatial info is present, point_id and x,y arguments are ignored
+          making f a function of time only.
+
+          
+          FIXME: point_id could also be a slice
+          FIXME: What if x and y are vectors?
+          FIXME: What about f(x,y) without t?
+        """
+
+        from math import pi, cos, sin, sqrt
+        from Numeric import zeros, Float
+        from util import mean        
+
+        if self.spatial is True:
+            if point_id is None:
+                if x is None or y is None:
+                    msg = 'Either point_id or x and y must be specified'
+                    raise msg
+            else:
+                if self.interpolation_points is None:
+                    msg = 'Interpolation_function must be instantiated ' +\
+                          'with a list of interpolation points before parameter ' +\
+                          'point_id can be used'
+                    raise msg
+
+
+        msg = 'Time interval [%s:%s]' %(self.T[0], self.T[1])
+        msg += ' does not match model time: %s\n' %t
+        if t < self.T[0]: raise msg
+        if t > self.T[-1]: raise msg
+
+        oldindex = self.index #Time index
+
+        #Find current time slot
+        while t > self.T[self.index]: self.index += 1
+        while t < self.T[self.index]: self.index -= 1
+
+        if t == self.T[self.index]:
+            #Protect against case where t == T[-1] (last time)
+            # - also works in general when t == T[i]
+            ratio = 0
+        else:
+            #t is now between index and index+1
+            ratio = (t - self.T[self.index])/\
+                    (self.T[self.index+1] - self.T[self.index])
+
+        #Compute interpolated values
+        q = zeros(len(self.quantity_names), Float)
+
+        for i, name in enumerate(self.quantity_names):
+            Q = self.precomputed_values[name]
+
+            if self.spatial is False:
+                #If there is no spatial info                
+                assert len(Q.shape) == 1
+
+                Q0 = Q[self.index]
+                if ratio > 0: Q1 = Q[self.index+1]
+
+            else:
+                if x is not None and y is not None:
+                    #Interpolate to x, y
+                    
+                    raise 'x,y interpolation not yet implemented'
+                else:
+                    #Use precomputed point
+                    Q0 = Q[self.index, point_id]
+                    if ratio > 0: Q1 = Q[self.index+1, point_id]
+
+            #Linear temporal interpolation    
+            if ratio > 0:
+                q[i] = Q0 + ratio*(Q1 - Q0)
+            else:
+                q[i] = Q0
+
+            #if ratio > 0:
+            #    q[i] = Q[self.index, point_id] +\
+            #           ratio*(Q[self.index+1, point_id] -\
+            #                  Q[self.index, point_id])
+            #else:
+            #    q[i] = Q[self.index, point_id]
+
+        #Return vector of interpolated values
+        if len(q) == 1:
+            return q[0]
+        else:
+            return q
+
+
+
+
+
 #-------------------------------------------------------------
 if __name__ == "__main__":

inundation/ga/storm_surge/pyvolution/mesh.py

@@ -416 +416 @@
     def get_boundary_polygon(self):
         """Return bounding polygon as a list of points
+
+        FIXME: If triangles are listed as discontinuous
+        (e.g vertex coordinates listed multiple times),
+        this may not work as expected.
         """
         from Numeric import allclose, sqrt, array, minimum, maximum
@@ -421 +425 @@
 
 
-        V = self.get_vertex_coordinates()
+        #V = self.get_vertex_coordinates()
         segments = {}
 

inundation/ga/storm_surge/pyvolution/test_least_squares.py

@@ -989 +989 @@
 
 
+    def test_interpolation_function_time_only(self):
+        """Test spatio-temporal interpolation
+        Test that spatio temporal function performs the correct
+        interpolations in both time and space
+        """
+
+
+        #Three timesteps
+        time = [1.0, 5.0, 6.0]
+        
+
+        #One quantity
+        Q = zeros( (3,6), Float )
+
+        #Linear in time and space
+        a = [0.0, 0.0]
+        b = [0.0, 2.0]
+        c = [2.0, 0.0]
+        d = [0.0, 4.0]
+        e = [2.0, 2.0]
+        f = [4.0, 0.0]
+
+        points = [a, b, c, d, e, f]
+        
+        for i, t in enumerate(time):
+            Q[i, :] = t*linear_function(points)
+
+            
+        #Check basic interpolation of one quantity using averaging
+        #(no interpolation points or spatial info)
+        from util import mean        
+        I = Interpolation_function(time, [mean(Q[0,:]),
+                                          mean(Q[1,:]),
+                                          mean(Q[2,:])])
+
+
+
+        #Check temporal interpolation
+        for i in [0,1,2]:
+            assert allclose(I(time[i]), mean(Q[i,:]))
+
+        #Midway    
+        assert allclose(I( (time[0] + time[1])/2 ),
+                        (I(time[0]) + I(time[1]))/2 )
+
+        assert allclose(I( (time[1] + time[2])/2 ),
+                        (I(time[1]) + I(time[2]))/2 )
+
+        assert allclose(I( (time[0] + time[2])/2 ),
+                        (I(time[0]) + I(time[2]))/2 )                 
+
+        #1/3
+        assert allclose(I( (time[0] + time[2])/3 ),
+                        (I(time[0]) + I(time[2]))/3 )                         
+
+
+        #Out of bounds checks
+        try:
+            I(time[0]-1) 
+        except:
+            pass
+        else:
+            raise 'Should raise exception'
+
+        try:
+            I(time[-1]+1) 
+        except:
+            pass
+        else:
+            raise 'Should raise exception'        
+
+
+        
+
+    def test_interpolation_function(self):
+        """Test spatio-temporal interpolation
+        Test that spatio temporal function performs the correct
+        interpolations in both time and space
+        """
+
+
+        #Three timesteps
+        time = [1.0, 5.0, 6.0]
+        
+
+        #Setup mesh used to represent fitted function
+        a = [0.0, 0.0]
+        b = [0.0, 2.0]
+        c = [2.0, 0.0]
+        d = [0.0, 4.0]
+        e = [2.0, 2.0]
+        f = [4.0, 0.0]
+
+        points = [a, b, c, d, e, f]
+        #bac, bce, ecf, dbe
+        triangles = [[1,0,2], [1,2,4], [4,2,5], [3,1,4]]
+
+
+        #New datapoints where interpolated values are sought
+        interpolation_points = [[ 0.0, 0.0],
+                                [ 0.5, 0.5],
+                                [ 0.7, 0.7],
+                                [ 1.0, 0.5],
+                                [ 2.0, 0.4],
+                                [ 2.8, 1.2]]
+
+
+        #One quantity
+        Q = zeros( (3,6), Float )
+
+        #Linear in time and space
+        for i, t in enumerate(time):
+            Q[i, :] = t*linear_function(points)
+
+
+        #Check interpolation of one quantity using interpolaton points)
+        I = Interpolation_function(time, Q,
+                                   vertex_coordinates = points,
+                                   triangles = triangles, 
+                                   interpolation_points = interpolation_points,
+                                   verbose = False)
+
+
+        answer = linear_function(interpolation_points)
+
+        t = time[0]
+        for j in range(50): #t in [1, 6]
+            for id in range(len(interpolation_points)):
+                assert allclose(I(t, id), t*answer[id])
+
+            t += 0.1    
+
+
+        try:    
+            I(1)
+        except:
+            pass
+        else:
+            raise 'Should raise exception'
+            
+        #
+        #interpolation_points = [[ 0.0, 0.0],
+        #                        [ 0.5, 0.5],
+        #                        [ 0.7, 0.7],
+        #                        [-13, 65],  #Outside
+        #                        [ 1.0, 0.5],
+        #                        [ 2.0, 0.4],
+        #                        [ 2.8, 1.2]]
+        #
+        #try:
+        #    I = Interpolation_function(time, Q,
+        #                               vertex_coordinates = points,
+        #                               triangles = triangles, 
+        #                               interpolation_points = interpolation_points,
+        #                               verbose = False)
+        #except:
+        #    pass
+        #else:
+        #    raise 'Should raise exception'
 
 

inundation/ga/storm_surge/pyvolution/test_mesh.py

@@ -800 +800 @@
 
 
+    def test_boundary_polygon_IV(self):
+        """Reproduce test test_spatio_temporal_file_function_time
+        from test_util.py that looked as if it produced the wrong boundary
+        """
+
+        from mesh import Mesh
+        from Numeric import zeros, Float
+        from util import inside_polygon
+        from mesh_factory import rectangular        
+
+        #Create a domain to hold test grid
+        #(0:15, -20:10)
+        points, vertices, boundary =\
+                rectangular(4, 4, 15, 30, origin = (0, -20))        
+
+        #####
+        mesh = Mesh(points, vertices)
+        mesh.check_integrity()
+
+        P = mesh.get_boundary_polygon()
+
+        #print P
+        assert len(P) == 16
+        for p in points:
+            assert inside_polygon(p, P)
+
+
+
+        #####
+        mesh = Mesh(points, vertices, boundary)
+        mesh.check_integrity()
+
+        P = mesh.get_boundary_polygon()
+
+        
+        #print P, len(P)
+        assert len(P) == 16
+
+        for p in points:
+            assert inside_polygon(p, P)            
+
+
+
+
 
 #-------------------------------------------------------------

inundation/ga/storm_surge/pyvolution/test_util.py

@@ -495 +495 @@
 
         #Create a domain to hold test grid
-
+        #(0:15, -20:10)
         points, vertices, boundary =\
                 rectangular(4, 4, 15, 30, origin = (0, -20))
@@ -566 +566 @@
                           quantities = domain.conserved_quantities,
                           interpolation_points = interpolation_points)                
-        #for id in range(len(interpolation_points)):
-        for id in [1]:
+        for id in range(len(interpolation_points)):
             x = interpolation_points[id][0]
             y = interpolation_points[id][1]
@@ -617 +616 @@
 
         #Now try interpolation with delta offset
-        for id in [1]:
+        for id in range(len(interpolation_points)):            
             x = interpolation_points[id][0]
             y = interpolation_points[id][1]

inundation/ga/storm_surge/pyvolution/util.py

@@ -291 +291 @@
     
 
-    return Interpolation_function(vertex_coordinates,
+    from least_squares import Interpolation_function
+    return Interpolation_function(time,
+                                  quantities,
+                                  quantity_names,
+                                  vertex_coordinates,
                                   triangles,
-                                  time,
-                                  quantities,
-                                  quantity_names,  
                                   interpolation_points,
-                                  alpha = 0,
                                   verbose = verbose)
 
@@ -303 +303 @@
 
 
-
-class Interpolation_function:
-    """Interpolation_function - creates callable object f(t, id) or f(t,x,y)
-    which is interpolated from time series defined at vertices of
-    triangular mesh 
-
-    Input
-        vertex_coordinates: mx2 array of coordinates (Float)
-        triangles:          nx3 array of indices into vertex_coordinates (Int)
-        time:               px1 array of times (Float)               
-        quantities:         Dictionary of pxm values
-        quantity_names:
-        interpolation_points:
-        alpha:
-        verbose:
-    
-
-    
-    The quantities returned by the callable object are specified by
-    the list quantities which must contain the names of the
-    quantities to be returned and also reflect the order, e.g. for
-    the shallow water wave equation, on would have
-    quantities = ['stage', 'xmomentum', 'ymomentum']
-
-    The parameter interpolation_points decides at which points interpolated
-    quantities are to be computed whenever object is called.
-    If None, return average value
-    """
-
-    
-
-    
-    def __init__(self,
-                 vertex_coordinates,
-                 triangles,
-                 time,
-                 quantities,
-                 quantity_names = None,  
-                 interpolation_points = None,
-                 alpha = 0,
-                 verbose = False):
-        """
-        """
-
-        from Numeric import array, zeros, Float, alltrue, concatenate,\
-             reshape, ArrayType
-
-        from config import time_format
-        from least_squares import Interpolation
-
-        #Check
-        msg = 'Time must be a monotonuosly '
-        msg += 'increasing sequence %s' %time
-        assert alltrue(time[1:] - time[:-1] > 0 ), msg
-
-
-        #Check if quantities is a single array only
-        if type(quantities) == ArrayType:
-            quantity_names = 'Attribute'
-            quantities = {quantity_names: quantities}            
-
-        #Use keys if no names specified
-        if quantity_names is not None:
-            self.quantity_names = quantity_names
-        else:
-            self.quantity_names = quantities.keys()
-
-
-        self.interpolation_points = interpolation_points
-        self.T = time[:]  #Time assumed to be relative to starttime
-        self.index = 0    #Initial time index
-        self.values = {}
-            
-
-        if interpolation_points is not None:
-            #Spatial interpolation
-
-            try:
-                interpolation_points = ensure_numeric(interpolation_points)
-            except:
-                msg = 'Interpolation points must be an N x 2 Numeric array or a list of points\n'
-                msg += 'I got: %s.' %( str(interpolation_points)[:60] + '...')
-                raise msg
-
-
-            for name in quantity_names:
-                self.values[name] = zeros( (len(self.T),
-                                            len(interpolation_points)),
-                                            Float)
-
-            #Build interpolator
-            interpol = Interpolation(vertex_coordinates,
-                                     triangles,
-                                     point_coordinates = interpolation_points,
-                                     alpha = 0,
-                                     verbose = verbose)
-
-
-            if verbose: print 'Interpolate'
-            for i, t in enumerate(self.T):
-                #Interpolate quantities at this timestep
-                if verbose: print ' time step %d of %d' %(i, len(self.T))
-                for name in quantity_names:
-                    self.values[name][i, :] =\
-                    interpol.interpolate(quantities[name][i,:])
-
-            #Report
-            if verbose:
-                x = vertex_coordinates[:,0]
-                y = vertex_coordinates[:,1]                
-                
-                print '------------------------------------------------'
-                print 'Interpolation_function statistics:'
-                print '  Extent:'
-                print '    x in [%f, %f], len(x) == %d'\
-                      %(min(x.flat), max(x.flat), len(x.flat))
-                print '    y in [%f, %f], len(y) == %d'\
-                      %(min(y.flat), max(y.flat), len(y.flat))
-                print '    t in [%f, %f], len(t) == %d'\
-                      %(min(self.T), max(self.T), len(self.T))
-                print '  Quantities:'
-                for name in quantity_names:
-                    q = quantities[name][:].flat
-                    print '    %s in [%f, %f]' %(name, min(q), max(q))
-
-
-                print '  Interpolation points (xi, eta):'\
-                      'number of points == %d ' %interpolation_points.shape[0]
-                print '    xi in [%f, %f]' %(min(interpolation_points[:,0]),
-                                             max(interpolation_points[:,0]))
-                print '    eta in [%f, %f]' %(min(interpolation_points[:,1]),
-                                              max(interpolation_points[:,1]))
-                print '  Interpolated quantities (over all timesteps):'
-                for name in quantity_names:
-                    q = self.values[name][:].flat
-                    print '    %s at interpolation points in [%f, %f]'\
-                          %(name, min(q), max(q))
-                print '------------------------------------------------'
-        else:
-            #Return an average, making this a time series
-            #FIXME: Needs a test
-            for name in quantity_names:
-                self.values[name] = avg(quantities[name][i,:])
-
-
-    def __repr__(self):
-        return 'Interpolation function'
-
-    def __call__(self, t, point_id = None, x = None, y = None):
-        """Evaluate f(t), f(t, point_id) or f(t, x, y)
-
-        Inputs:
-          t: time - Model time. Must lie within existing timesteps
-          point_id: index of one of the preprocessed points.
-                    If point_id is None all preprocessed points are computed
-                    
-          If no x,y info is present, point_id and x,y arguments are ignored
-          making f a function of time only.
-          
-          FIXME: point_id could also be a slice
-          FIXME: One could allow arbitrary x, y coordinates
-          (requires computation of a new interpolator)
-          FIXME: What if x and y are vectors?
-          FIXME: Does point_id or x,y take precedence?
-          FIXME: What about f(x,y) without t
-          
-        """
-
-        from math import pi, cos, sin, sqrt
-        from Numeric import zeros, Float
-
-
-        msg = 'Time interval [%s:%s]' %(self.T[0], self.T[1])
-        msg += ' does not match model time: %s\n' %t
-        if t < self.T[0]: raise msg
-        if t > self.T[-1]: raise msg
-
-        oldindex = self.index #Time index
-
-        #Find current time slot
-        while t > self.T[self.index]: self.index += 1
-        while t < self.T[self.index]: self.index -= 1
-
-        if t == self.T[self.index]:
-            #Protect against case where t == T[-1] (last time)
-            # - also works in general when t == T[i]
-            ratio = 0
-        else:
-            #t is now between index and index+1
-            ratio = (t - self.T[self.index])/\
-                    (self.T[self.index+1] - self.T[self.index])
-
-        #Compute interpolated values
-        q = zeros(len(self.quantity_names), Float)
-
-        for i, name in enumerate(self.quantity_names):
-            Q = self.values[name]
-
-            if point_id is None:
-                Q0 = Q[self.index]
-
-                try:
-                    Q1 = Q[self.index+1]
-                except:
-                    pass                
-            else:    
-                Q0 = Q[self.index, point_id]
-
-                try:
-                    Q1 = Q[self.index+1, point_id]
-                except:
-                    pass
-                
-
-            #Linear temporal interpolation    
-            if ratio > 0:
-                q[i] = Q0 + ratio*(Q1 - Q0)
-            else:
-                q[i] = Q0
-            
-
-            #if ratio > 0:
-            #    q[i] = Q[self.index, point_id] +\
-            #           ratio*(Q[self.index+1, point_id] -\
-            #                  Q[self.index, point_id])
-            #else:
-            #    q[i] = Q[self.index, point_id]
-
-        #Return vector of interpolated values
-        return q
 
 

inundation/ga/storm_surge/pyvolution/wiki/issues.txt

@@ -8 +8 @@
  are outside the mesh.
 Importance: Mid
+Comment (Ole): If points and triangles are listed as discontinuous e.g. vertex coordinates replicated for each triangle, then the bounding polygon gets mangled and some points appear to be outside.
 Who: DSG