Changeset 7685


Ignore:
Timestamp:
Apr 14, 2010, 6:53:12 PM (10 years ago)
Author:
hudson
Message:

Ticket 344 completed and unit tests pass.
Stripped out a heap of tabs and replaced them with 4 spaces.
Removed some redundant imports and deprecated functions.

Location:
anuga_core/source/anuga
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/abstract_2d_finite_volumes/gauge.py

    r7679 r7685  
    9292    import string
    9393    from anuga.shallow_water.data_manager import get_all_swwfiles
    94     from anuga.abstract_2d_finite_volumes.util import file_function     
     94    from anuga.abstract_2d_finite_volumes.util import file_function   
    9595
    9696    assert type(gauge_file) == type(''), 'Gauge filename must be a string'
     
    246246                        if quantity == 'ycentroid':
    247247                            points_list.append(callable_sww.centroids[point_i][1])
    248                                                        
     248                           
    249249                points_writer[point_i].writerow(points_list)
    250250
     
    281281                   use_cache=False,
    282282                   verbose=False,
    283                                    output_centroids=False):
     283                   output_centroids=False):
    284284    """ Read sww file and plot the time series for the
    285285    prescribed quantities at defined gauge locations and
     
    399399                        use_cache,
    400400                        verbose,
    401                                                 output_centroids = output_centroids)
     401                        output_centroids = output_centroids)
    402402    return k
    403403
     
    435435                    use_cache = False,
    436436                    verbose = False,
    437                                         output_centroids = False):   
     437                    output_centroids = False):   
    438438       
    439439    # FIXME(Ole): Shouldn't print statements here be governed by verbose?
     
    499499                          verbose = verbose,
    500500                          use_cache = use_cache,
    501                                                   output_centroids = output_centroids)
     501                          output_centroids = output_centroids)
    502502
    503503        # determine which gauges are contained in sww file
  • anuga_core/source/anuga/abstract_2d_finite_volumes/test_gauge.py

    r7679 r7685  
    389389        os.remove(point1_filename)
    390390        os.remove(point2_filename)
    391                
     391       
    392392
    393393    def test_sww2csv_centroid(self):
     
    399399       
    400400        Test the ability to get a timeseries at the centroid of a triangle, rather
    401                 than the given gauge point.
     401        than the given gauge point.
    402402        """
    403403       
     
    453453        points_file = tempfile.mktemp(".csv")
    454454        file_id = open(points_file,"w")
    455 # These values are where the centroids should be               
     455# These values are where the centroids should be       
    456456#        file_id.write("name, easting, northing, elevation \n\
    457457#point1, 2.0, 2.0, 3.0\n\
     
    501501        point1_handle.close()
    502502        point2_handle.close()
    503         os.remove(mesh_file)           
     503        os.remove(mesh_file)       
    504504        os.remove(sww.filename)
    505505        os.remove(points_file)
     
    516516       
    517517        Test the ability to get a timeseries at the centroid of a triangle, rather
    518                 than the given gauge point, then output the results.
     518        than the given gauge point, then output the results.
    519519        """
    520520       
     
    596596        # clean up
    597597        point1_handle.close()
    598         os.remove(mesh_file)           
     598        os.remove(mesh_file)       
    599599        os.remove(sww.filename)
    600600        os.remove(points_file)
    601601        os.remove(point1_filename)
    602                
    603                
     602       
     603       
    604604
    605605#-------------------------------------------------------------
  • anuga_core/source/anuga/abstract_2d_finite_volumes/util.py

    r7679 r7685  
    158158    """Temporary Interface to new location"""
    159159
    160     import anuga.utilities.numerical_tools as NT       
     160    import anuga.utilities.numerical_tools as NT   
    161161   
    162162    msg = 'angle has moved from util.py.  '
     
    187187    warn(msg, DeprecationWarning)
    188188
    189     return utilities.polygon.point_on_line(*args, **kwargs)     
     189    return utilities.polygon.point_on_line(*args, **kwargs)   
    190190
    191191   
     
    211211    return utilities.polygon.outside_polygon(*args, **kwargs)   
    212212
    213 
    214 # @note TEMP
    215 def separate_points_by_polygon(*args, **kwargs):
    216     """Temporary Interface to new location"""
    217 
    218     log.critical('separate_points_by_polygon has moved from util.py.')
    219     log.critical('Please use "from anuga.utilities.polygon import '
    220                  'separate_points_by_polygon"')
    221 
    222     return utilities.polygon.separate_points_by_polygon(*args, **kwargs)   
    223 
    224 
     213   
    225214# @note TEMP
    226215def read_polygon(*args, **kwargs):
     
    947936        # verts=num.take(verts,X,axis=0)  to Remove the loners from verts
    948937        # but I think it would use more memory
    949         new_i = lone_start      # point at first loner - 'shuffle down' target
     938        new_i = lone_start    # point at first loner - 'shuffle down' target
    950939        for i in range(lone_start, N):
    951             if loners[i] >= N:  # [i] is a loner, leave alone
     940            if loners[i] >= N:    # [i] is a loner, leave alone
    952941                pass
    953             else:               # a non-loner, move down
     942            else:        # a non-loner, move down
    954943                loners[i] = new_i
    955944                verts[new_i] = verts[i]
     
    15751564    """
    15761565    from gauge import sww2csv_gauges as sww2csv
    1577        
     1566   
    15781567    return sww2csv(sww_file, gauge_file, out_name, quantities, verbose, use_cache)
    1579        
     1568   
    15801569def sww2timeseries(swwfiles,
    15811570                   gauge_filename,
     
    16101599                   use_cache,
    16111600                   verbose,
    1612                    output_centroids)                               
    1613        
     1601                   output_centroids)                   
     1602   
    16141603##
    16151604# @brief Get a wave height at a certain depth given wave height at another depth.
  • anuga_core/source/anuga/fit_interpolate/general_fit_interpolate.py

    r7317 r7685  
    2727from anuga.utilities.cg_solve import conjugate_gradient, VectorShapeError
    2828from anuga.utilities.numerical_tools import ensure_numeric
    29 from anuga.utilities.polygon import in_and_outside_polygon
    3029from anuga.utilities.quad import build_quadtree
    3130
  • anuga_core/source/anuga/fit_interpolate/interpolate.py

    r7682 r7685  
    155155                                 use_cache=use_cache,
    156156                                 verbose=verbose,
    157                                                                 output_centroids=output_centroids)
     157                                output_centroids=output_centroids)
    158158
    159159    return result
     
    230230                    start_blocking_len=500000,
    231231                    verbose=False,
    232                                         output_centroids=False):
     232                    output_centroids=False):
    233233        """Interpolate mesh data f to determine values, z, at points.
    234234
     
    324324    # @param use_cache True if caching should be used to accelerate the calculations
    325325    # @param verbose True if this function is verbose.
    326     # @return ??
     326    # @return interpolated f
    327327    def interpolate_block(self, f, point_coordinates,
    328328                          use_cache=False, verbose=False, output_centroids=False):
     
    407407    ##
    408408    # @brief Get interpolated data at given points.
    409         #        Applies a transform to all points to calculate the
    410         #        interpolated values. Points outside the mesh are returned as NaN.
    411         # @note self._A matrix must be valid
     409    #        Applies a transform to all points to calculate the
     410    #        interpolated values. Points outside the mesh are returned as NaN.
     411    # @note self._A matrix must be valid
    412412    # @param f Array of arbitrary data
    413413    # @param verbose True if this function is to be verbose.
     
    431431    # @brief Build NxM interpolation matrix.
    432432    # @param point_coordinates Points to sample at
    433         # @param output_centroids set to True to always sample from the centre
    434         #                         of the intersected triangle, instead of the intersection
    435         #                         point.
     433    # @param output_centroids set to True to always sample from the centre
     434    #                         of the intersected triangle, instead of the intersection
     435    #                         point.
    436436    # @param verbose True if this function is to be verbose.
    437437    # @return Interpolation matrix A, plus lists of the points inside and outside the mesh
    438         #         and the list of centroids, if requested.
     438    #         and the list of centroids, if requested.
    439439    def _build_interpolation_matrix_A(self,
    440440                                      point_coordinates,
     
    464464        if verbose: log.critical('Getting indices inside mesh boundary')
    465465
    466         inside_poly_indices, outside_poly_indices = \
     466        # Quick test against boundary, but will not deal with holes in the mesh
     467        inside_boundary_indices, outside_poly_indices = \
    467468            in_and_outside_polygon(point_coordinates,
    468469                                   self.mesh.get_boundary_polygon(),
     
    474475
    475476        # Since you can block, throw a warning, not an error.
    476         if verbose and 0 == len(inside_poly_indices):
     477        if verbose and 0 == len(inside_boundary_indices):
    477478            log.critical('WARNING: No points within the mesh!')
    478479
     
    485486        A = Sparse(n,m)
    486487
    487         n = len(inside_poly_indices)
     488        n = len(inside_boundary_indices)
    488489
    489490        centroids = []
    490                
     491       
     492        inside_poly_indices = []
     493       
    491494        # Compute matrix elements for points inside the mesh
    492495        if verbose: log.critical('Building interpolation matrix from %d points'
    493496                                 % n)
    494497
    495         for d, i in enumerate(inside_poly_indices):
     498        for d, i in enumerate(inside_boundary_indices):
    496499            # For each data_coordinate point
    497500            if verbose and d%((n+10)/10)==0: log.critical('Doing %d of %d'
     
    501504            element_found, sigma0, sigma1, sigma2, k = \
    502505                           search_tree_of_vertices(self.root, x)
     506                       
     507        # Update interpolation matrix A if necessary
     508            if element_found is True:
     509
     510                if verbose:
     511                    print 'Point is within mesh:', d, i           
    503512           
    504             # Update interpolation matrix A if necessary
    505             if element_found is True:
     513                inside_poly_indices.append(i)
     514               
    506515                # Assign values to matrix A
    507516                j0 = self.mesh.triangles[k,0] # Global vertex id for sigma0
     
    509518                j2 = self.mesh.triangles[k,2] # Global vertex id for sigma2
    510519                js = [j0, j1, j2]
    511                                
     520               
    512521                if output_centroids is False:
    513522                    # Weight each vertex according to its distance from x
     
    516525                        A[i, j] = sigmas[j]
    517526                else:
    518                                     # If centroids are needed, weight all 3 vertices equally
     527                    # If centroids are needed, weight all 3 vertices equally
    519528                    for j in js:
    520529                        A[i, j] = 1.0/3.0
    521                     centroids.append(self.mesh.centroid_coordinates[k])                                         
     530                    centroids.append(self.mesh.centroid_coordinates[k])                       
    522531            else:
    523                  # Mesh has a hole - move this point to the outside list
    524 #                 num.delete(inside_poly_indices,[i], axis=0)
    525                  num.append(outside_poly_indices, [i], axis=0)
    526 #                msg = 'Could not find triangle for point', x
    527 #                raise Exception(msg)
     532                if verbose:
     533                    log.critical('Mesh has a hole - moving this point to outside list')
     534
     535                # This is a numpy arrays, so we need to do a slow transfer
     536                outside_poly_indices = num.append(outside_poly_indices, [i], axis=0)
    528537
    529538        return A, inside_poly_indices, outside_poly_indices, centroids
     
    775784                 verbose=False,
    776785                 gauge_neighbour_id=None,
    777                                 output_centroids=False):
     786                output_centroids=False):
    778787        """Initialise object and build spatial interpolation if required
    779788
     
    10031012                                                      verbose=False,
    10041013                                                      output_centroids=output_centroids)
    1005                         self.centroids = interpol.centroids                                                                                                              
     1014                        self.centroids = interpol.centroids                                                         
    10061015                    elif triangles is None and vertex_coordinates is not None:
    10071016                        result = interpolate_polyline(Q,
     
    10121021
    10131022                    #assert len(result), len(interpolation_points)
    1014                     self.precomputed_values[name][i, :] = result                                                                       
    1015                                        
     1023                    self.precomputed_values[name][i, :] = result                                   
     1024                   
    10161025            # Report
    10171026            if verbose:
    1018                 log.critical(self.statistics())                 
     1027                log.critical(self.statistics())           
    10191028        else:
    10201029            # Store quantitites as is
  • anuga_core/source/anuga/fit_interpolate/test_interpolate.py

    r7675 r7685  
    107107        assert num.allclose(A.todense(), [[1./3, 1./3, 1./3]])
    108108
    109 
     109    def test_datapoint_in_hole(self):
     110        # create 3 right-angled triangles arranged in a bigger triangle       
     111        a = [0.0, 0.0] #0
     112        b = [0.0, 2.0] #1
     113        c = [2.0,0.0]  #2
     114        d = [0.0,4.0]  #3
     115        e = [2.0,2.0]  #4
     116        f = [4.0,0.0]  #5
     117        points = [a, b, c, d, e, f]
     118        vertices = [ [1,0,2], [3,1,4], [4,2,5] ]   #bac dbe ecf
     119
     120        point_in_hole = [1.5, 1.5] # a point in the hole
     121        data = [ [20, 20], [0.3, 0.3], point_in_hole, [2.5, 0.3], [30, 30] ] #some points inside and outside the hole
     122
     123        # any function for the vertices, we don't care about the result
     124        f = num.array([linear_function(points), 2*linear_function(points)])
     125        f = num.transpose(f)
     126
     127        interp = Interpolate(points, vertices)
     128        interp.interpolate(f, data)
     129       
     130        assert not set(interp.inside_poly_indices).intersection(set(interp.outside_poly_indices)), \
     131               'Some points are in both lists!'   
     132        assert len(interp.inside_poly_indices) == 2
     133        assert len(interp.outside_poly_indices) == 3
     134       
     135        interp.outside_poly_indices.sort()
     136        assert interp.outside_poly_indices[1] == 2, \
     137               'third outside point should be inside the hole!'
    110138
    111139    def test_simple_interpolation_example(self):
     
    391419        assert num.allclose(num.sum(results, axis=1), 1.0)
    392420
    393                
     421       
    394422    def test_arbitrary_datapoints_return_centroids(self):
    395423        #Try arbitrary datapoints, making sure they return
     
    409437        third = [1.0/3.0, 1.0/3.0, 1.0/3.0]
    410438        answer = [third, third, third]
    411                
     439       
    412440        A,_,_,_ = interp._build_interpolation_matrix_A(data, output_centroids=True)
    413441        results = A.todense()
    414         assert num.allclose(results, answer)           
    415                
    416                
     442        assert num.allclose(results, answer)       
     443       
     444       
    417445    def test_arbitrary_datapoints_some_outside(self):
    418446        #Try arbitrary datapoints one outside the triangle.
  • anuga_core/source/anuga/shallow_water/data_manager.py

    r7641 r7685  
    22582258    import types
    22592259
    2260     from anuga.utilities.polygon import inside_polygon, outside_polygon, \
    2261          separate_points_by_polygon
     2260    from anuga.utilities.polygon import inside_polygon, outside_polygon
    22622261    from anuga.abstract_2d_finite_volumes.util import \
    22632262         apply_expression_to_dictionary
     
    26902689
    26912690    import sys
    2692     from anuga.utilities.polygon import inside_polygon, outside_polygon, \
    2693              separate_points_by_polygon
     2691    from anuga.utilities.polygon import inside_polygon, outside_polygon
    26942692    from anuga.abstract_2d_finite_volumes.util import \
    26952693             apply_expression_to_dictionary
  • anuga_core/source/anuga/shallow_water/data_manager_joaquims_patch.py

    r7340 r7685  
    18471847    from Numeric import array2string
    18481848
    1849     from anuga.utilities.polygon import inside_polygon, outside_polygon, \
    1850          separate_points_by_polygon
     1849    from anuga.utilities.polygon import inside_polygon, outside_polygon
    18511850    from anuga.abstract_2d_finite_volumes.util import \
    18521851         apply_expression_to_dictionary
     
    22322231    from Numeric import array2string
    22332232
    2234     from anuga.utilities.polygon import inside_polygon, outside_polygon, separate_points_by_polygon
     2233    from anuga.utilities.polygon import inside_polygon, outside_polygon
    22352234    from anuga.abstract_2d_finite_volumes.util import apply_expression_to_dictionary
    22362235
Note: See TracChangeset for help on using the changeset viewer.