Changeset 5254


Ignore:
Timestamp:
Apr 30, 2008, 4:34:19 PM (16 years ago)
Author:
duncan
Message:

Adding northern hemisphere tests for URS_points_needed.

Location:
anuga_core/source/anuga/shallow_water
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/shallow_water/data_manager.py

    r5253 r5254  
    44394439    ll-long - lower left longitude, in decimal degrees
    44404440    grid_spacing - in deciamal degrees
    4441     lat_amount - number of latitudes
    4442     long_amount- number of longs
     4441    lat_amount - number of latitudes 
     4442    long_amount- number of longs 
    44434443
    44444444
     
    44564456            file_name = file_name[:-4] + ".csv"
    44574457        geo.export_points_file(file_name)
     4458    return geo
    44584459
    44594460def URS_points_needed(boundary_polygon, zone, ll_lat,
     
    45314532                  isSouthHemisphere):
    45324533    """
    4533     seg is one point, in UTM
     4534    seg is two points, in UTM
    45344535    return a list of the points, in lats and longs that are needed to
    45354536    interpolate any point on the segment.
     
    45384539    #print "zone",zone
    45394540    geo_reference = Geo_reference(zone=zone)
     4541    #print "seg", seg
    45404542    geo = Geospatial_data(seg,geo_reference=geo_reference)
    45414543    seg_lat_long = geo.get_data_points(as_lat_long=True,
    45424544                                       isSouthHemisphere=isSouthHemisphere)
     4545    #print "seg_lat_long", seg_lat_long
    45434546    # 1.415 = 2^0.5, rounded up....
    45444547    sqrt_2_rounded_up = 1.415
     
    45504553    min_long = min(seg_lat_long[0][1], seg_lat_long[1][1]) - buffer
    45514554
     4555    #print "min_long", min_long
     4556    #print "ll_long", ll_long
     4557    #print "grid_spacing", grid_spacing
    45524558    first_row = (min_long - ll_long)/grid_spacing
    45534559    # To round up
     
    45774583    for index_lat in range(first_row_lat, last_row_lat + 1):
    45784584        for index_long in range(first_row_long, last_row_long + 1):
     4585            #print "index_lat", index_lat
     4586            #print "index_long", index_long
    45794587            lat = ll_lat + index_lat*grid_spacing
    45804588            long = ll_long + index_long*grid_spacing
     
    45984606    x2 = seg[1][0]
    45994607    y2 = seg[1][1]
    4600 
    46014608    x2_1 = x2-x1
    46024609    y2_1 = y2-y1
    4603     d = abs((x2_1)*(y1-y0)-(x1-x0)*(y2_1))/sqrt((x2_1)*(x2_1)+(y2_1)*(y2_1))
     4610    try:
     4611        d = abs((x2_1)*(y1-y0)-(x1-x0)*(y2_1))/sqrt( \
     4612            (x2_1)*(x2_1)+(y2_1)*(y2_1))
     4613    except ZeroDivisionError:
     4614        #print "seg", seg
     4615        #print "x0", x0
     4616        #print "y0", y0
     4617        if sqrt((x2_1)*(x2_1)+(y2_1)*(y2_1)) == 0 and \
     4618           abs((x2_1)*(y1-y0)-(x1-x0)*(y2_1)) == 0:
     4619            return True
     4620        else:
     4621            return False
     4622   
    46044623    if d <= max_distance:
    46054624        return True
  • anuga_core/source/anuga/shallow_water/test_data_manager.py

    r5253 r5254  
    61176117                          lat_amount, long_amount,
    61186118                          verbose=self.verbose)
    6119 
    6120     def test_URS_points_needed_poly1(self):
     6119       
     6120    def test_URS_points_northern_hemisphere(self):
     6121               
     6122        LL_LAT = 8.0
     6123        LL_LONG = 97.0
     6124        GRID_SPACING = 2.0/60.0
     6125        LAT_AMOUNT = 2
     6126        LONG_AMOUNT = 2
     6127        ZONE = 47
     6128
     6129        #
     6130        points = []
     6131        for i in range(2):
     6132            for j in range(2):
     6133                points.append((degminsec2decimal_degrees(8,1+i*2,0),
     6134                               degminsec2decimal_degrees(97,1+i*2,0)))
     6135        #print "points", points
     6136        geo_poly = Geospatial_data(data_points=points,
     6137                                     points_are_lats_longs=True)
     6138        poly_lat_long = geo_poly.get_data_points(as_lat_long=False,
     6139                                       isSouthHemisphere=False)
     6140        #print "seg_lat_long",  poly_lat_long
     6141       
     6142        #geo=URS_points_needed_to_file('test_example_poly3', poly_lat_long,
     6143        geo=URS_points_needed(poly_lat_long,
     6144                                  ZONE,
     6145                                  LL_LAT, LL_LONG,
     6146                                  GRID_SPACING,
     6147                                  LAT_AMOUNT, LONG_AMOUNT,
     6148                                  isSouthHemisphere=False,
     6149                                  verbose=self.verbose)
     6150        results = ImmutableSet(geo.get_data_points(as_lat_long=True,
     6151                                  isSouthHemisphere=False))
     6152        #print 'results',results
     6153
     6154        # These are a set of points that have to be in results
     6155        points = []
     6156        for i in range(2):
     6157            for j in range(2):
     6158                points.append((degminsec2decimal_degrees(8,i*2,0),
     6159                               degminsec2decimal_degrees(97,i*2,0)))
     6160        #print "answer points", points
     6161        answer = ImmutableSet(points)
     6162       
     6163        for point in points:
     6164            found = False
     6165            for result in results:
     6166                if allclose(point, result):
     6167                    found = True
     6168                    break
     6169            if not found:
     6170                assert False
     6171       
     6172
     6173    def covered_in_other_tests_test_URS_points_needed_poly1(self):
    61216174        # Values used for FESA 2007 results
    61226175        # domain in southern hemisphere zone 51       
     
    61416194                                  GRID_SPACING,
    61426195                                  LAT_AMOUNT, LONG_AMOUNT,
    6143                                   verbose=self.verbose)         
    6144 
    6145 
    6146     def test_URS_points_needed_poly2(self):
     6196                                  verbose=self.verbose)
     6197       
     6198
     6199
     6200    def covered_in_other_tests_test_URS_points_needed_poly2(self):
    61476201        # Values used for 2004 validation work
    61486202        # domain in northern hemisphere zone 47       
     
    75627616    suite = unittest.makeSuite(Test_Data_Manager,'test')
    75637617    #suite = unittest.makeSuite(Test_Data_Manager,'test_get_flow_through_cross_section')
    7564     #suite = unittest.makeSuite(Test_Data_Manager,'Xtest')
     7618    #suite = unittest.makeSuite(Test_Data_Manager,'covered_')
    75657619
    75667620   
Note: See TracChangeset for help on using the changeset viewer.