Changeset 6149
- Timestamp:
- Jan 13, 2009, 1:17:54 PM (16 years ago)
- Location:
- anuga_core/source/anuga/coordinate_transforms
- Files:
-
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/coordinate_transforms/geo_reference.py
r5736 r6149 8 8 9 9 import types, sys 10 from Numeric import array, Float, ArrayType, reshape, allclose11 10 from anuga.utilities.numerical_tools import ensure_numeric 12 11 from anuga.utilities.anuga_exceptions import ANUGAError, TitleError, ParsingError, \ 13 12 ShapeError 13 14 import Numeric as num 14 15 15 16 … … 98 99 99 100 # Fix some assertion failures 100 if type(self.zone) == ArrayType and self.zone.shape == ():101 if type(self.zone) == num.ArrayType and self.zone.shape == (): 101 102 self.zone = self.zone[0] 102 if type(self.xllcorner) == ArrayType and self.xllcorner.shape == ():103 if type(self.xllcorner) == num.ArrayType and self.xllcorner.shape == (): 103 104 self.xllcorner = self.xllcorner[0] 104 if type(self.yllcorner) == ArrayType and self.yllcorner.shape == ():105 if type(self.yllcorner) == num.ArrayType and self.yllcorner.shape == (): 105 106 self.yllcorner = self.yllcorner[0] 106 107 … … 172 173 173 174 # Fix some assertion failures 174 if(type(self.zone) == ArrayType and self.zone.shape == ()):175 if(type(self.zone) == num.ArrayType and self.zone.shape == ()): 175 176 self.zone = self.zone[0] 176 if type(self.xllcorner) == ArrayType and self.xllcorner.shape == ():177 if type(self.xllcorner) == num.ArrayType and self.xllcorner.shape == (): 177 178 self.xllcorner = self.xllcorner[0] 178 if type(self.yllcorner) == ArrayType and self.yllcorner.shape == ():179 if type(self.yllcorner) == num.ArrayType and self.yllcorner.shape == (): 179 180 self.yllcorner = self.yllcorner[0] 180 181 … … 196 197 is_list = True 197 198 198 points = ensure_numeric(points, Float)199 points = ensure_numeric(points, num.Float) 199 200 200 201 if len(points.shape) == 1: … … 202 203 msg = 'Single point must have two elements' 203 204 assert len(points) == 2, msg 204 points = reshape(points, (1,2))205 points = num.reshape(points, (1,2)) 205 206 206 207 msg = 'Points array must be two dimensional.\n' … … 239 240 """ 240 241 241 return allclose([self.xllcorner, self.yllcorner], 0)242 return num.allclose([self.xllcorner, self.yllcorner], 0) 242 243 243 244 … … 255 256 is_list = True 256 257 257 points = ensure_numeric(points, Float)258 points = ensure_numeric(points, num.Float) 258 259 if len(points.shape) == 1: 259 260 #One point has been passed … … 295 296 is_list = True 296 297 297 points = ensure_numeric(points, Float)298 points = ensure_numeric(points, num.Float) 298 299 if len(points.shape) == 1: 299 300 #One point has been passed -
anuga_core/source/anuga/coordinate_transforms/redfearn.py
r5662 r6149 8 8 """ 9 9 from anuga.coordinate_transforms.geo_reference import Geo_reference, DEFAULT_ZONE 10 from Numeric import array 10 11 11 12 12 def degminsec2decimal_degrees(dd,mm,ss): -
anuga_core/source/anuga/coordinate_transforms/test_geo_reference.py
r6086 r6149 7 7 8 8 from geo_reference import * 9 from Numeric import allclose,array10 9 from anuga.config import netcdf_mode_r, netcdf_mode_w, netcdf_mode_a 10 11 import Numeric as num 11 12 12 13 … … 165 166 y = 443.0 166 167 g = Geo_reference(56,x,y) 167 lofl = array([[3.0,323.0], [6.0,645.0]])168 lofl = num.array([[3.0,323.0], [6.0,645.0]]) 168 169 new_lofl = g.change_points_geo_ref(lofl) 169 170 #print "4 lofl",lofl 170 171 #print "4 new_lofl",new_lofl 171 172 172 self.failUnless(type(new_lofl) == ArrayType, ' failed')173 self.failUnless(type(new_lofl) == num.ArrayType, ' failed') 173 174 self.failUnless(type(new_lofl) == type(lofl), ' failed') 174 175 lofl[:,0] -= x 175 176 lofl[:,1] -= y 176 assert allclose(lofl,new_lofl)177 assert num.allclose(lofl,new_lofl) 177 178 178 179 def test_change_points_geo_ref5(self): … … 180 181 y = 3.0 181 182 g = Geo_reference(56,x,y) 182 lofl = array([[3.0,323.0]])183 lofl = num.array([[3.0,323.0]]) 183 184 184 185 … … 188 189 #print "5 new_lofl",new_lofl 189 190 190 self.failUnless(type(new_lofl) == ArrayType, ' failed')191 self.failUnless(type(new_lofl) == num.ArrayType, ' failed') 191 192 self.failUnless(type(new_lofl) == type(lofl), ' failed') 192 193 … … 200 201 y = 3.0 201 202 g = Geo_reference(56,x,y) 202 lofl = array([355.0,3.0])203 lofl = num.array([355.0,3.0]) 203 204 new_lofl = g.change_points_geo_ref(lofl.copy()) 204 205 #print "lofl",lofl 205 206 #print "new_lofl",new_lofl 206 207 207 self.failUnless(type(new_lofl) == ArrayType, ' failed')208 self.failUnless(type(new_lofl) == num.ArrayType, ' failed') 208 209 self.failUnless(type(new_lofl) == type(lofl), ' failed') 209 210 for point,new_point in map(None,[lofl],new_lofl): -
anuga_core/source/anuga/coordinate_transforms/test_lat_long_UTM_conversion.py
r4197 r6149 7 7 8 8 import unittest 9 from Numeric import allclose10 9 11 10 from lat_long_UTM_conversion import * 12 11 from redfearn import degminsec2decimal_degrees, decimal_degrees2degminsec 13 12 from anuga.utilities.anuga_exceptions import ANUGAError 13 14 import Numeric as num 15 14 16 15 17 #------------------------------------------------------------- … … 28 30 lat = degminsec2decimal_degrees(-37,39,10.15610) 29 31 lon = degminsec2decimal_degrees(143,55,35.38390) 30 assert allclose(lat, -37.65282114)31 assert allclose(lon, 143.9264955)32 assert num.allclose(lat, -37.65282114) 33 assert num.allclose(lon, 143.9264955) 32 34 33 35 … … 35 37 36 38 assert zone == 54 37 assert allclose(easting, 758173.797)38 assert allclose(northing, 5828674.340)39 assert num.allclose(easting, 758173.797) 40 assert num.allclose(northing, 5828674.340) 39 41 40 42 lat_calced, long_calced = UTMtoLL(northing, easting, zone) 41 assert allclose(lat, lat_calced)42 assert allclose(lon, long_calced)43 assert num.allclose(lat, lat_calced) 44 assert num.allclose(lon, long_calced) 43 45 44 46 … … 61 63 62 64 assert zone == 55 63 assert allclose(easting, 273741.297)64 assert allclose(northing, 5796489.777)65 assert num.allclose(easting, 273741.297) 66 assert num.allclose(northing, 5796489.777) 65 67 66 68 lat_calced, long_calced = UTMtoLL(northing, easting, zone) 67 assert allclose(lat, lat_calced)68 assert allclose(lon, long_calced)69 assert num.allclose(lat, lat_calced) 70 assert num.allclose(lon, long_calced) 69 71 70 72 … … 78 80 79 81 assert zone == 52 80 assert allclose(easting, 555776.267)81 assert allclose(northing, 3348167.264)82 assert num.allclose(easting, 555776.267) 83 assert num.allclose(northing, 3348167.264) 82 84 83 85 Lat, Long = UTMtoLL(northing, easting, zone) … … 94 96 95 97 assert zone == 33 96 assert allclose(easting, 348157.631)97 assert allclose(northing, 6175612.993)98 assert num.allclose(easting, 348157.631) 99 assert num.allclose(northing, 6175612.993) 98 100 99 101 lat_calced, long_calced = UTMtoLL(northing, easting, zone, 100 102 isSouthernHemisphere=False) 101 assert allclose(lat, lat_calced)102 assert allclose(lon, long_calced)103 assert num.allclose(lat, lat_calced) 104 assert num.allclose(lon, long_calced) 103 105 104 106 def test_UTM_5(self): … … 113 115 114 116 assert zone == 56 115 assert allclose(easting, 308728.009)116 assert allclose(northing, 6180432.601)117 assert num.allclose(easting, 308728.009) 118 assert num.allclose(northing, 6180432.601) 117 119 118 120 lat_calced, long_calced = UTMtoLL(northing, easting, zone) 119 assert allclose(lat, lat_calced)120 assert allclose(lon, long_calced)121 assert num.allclose(lat, lat_calced) 122 assert num.allclose(lon, long_calced) 121 123 #------------------------------------------------------------- 122 124 if __name__ == "__main__": -
anuga_core/source/anuga/coordinate_transforms/test_redfearn.py
r5663 r6149 7 7 8 8 import unittest 9 from Numeric import allclose10 9 11 10 from redfearn import * 12 11 from anuga.utilities.anuga_exceptions import ANUGAError 12 13 import Numeric as num 14 13 15 14 16 #------------------------------------------------------------- … … 19 21 lat = degminsec2decimal_degrees(-37,39,10.15610) 20 22 lon = degminsec2decimal_degrees(143,55,35.38390) 21 assert allclose(lat, -37.65282114)22 assert allclose(lon, 143.9264955)23 assert num.allclose(lat, -37.65282114) 24 assert num.allclose(lon, 143.9264955) 23 25 24 26 dd,mm,ss = decimal_degrees2degminsec(-37.65282114) 25 27 assert dd==-37 26 28 assert mm==39 27 assert allclose(ss, 10.15610)29 assert num.allclose(ss, 10.15610) 28 30 29 31 dd,mm,ss = decimal_degrees2degminsec(143.9264955) 30 32 assert dd==143 31 33 assert mm==55 32 assert allclose(ss, 35.38390)34 assert num.allclose(ss, 35.38390) 33 35 34 36 … … 44 46 lat = degminsec2decimal_degrees(-37,39,10.15610) 45 47 lon = degminsec2decimal_degrees(143,55,35.38390) 46 assert allclose(lat, -37.65282114)47 assert allclose(lon, 143.9264955)48 assert num.allclose(lat, -37.65282114) 49 assert num.allclose(lon, 143.9264955) 48 50 49 51 … … 51 53 52 54 assert zone == 54 53 assert allclose(easting, 758173.797)54 assert allclose(northing, 5828674.340)55 assert num.allclose(easting, 758173.797) 56 assert num.allclose(northing, 5828674.340) 55 57 56 58 … … 72 74 73 75 assert zone == 55 74 assert allclose(easting, 273741.297)75 assert allclose(northing, 5796489.777)76 assert num.allclose(easting, 273741.297) 77 assert num.allclose(northing, 5796489.777) 76 78 77 79 … … 85 87 86 88 assert zone == 52 87 assert allclose(easting, 555776.267)88 assert allclose(northing, 3348167.264)89 assert num.allclose(easting, 555776.267) 90 assert num.allclose(northing, 3348167.264) 89 91 90 92 … … 101 103 102 104 assert zone == 33 103 assert allclose(easting, 348157.631)104 assert allclose(northing, 6175612.993)105 assert num.allclose(easting, 348157.631) 106 assert num.allclose(northing, 6175612.993) 105 107 106 108 … … 116 118 117 119 assert zone == 56 118 assert allclose(easting, 308728.009)119 assert allclose(northing, 6180432.601)120 assert num.allclose(easting, 308728.009) 121 assert num.allclose(northing, 6180432.601) 120 122 121 123 def test_UTM_6_nonstandard_projection(self): … … 126 128 127 129 assert zone == 50 128 assert allclose(easting, 213251.040253)129 assert allclose(northing, 6762559.15978)130 assert num.allclose(easting, 213251.040253) 131 assert num.allclose(northing, 6762559.15978) 130 132 131 133 #Testing using the native zone … … 133 135 134 136 assert zone == 50 135 assert allclose(easting, 213251.040253)136 assert allclose(northing, 6762559.15978)137 assert num.allclose(easting, 213251.040253) 138 assert num.allclose(northing, 6762559.15978) 137 139 138 140 #Then project to zone 49 … … 140 142 141 143 assert zone == 49 142 assert allclose(easting, 796474.020057)143 assert allclose(northing, 6762310.25162)144 assert num.allclose(easting, 796474.020057) 145 assert num.allclose(northing, 6762310.25162) 144 146 145 147 … … 151 153 152 154 assert zone == 49 153 assert allclose(easting, 788653.192779)154 assert allclose(northing, 6773605.46384)155 assert num.allclose(easting, 788653.192779) 156 assert num.allclose(northing, 6773605.46384) 155 157 156 158 #Then project to zone 50 … … 158 160 159 161 assert zone == 50 160 assert allclose(easting, 204863.606467)161 assert allclose(northing, 6773440.04726)162 assert num.allclose(easting, 204863.606467) 163 assert num.allclose(northing, 6773440.04726) 162 164 163 165 #Testing point on zone boundary … … 166 168 167 169 assert zone == 50 168 assert allclose(easting, 208199.768268)169 assert allclose(northing, 6769820.01453)170 assert num.allclose(easting, 208199.768268) 171 assert num.allclose(northing, 6769820.01453) 170 172 171 173 #Then project to zone 49 … … 173 175 174 176 assert zone == 49 175 assert allclose(easting, 791800.231817)176 assert allclose(northing, 6769820.01453)177 assert num.allclose(easting, 791800.231817) 178 assert num.allclose(northing, 6769820.01453) 177 179 178 180 #Testing furthest point in Geraldton scenario) … … 181 183 182 184 assert zone == 49 183 assert allclose(easting, 737178.16131)184 assert allclose(northing, 6876426.38578)185 assert num.allclose(easting, 737178.16131) 186 assert num.allclose(northing, 6876426.38578) 185 187 186 188 #Then project to zone 50 … … 188 190 189 191 assert zone == 50 190 assert allclose(easting, 148260.567427)191 assert allclose(northing, 6873587.50926)192 assert num.allclose(easting, 148260.567427) 193 assert num.allclose(northing, 6873587.50926) 192 194 193 195 #Testing outside GDA zone (New Zeland) … … 196 198 197 199 assert zone == 60 198 assert allclose(easting, 580174.259843)199 assert allclose(northing, 5127641.114461)200 assert num.allclose(easting, 580174.259843) 201 assert num.allclose(northing, 5127641.114461) 200 202 201 203 #Then project to zone 59 … … 203 205 204 206 assert zone == 59 205 assert allclose(easting, 1061266.922118)206 assert allclose(northing, 5104249.395469)207 assert num.allclose(easting, 1061266.922118) 208 assert num.allclose(northing, 5104249.395469) 207 209 208 210 #Then skip three zones 57 (native 60) … … 210 212 211 213 assert zone == 57 212 assert allclose(easting, 2023865.527497)213 assert allclose(northing, 4949253.934967)214 assert num.allclose(easting, 2023865.527497) 215 assert num.allclose(northing, 4949253.934967) 214 216 215 217 # Note Projecting into the Northern Hemisphere Does not coincide … … 227 229 # Google Earth interpretation 228 230 assert zone == 57 229 assert allclose(easting, 259473.69)230 assert allclose(northing, 4876249.13)231 assert num.allclose(easting, 259473.69) 232 assert num.allclose(northing, 4876249.13) 231 233 232 234 # ArcMap's interpretation 233 235 #assert zone == 57 234 #assert allclose(easting, 259473.678944)235 #assert allclose(northing, 14876249.1268)236 #assert num.allclose(easting, 259473.678944) 237 #assert num.allclose(northing, 14876249.1268) 236 238 237 239 #Then project to zone 56 … … 239 241 240 242 assert zone == 56 241 assert allclose(easting, 740526.321055)242 assert allclose(northing, 4876249.13)243 assert num.allclose(easting, 740526.321055) 244 assert num.allclose(northing, 4876249.13) 243 245 244 246 … … 283 285 points, zone = convert_from_latlon_to_utm(latitudes=lats, longitudes=longs) 284 286 285 assert allclose(points[0][0], 308728.009)286 assert allclose(points[0][1], 6180432.601)287 assert allclose(points[1][0], 222908.705)288 assert allclose(points[1][1], 6233785.284)287 assert num.allclose(points[0][0], 308728.009) 288 assert num.allclose(points[0][1], 6180432.601) 289 assert num.allclose(points[1][0], 222908.705) 290 assert num.allclose(points[1][1], 6233785.284) 289 291 self.failUnless(zone == 56, 290 292 'Bad zone error!') … … 355 357 points, zone = convert_from_latlon_to_utm(points=points) 356 358 #print "points",points 357 assert allclose(points[0][0], 308728.009)358 assert allclose(points[0][1], 6180432.601)359 assert allclose(points[1][0], 222908.705)360 assert allclose(points[1][1], 6233785.284)359 assert num.allclose(points[0][0], 308728.009) 360 assert num.allclose(points[0][1], 6180432.601) 361 assert num.allclose(points[1][0], 222908.705) 362 assert num.allclose(points[1][1], 6233785.284) 361 363 self.failUnless(zone == 56, 362 364 'Bad zone error!')
Note: See TracChangeset
for help on using the changeset viewer.