Changeset 4238

Timestamp:

Feb 7, 2007, 5:14:03 PM (18 years ago)

Author:

duncan

Message:

more for updated urs2sww

Files:

• anuga_core/source/anuga/geospatial_data/geospatial_data.py (modified) (2 diffs)
• anuga_core/source/anuga/geospatial_data/test_geospatial_data.py (modified) (3 diffs)
• anuga_core/source/anuga/shallow_water/data_manager.py (modified) (2 diffs)
• anuga_core/source/anuga/shallow_water/test_data_manager.py (modified) (2 diffs)
• anuga_work/development/friction_dam_2006/run_dam.py (modified) (1 diff)

anuga_core/source/anuga/geospatial_data/geospatial_data.py

@@ -608 +608 @@
                             self.get_all_attributes(),
                             as_lat_long=as_lat_long)
-                                    
+              
+        elif file_name[-4:] == ".urs" :
+            msg = "ERROR: Can not write a .urs file as a relative file."
+            assert absolute, msg
+            _write_urs_file(file_name,
+                            self.get_data_points(as_lat_long=True))
+                                                          
         else:
             msg = 'Unknown file type %s ' %file_name
@@ -1223 +1229 @@
 
     fd.close()
+  
+def _write_urs_file(file_name,
+                    points,
+                    delimiter=' '):
+    """
+    export a file, file_name, with the urs format
+    the data points are in lats and longs
+    
+    """    
+    fd = open(file_name,'w')   
+    fd.write(str(len(points))+"\n")   
+    # <lat> <long> <id#>
+    for i, vert in enumerate( points):
+        fd.write(str(round(vert[0],7)) + delimiter + \
+                 str(round(vert[1],7)) + delimiter +str(i+1)+ "\n")
+    fd.close()
     
 def _point_atts2array(point_atts):

anuga_core/source/anuga/geospatial_data/test_geospatial_data.py

@@ -2035 +2035 @@
 
 
+    def test_write_urs_file(self):
+        lat_gong = degminsec2decimal_degrees(-34,00,0)
+        lon_gong = degminsec2decimal_degrees(150,30,0.)
+        
+        lat_2 = degminsec2decimal_degrees(-34,00,1)
+        lon_2 = degminsec2decimal_degrees(150,00,0.)
+        p1 = (lat_gong, lon_gong)
+        p2 = (lat_2, lon_2)
+        points = ImmutableSet([p1, p2, p1])
+        gsd = Geospatial_data(data_points=list(points),
+                              points_are_lats_longs=True)
+        
+        fn = tempfile.mktemp(".urs")
+        gsd.export_points_file(fn)
+        #print "fn", fn 
+        handle = open(fn)
+        lines = handle.readlines()
+        assert lines[0],'2'
+        assert lines[1],'-34.0002778 150.0 1'
+        assert lines[2],'-34.0 150.5 2'
+        handle.close()
+        os.remove(fn)
+        
     def test_lat_long_set(self):
         lat_gong = degminsec2decimal_degrees(-34,30,0.)
@@ -2096 +2119 @@
                   [5.0, 1.0], [5.0, 2.0],[5.0, 3.0], [5.0, 4.0], [5.0, 5.0]]
         attributes = {'depth':[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 
-                      14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25],'speed':[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 
+                      14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25],
+                      'speed':[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 
                       14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25]}
         G = Geospatial_data(points, attributes)
@@ -2104 +2128 @@
         factor = 0.21
 
-        G1, G2  = G.split(factor) #will return G1 with 10% of points and G2 with 90%
-        
-#        print 'len(G): %s  len(G1): %s len(G2): %s' %(len(G), len(G1), len(G2))
+        #will return G1 with 10% of points and G2 with 90%
+        G1, G2  = G.split(factor) 
+        
+#        print 'len(G): %s len(G1): %s len(G2): %s' %(len(G), len(G1), len(G2))
 #        print 'G: ', len(G),'G1: ', len(G1), 'G2: ', len(G2)
 

anuga_core/source/anuga/shallow_water/data_manager.py

@@ -4409 +4409 @@
     #### URS UNGRIDDED 2 SWW ###
 
-def URS_points_needed_to_file(boundary_polygon, ll_lat, ll_long,
+def URS_points_needed_to_file(file_name, boundary_polygon, ll_lat, ll_long,
                               grid_spacing, 
                       lat_amount, long_amount, zone=None):
-    
+    """
+    file_name - name of the urs file produced for David.
+    boundary_polygon - a list of points that describes a polygon.
+                      The last point is assumed ot join the first point.
+                      This is in UTM (lat long would be better though)
+
+    ll_lat - lower left latitude, in decimal degrees
+    ll-long - lower left longitude, in decimal degrees
+    grid_spacing - in deciamal degrees
+
+
+    Don't add the file extension.  It will be added.
+    """
     geo = URS_points_needed(boundary_polygon, ll_lat, ll_long, grid_spacing, 
                       lat_amount, long_amount, zone=zone)
+    if not file[-4:] == ".urs":
+        file += ".urs"
+    geo.export_points_file(file_name)
     
 def URS_points_needed(boundary_polygon, ll_lat, ll_long, grid_spacing, 
@@ -4471 +4486 @@
     max_lat = max(seg_lat_long[0][0], seg_lat_long[1][0]) + buffer
     max_long = max(seg_lat_long[0][1], seg_lat_long[1][1]) + buffer
-    min_lat = max(seg_lat_long[0][0], seg_lat_long[1][0]) - buffer
-    min_long = max(seg_lat_long[0][1], seg_lat_long[1][1]) - buffer
+    min_lat = min(seg_lat_long[0][0], seg_lat_long[1][0]) - buffer
+    min_long = min(seg_lat_long[0][1], seg_lat_long[1][1]) - buffer
 
     #print "ll_lat", ll_lat

anuga_core/source/anuga/shallow_water/test_data_manager.py

@@ -5040 +5040 @@
 
     #### TESTS URS UNGRIDDED 2 SWW ###
-    def X_test_URS_points_needed(self):
+    def test_URS_points_needed(self):
         
         ll_lat = -21.5
@@ -5057 +5057 @@
         #Maybe see if I can fit the data to the polygon - have to represent
         # the poly as points though.
-        geo.export_points_file("results")
+        #geo.export_points_file("results.txt", as_lat_long=True)
         results = ImmutableSet(geo.get_data_points(as_lat_long=True))
-        print 'results',results
+        #print 'results',results
 
         # These are a set of points that have to be in results
         points = []
         for i in range(18):
-            lat = -21.0 - 8/60 - grid_spacing * i
+            lat = -21.0 - 8./60 - grid_spacing * i
             points.append((lat,degminsec2decimal_degrees(114,35,0))) 
-            #points.append((lat,degminsec2decimal_degrees(114,36,0))) 
-            #points.append((lat,degminsec2decimal_degrees(114,52,0))) 
-            #points.append((lat,degminsec2decimal_degrees(114,53,0)))
-            
+            points.append((lat,degminsec2decimal_degrees(114,36,0))) 
+            points.append((lat,degminsec2decimal_degrees(114,52,0))) 
+            points.append((lat,degminsec2decimal_degrees(114,53,0)))
+        geo_answer = Geospatial_data(data_points=points,
+                                     points_are_lats_longs=True) 
+        #geo_answer.export_points_file("answer.txt", as_lat_long=True)  
         answer = ImmutableSet(points)
-        print "answer", answer
+        
         outs = answer.difference(results)
-        print "outs", outs 
-        assert answer.issubset(results)
+        #print "outs", outs
+        # This doesn't work.  Though visualising the results shows that
+        # it is correct.
+        #assert answer.issubset(results)
+        # this is why;
+        #point (-21.133333333333333, 114.58333333333333)
+        #result (-21.133333332232368, 114.58333333300342)
+        
+        for point in points:
+            found = False
+            for result in results:
+                if allclose(point, result):
+                    found = True
+                    break
+            if not found:
+                assert False
+        
         
     def X_test_URS_points_neededII(self):

anuga_work/development/friction_dam_2006/run_dam.py

@@ -153 +153 @@
 if __name__ == "__main__":
     for friction in [ 0.05]:
-        main(friction, is_trial_run = False, outputdir_name='friction_set_A')
+        main(friction, is_trial_run = True, outputdir_name='friction_set_A')
     
     #    to do