Changeset 3761


Ignore:
Timestamp:
Oct 11, 2006, 5:46:10 PM (18 years ago)
Author:
ole
Message:

Removed old tests in data_manager as per ticket:86 - they were old versions of dem2pts tests that had been superseded. All good.

Location:
anuga_core/source
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • anuga_core/source/anuga/shallow_water/test_data_manager.py

    r3750 r3761  
    42694269
    42704270
    4271     def NOT_test_dem2pts(self):
    4272         """Test conversion from dem in ascii format to native NetCDF xya format
    4273         """
    4274 
    4275         import time, os
    4276         from Numeric import array, zeros, allclose, Float, concatenate
    4277         from Scientific.IO.NetCDF import NetCDFFile
    4278 
    4279         #Write test asc file
    4280         root = 'demtest'
    4281 
    4282         filename = root+'.asc'
    4283         fid = open(filename, 'w')
    4284         fid.write("""ncols         5
    4285 nrows         6
    4286 xllcorner     2000.5
    4287 yllcorner     3000.5
    4288 cellsize      25
    4289 NODATA_value  -9999
    4290 """)
    4291         #Create linear function
    4292 
    4293         ref_points = []
    4294         ref_elevation = []
    4295         for i in range(6):
    4296             y = (6-i)*25.0
    4297             for j in range(5):
    4298                 x = j*25.0
    4299                 z = x+2*y
    4300 
    4301                 ref_points.append( [x,y] )
    4302                 ref_elevation.append(z)
    4303                 fid.write('%f ' %z)
    4304             fid.write('\n')
    4305 
    4306         fid.close()
    4307 
    4308         #Write prj file with metadata
    4309         metafilename = root+'.prj'
    4310         fid = open(metafilename, 'w')
    4311 
    4312 
    4313         fid.write("""Projection UTM
    4314 Zone 56
    4315 Datum WGS84
    4316 Zunits NO
    4317 Units METERS
    4318 Spheroid WGS84
    4319 Xshift 0.0000000000
    4320 Yshift 10000000.0000000000
    4321 Parameters
    4322 """)
    4323         fid.close()
    4324 
    4325         #Convert to NetCDF pts
    4326         convert_dem_from_ascii2netcdf(root)
    4327         dem2pts(root)
    4328 
    4329         #Check contents
    4330         #Get NetCDF
    4331         fid = NetCDFFile(root+'.pts', 'r')
    4332 
    4333         # Get the variables
    4334         #print fid.variables.keys()
    4335         points = fid.variables['points']
    4336         elevation = fid.variables['elevation']
    4337 
    4338         #Check values
    4339 
    4340         #print points[:]
    4341         #print ref_points
    4342         assert allclose(points, ref_points)
    4343 
    4344         #print attributes[:]
    4345         #print ref_elevation
    4346         assert allclose(elevation, ref_elevation)
    4347 
    4348         #Cleanup
    4349         fid.close()
    4350 
    4351 
    4352         os.remove(root + '.pts')
    4353         os.remove(root + '.dem')
    4354         os.remove(root + '.asc')
    4355         os.remove(root + '.prj')
    4356 
    4357 
    4358 
    4359     def NOT_test_dem2pts_bounding_box(self):
    4360         """Test conversion from dem in ascii format to native NetCDF xya format
    4361         """
    4362 
    4363         import time, os
    4364         from Numeric import array, zeros, allclose, Float, concatenate
    4365         from Scientific.IO.NetCDF import NetCDFFile
    4366 
    4367         #Write test asc file
    4368         root = 'demtest'
    4369 
    4370         filename = root+'.asc'
    4371         fid = open(filename, 'w')
    4372         fid.write("""ncols         5
    4373 nrows         6
    4374 xllcorner     2000.5
    4375 yllcorner     3000.5
    4376 cellsize      25
    4377 NODATA_value  -9999
    4378 """)
    4379         #Create linear function
    4380 
    4381         ref_points = []
    4382         ref_elevation = []
    4383         for i in range(6):
    4384             y = (6-i)*25.0
    4385             for j in range(5):
    4386                 x = j*25.0
    4387                 z = x+2*y
    4388 
    4389                 ref_points.append( [x,y] )
    4390                 ref_elevation.append(z)
    4391                 fid.write('%f ' %z)
    4392             fid.write('\n')
    4393 
    4394         fid.close()
    4395 
    4396         #Write prj file with metadata
    4397         metafilename = root+'.prj'
    4398         fid = open(metafilename, 'w')
    4399 
    4400 
    4401         fid.write("""Projection UTM
    4402 Zone 56
    4403 Datum WGS84
    4404 Zunits NO
    4405 Units METERS
    4406 Spheroid WGS84
    4407 Xshift 0.0000000000
    4408 Yshift 10000000.0000000000
    4409 Parameters
    4410 """)
    4411         fid.close()
    4412 
    4413         #Convert to NetCDF pts
    4414         convert_dem_from_ascii2netcdf(root)
    4415         dem2pts(root, easting_min=2010.0, easting_max=2110.0,
    4416                 northing_min=3035.0, northing_max=3125.5)
    4417 
    4418         #Check contents
    4419         #Get NetCDF
    4420         fid = NetCDFFile(root+'.pts', 'r')
    4421 
    4422         # Get the variables
    4423         #print fid.variables.keys()
    4424         points = fid.variables['points']
    4425         elevation = fid.variables['elevation']
    4426 
    4427         #Check values
    4428         assert fid.xllcorner[0] == 2010.0
    4429         assert fid.yllcorner[0] == 3035.0
    4430 
    4431         #create new reference points
    4432         ref_points = []
    4433         ref_elevation = []
    4434         for i in range(4):
    4435             y = (4-i)*25.0 + 25.0
    4436             y_new = y + 3000.5 - 3035.0
    4437             for j in range(4):
    4438                 x = j*25.0 + 25.0
    4439                 x_new = x + 2000.5 - 2010.0
    4440                 z = x+2*y
    4441 
    4442                 ref_points.append( [x_new,y_new] )
    4443                 ref_elevation.append(z)
    4444 
    4445         #print points[:]
    4446         #print ref_points
    4447         assert allclose(points, ref_points)
    4448 
    4449         #print attributes[:]
    4450         #print ref_elevation
    4451         assert allclose(elevation, ref_elevation)
    4452 
    4453         #Cleanup
    4454         fid.close()
    4455 
    4456 
    4457         os.remove(root + '.pts')
    4458         os.remove(root + '.dem')
    4459         os.remove(root + '.asc')
    4460         os.remove(root + '.prj')
    4461 
    4462 
    4463 
    4464     def NOT_test_dem2pts_remove_Nullvalues(self):
    4465         """Test conversion from dem in ascii format to native NetCDF xya format
    4466         """
    4467 
    4468         import time, os
    4469         from Numeric import array, zeros, allclose, Float, concatenate
    4470         from Scientific.IO.NetCDF import NetCDFFile
    4471 
    4472         #Write test asc file
    4473         root = 'demtest'
    4474 
    4475         filename = root+'.asc'
    4476         fid = open(filename, 'w')
    4477         fid.write("""ncols         5
    4478 nrows         6
    4479 xllcorner     2000.5
    4480 yllcorner     3000.5
    4481 cellsize      25
    4482 NODATA_value  -9999
    4483 """)
    4484         #Create linear function
    4485         # ref_ will write all the values
    4486         # new_ref_ will write the values except for NODATA_values
    4487         ref_points = []
    4488         ref_elevation = []
    4489         new_ref_pts = []
    4490         new_ref_elev = []
    4491         NODATA_value = -9999
    4492         for i in range(6):
    4493             y = (6-i)*25.0
    4494             for j in range(5):
    4495                 x = j*25.0
    4496                 z = x+2*y
    4497                 if j == 4: z = NODATA_value # column
    4498                 if i == 2 and j == 2: z = NODATA_value # random
    4499                 if i == 5 and j == 1: z = NODATA_value
    4500                 if i == 1: z = NODATA_value # row
    4501                 if i == 3 and j == 1: z = NODATA_value # two pts/row
    4502                 if i == 3 and j == 3: z = NODATA_value
    4503 
    4504 
    4505                 if z <> NODATA_value:
    4506                     new_ref_elev.append(z)
    4507                     new_ref_pts.append( [x,y] )
    4508 
    4509                 ref_points.append( [x,y] )
    4510                 ref_elevation.append(z)
    4511 
    4512                 fid.write('%f ' %z)
    4513             fid.write('\n')
    4514 
    4515         fid.close()
    4516 
    4517 
    4518         #Write prj file with metadata
    4519         metafilename = root+'.prj'
    4520         fid = open(metafilename, 'w')
    4521 
    4522 
    4523         fid.write("""Projection UTM
    4524 Zone 56
    4525 Datum WGS84
    4526 Zunits NO
    4527 Units METERS
    4528 Spheroid WGS84
    4529 Xshift 0.0000000000
    4530 Yshift 10000000.0000000000
    4531 Parameters
    4532 """)
    4533         fid.close()
    4534 
    4535         #Convert to NetCDF pts
    4536         convert_dem_from_ascii2netcdf(root)
    4537         dem2pts(root)
    4538 
    4539         #Check contents
    4540         #Get NetCDF
    4541         fid = NetCDFFile(root+'.pts', 'r')
    4542 
    4543         # Get the variables
    4544         #print fid.variables.keys()
    4545         points = fid.variables['points']
    4546         elevation = fid.variables['elevation']
    4547 
    4548         #Check values
    4549         #print 'points', points[:]
    4550         assert len(points) == len(new_ref_pts), 'length of returned points not correct'
    4551         assert allclose(points, new_ref_pts), 'points do not align'
    4552 
    4553         #print 'elevation', elevation[:]
    4554         assert len(elevation) == len(new_ref_elev), 'length of returned elevation not correct'
    4555         assert allclose(elevation, new_ref_elev), 'elevations do not align'
    4556 
    4557         #Cleanup
    4558         fid.close()
    4559 
    4560 
    4561         os.remove(root + '.pts')
    4562         os.remove(root + '.dem')
    4563         os.remove(root + '.asc')
    4564         os.remove(root + '.prj')
    4565 
    4566     def NOT_test_dem2pts_bounding_box_Nullvalues(self):
    4567         """Test conversion from dem in ascii format to native NetCDF xya format
    4568         """
    4569 
    4570         import time, os
    4571         from Numeric import array, zeros, allclose, Float, concatenate
    4572         from Scientific.IO.NetCDF import NetCDFFile
    4573 
    4574         #Write test asc file
    4575         root = 'demtest'
    4576 
    4577         filename = root+'.asc'
    4578         fid = open(filename, 'w')
    4579         fid.write("""ncols         5
    4580 nrows         6
    4581 xllcorner     2000.5
    4582 yllcorner     3000.5
    4583 cellsize      25
    4584 NODATA_value  -9999
    4585 """)
    4586         #Create linear function
    4587 
    4588         ref_points = []
    4589         ref_elevation = []
    4590         new_ref_pts1 = []
    4591         new_ref_elev1 = []
    4592         NODATA_value = -9999
    4593         for i in range(6):
    4594             y = (6-i)*25.0
    4595             for j in range(5):
    4596                 x = j*25.0
    4597                 z = x+2*y
    4598                 if j == 4: z = NODATA_value # column
    4599                 if i == 2 and j == 2: z = NODATA_value # random
    4600                 if i == 5 and j == 1: z = NODATA_value
    4601                 if i == 1: z = NODATA_value # row
    4602                 if i == 3 and j == 1: z = NODATA_value # two pts/row
    4603                 if i == 3 and j == 3: z = NODATA_value
    4604 
    4605                 if z <> NODATA_value:
    4606                     new_ref_elev1.append(z)
    4607                     new_ref_pts1.append( [x,y] )
    4608 
    4609                 ref_points.append( [x,y] )
    4610                 ref_elevation.append(z)
    4611                 fid.write('%f ' %z)
    4612             fid.write('\n')
    4613 
    4614         fid.close()
    4615 
    4616         #Write prj file with metadata
    4617         metafilename = root+'.prj'
    4618         fid = open(metafilename, 'w')
    4619 
    4620 
    4621         fid.write("""Projection UTM
    4622 Zone 56
    4623 Datum WGS84
    4624 Zunits NO
    4625 Units METERS
    4626 Spheroid WGS84
    4627 Xshift 0.0000000000
    4628 Yshift 10000000.0000000000
    4629 Parameters
    4630 """)
    4631         fid.close()
    4632 
    4633         #Convert to NetCDF pts
    4634         convert_dem_from_ascii2netcdf(root)
    4635         dem2pts(root, easting_min=2010.0, easting_max=2110.0,
    4636                 northing_min=3035.0, northing_max=3125.5)
    4637 
    4638         #Check contents
    4639         #Get NetCDF
    4640         fid = NetCDFFile(root+'.pts', 'r')
    4641 
    4642         # Get the variables
    4643         #print fid.variables.keys()
    4644         points = fid.variables['points']
    4645         elevation = fid.variables['elevation']
    4646 
    4647         #Check values
    4648         assert fid.xllcorner[0] == 2010.0
    4649         assert fid.yllcorner[0] == 3035.0
    4650 
    4651         #create new reference points
    4652         ref_points = []
    4653         ref_elevation = []
    4654         new_ref_pts2 = []
    4655         new_ref_elev2 = []
    4656         for i in range(4):
    4657             y = (4-i)*25.0 + 25.0
    4658             y_new = y + 3000.5 - 3035.0
    4659             for j in range(4):
    4660                 x = j*25.0 + 25.0
    4661                 x_new = x + 2000.5 - 2010.0
    4662                 z = x+2*y
    4663 
    4664                 if j == 3: z = NODATA_value # column
    4665                 if i == 1 and j == 1: z = NODATA_value # random
    4666                 if i == 4 and j == 0: z = NODATA_value
    4667                 if i == 0: z = NODATA_value # row
    4668                 if i == 2 and j == 0: z = NODATA_value # two pts/row
    4669                 if i == 2 and j == 2: z = NODATA_value
    4670 
    4671                 if z <> NODATA_value:
    4672                     new_ref_elev2.append(z)
    4673                     new_ref_pts2.append( [x_new,y_new] )
    4674 
    4675 
    4676                 ref_points.append( [x_new,y_new] )
    4677                 ref_elevation.append(z)
    4678 
    4679         #print points[:]
    4680         #print ref_points
    4681         #assert allclose(points, ref_points)
    4682 
    4683         #print attributes[:]
    4684         #print ref_elevation
    4685         #assert allclose(elevation, ref_elevation)
    4686 
    4687 
    4688         assert len(points) == len(new_ref_pts2), 'length of returned points not correct'
    4689         assert allclose(points, new_ref_pts2), 'points do not align'
    4690 
    4691         #print 'elevation', elevation[:]
    4692         assert len(elevation) == len(new_ref_elev2), 'length of returned elevation not correct'
    4693         assert allclose(elevation, new_ref_elev2), 'elevations do not align'
    4694         #Cleanup
    4695         fid.close()
    4696 
    4697 
    4698         os.remove(root + '.pts')
    4699         os.remove(root + '.dem')
    4700         os.remove(root + '.asc')
    4701         os.remove(root + '.prj')
    4702 
    47034271
    47044272########## testing nbed class ##################

  • anuga_core/source/obsolete_code/data_manager_obsolete_stuff.py

    r3565 r3761  
    391391    fid.close()
    392392
     393
     394
     395
     396    def NOT_test_dem2pts(self):
     397        """Test conversion from dem in ascii format to native NetCDF xya format
     398        """
     399
     400        import time, os
     401        from Numeric import array, zeros, allclose, Float, concatenate
     402        from Scientific.IO.NetCDF import NetCDFFile
     403
     404        #Write test asc file
     405        root = 'demtest'
     406
     407        filename = root+'.asc'
     408        fid = open(filename, 'w')
     409        fid.write("""ncols         5
     410nrows         6
     411xllcorner     2000.5
     412yllcorner     3000.5
     413cellsize      25
     414NODATA_value  -9999
     415""")
     416        #Create linear function
     417
     418        ref_points = []
     419        ref_elevation = []
     420        for i in range(6):
     421            y = (6-i)*25.0
     422            for j in range(5):
     423                x = j*25.0
     424                z = x+2*y
     425
     426                ref_points.append( [x,y] )
     427                ref_elevation.append(z)
     428                fid.write('%f ' %z)
     429            fid.write('\n')
     430
     431        fid.close()
     432
     433        #Write prj file with metadata
     434        metafilename = root+'.prj'
     435        fid = open(metafilename, 'w')
     436
     437
     438        fid.write("""Projection UTM
     439Zone 56
     440Datum WGS84
     441Zunits NO
     442Units METERS
     443Spheroid WGS84
     444Xshift 0.0000000000
     445Yshift 10000000.0000000000
     446Parameters
     447""")
     448        fid.close()
     449
     450        #Convert to NetCDF pts
     451        convert_dem_from_ascii2netcdf(root)
     452        dem2pts(root)
     453
     454        #Check contents
     455        #Get NetCDF
     456        fid = NetCDFFile(root+'.pts', 'r')
     457
     458        # Get the variables
     459        #print fid.variables.keys()
     460        points = fid.variables['points']
     461        elevation = fid.variables['elevation']
     462
     463        #Check values
     464
     465        #print points[:]
     466        #print ref_points
     467        assert allclose(points, ref_points)
     468
     469        #print attributes[:]
     470        #print ref_elevation
     471        assert allclose(elevation, ref_elevation)
     472
     473        #Cleanup
     474        fid.close()
     475
     476
     477        os.remove(root + '.pts')
     478        os.remove(root + '.dem')
     479        os.remove(root + '.asc')
     480        os.remove(root + '.prj')
     481
     482
     483
     484    def NOT_test_dem2pts_bounding_box(self):
     485        """Test conversion from dem in ascii format to native NetCDF xya format
     486        """
     487
     488        import time, os
     489        from Numeric import array, zeros, allclose, Float, concatenate
     490        from Scientific.IO.NetCDF import NetCDFFile
     491
     492        #Write test asc file
     493        root = 'demtest'
     494
     495        filename = root+'.asc'
     496        fid = open(filename, 'w')
     497        fid.write("""ncols         5
     498nrows         6
     499xllcorner     2000.5
     500yllcorner     3000.5
     501cellsize      25
     502NODATA_value  -9999
     503""")
     504        #Create linear function
     505
     506        ref_points = []
     507        ref_elevation = []
     508        for i in range(6):
     509            y = (6-i)*25.0
     510            for j in range(5):
     511                x = j*25.0
     512                z = x+2*y
     513
     514                ref_points.append( [x,y] )
     515                ref_elevation.append(z)
     516                fid.write('%f ' %z)
     517            fid.write('\n')
     518
     519        fid.close()
     520
     521        #Write prj file with metadata
     522        metafilename = root+'.prj'
     523        fid = open(metafilename, 'w')
     524
     525
     526        fid.write("""Projection UTM
     527Zone 56
     528Datum WGS84
     529Zunits NO
     530Units METERS
     531Spheroid WGS84
     532Xshift 0.0000000000
     533Yshift 10000000.0000000000
     534Parameters
     535""")
     536        fid.close()
     537
     538        #Convert to NetCDF pts
     539        convert_dem_from_ascii2netcdf(root)
     540        dem2pts(root, easting_min=2010.0, easting_max=2110.0,
     541                northing_min=3035.0, northing_max=3125.5)
     542
     543        #Check contents
     544        #Get NetCDF
     545        fid = NetCDFFile(root+'.pts', 'r')
     546
     547        # Get the variables
     548        #print fid.variables.keys()
     549        points = fid.variables['points']
     550        elevation = fid.variables['elevation']
     551
     552        #Check values
     553        assert fid.xllcorner[0] == 2010.0
     554        assert fid.yllcorner[0] == 3035.0
     555
     556        #create new reference points
     557        ref_points = []
     558        ref_elevation = []
     559        for i in range(4):
     560            y = (4-i)*25.0 + 25.0
     561            y_new = y + 3000.5 - 3035.0
     562            for j in range(4):
     563                x = j*25.0 + 25.0
     564                x_new = x + 2000.5 - 2010.0
     565                z = x+2*y
     566
     567                ref_points.append( [x_new,y_new] )
     568                ref_elevation.append(z)
     569
     570        #print points[:]
     571        #print ref_points
     572        assert allclose(points, ref_points)
     573
     574        #print attributes[:]
     575        #print ref_elevation
     576        assert allclose(elevation, ref_elevation)
     577
     578        #Cleanup
     579        fid.close()
     580
     581
     582        os.remove(root + '.pts')
     583        os.remove(root + '.dem')
     584        os.remove(root + '.asc')
     585        os.remove(root + '.prj')
     586
     587
     588
     589    def NOT_test_dem2pts_remove_Nullvalues(self):
     590        """Test conversion from dem in ascii format to native NetCDF xya format
     591        """
     592
     593        import time, os
     594        from Numeric import array, zeros, allclose, Float, concatenate
     595        from Scientific.IO.NetCDF import NetCDFFile
     596
     597        #Write test asc file
     598        root = 'demtest'
     599
     600        filename = root+'.asc'
     601        fid = open(filename, 'w')
     602        fid.write("""ncols         5
     603nrows         6
     604xllcorner     2000.5
     605yllcorner     3000.5
     606cellsize      25
     607NODATA_value  -9999
     608""")
     609        #Create linear function
     610        # ref_ will write all the values
     611        # new_ref_ will write the values except for NODATA_values
     612        ref_points = []
     613        ref_elevation = []
     614        new_ref_pts = []
     615        new_ref_elev = []
     616        NODATA_value = -9999
     617        for i in range(6):
     618            y = (6-i)*25.0
     619            for j in range(5):
     620                x = j*25.0
     621                z = x+2*y
     622                if j == 4: z = NODATA_value # column
     623                if i == 2 and j == 2: z = NODATA_value # random
     624                if i == 5 and j == 1: z = NODATA_value
     625                if i == 1: z = NODATA_value # row
     626                if i == 3 and j == 1: z = NODATA_value # two pts/row
     627                if i == 3 and j == 3: z = NODATA_value
     628
     629
     630                if z <> NODATA_value:
     631                    new_ref_elev.append(z)
     632                    new_ref_pts.append( [x,y] )
     633
     634                ref_points.append( [x,y] )
     635                ref_elevation.append(z)
     636
     637                fid.write('%f ' %z)
     638            fid.write('\n')
     639
     640        fid.close()
     641
     642
     643        #Write prj file with metadata
     644        metafilename = root+'.prj'
     645        fid = open(metafilename, 'w')
     646
     647
     648        fid.write("""Projection UTM
     649Zone 56
     650Datum WGS84
     651Zunits NO
     652Units METERS
     653Spheroid WGS84
     654Xshift 0.0000000000
     655Yshift 10000000.0000000000
     656Parameters
     657""")
     658        fid.close()
     659
     660        #Convert to NetCDF pts
     661        convert_dem_from_ascii2netcdf(root)
     662        dem2pts(root)
     663
     664        #Check contents
     665        #Get NetCDF
     666        fid = NetCDFFile(root+'.pts', 'r')
     667
     668        # Get the variables
     669        #print fid.variables.keys()
     670        points = fid.variables['points']
     671        elevation = fid.variables['elevation']
     672
     673        #Check values
     674        #print 'points', points[:]
     675        assert len(points) == len(new_ref_pts), 'length of returned points not correct'
     676        assert allclose(points, new_ref_pts), 'points do not align'
     677
     678        #print 'elevation', elevation[:]
     679        assert len(elevation) == len(new_ref_elev), 'length of returned elevation not correct'
     680        assert allclose(elevation, new_ref_elev), 'elevations do not align'
     681
     682        #Cleanup
     683        fid.close()
     684
     685
     686        os.remove(root + '.pts')
     687        os.remove(root + '.dem')
     688        os.remove(root + '.asc')
     689        os.remove(root + '.prj')
     690
     691    def NOT_test_dem2pts_bounding_box_Nullvalues(self):
     692        """Test conversion from dem in ascii format to native NetCDF xya format
     693        """
     694
     695        import time, os
     696        from Numeric import array, zeros, allclose, Float, concatenate
     697        from Scientific.IO.NetCDF import NetCDFFile
     698
     699        #Write test asc file
     700        root = 'demtest'
     701
     702        filename = root+'.asc'
     703        fid = open(filename, 'w')
     704        fid.write("""ncols         5
     705nrows         6
     706xllcorner     2000.5
     707yllcorner     3000.5
     708cellsize      25
     709NODATA_value  -9999
     710""")
     711        #Create linear function
     712
     713        ref_points = []
     714        ref_elevation = []
     715        new_ref_pts1 = []
     716        new_ref_elev1 = []
     717        NODATA_value = -9999
     718        for i in range(6):
     719            y = (6-i)*25.0
     720            for j in range(5):
     721                x = j*25.0
     722                z = x+2*y
     723                if j == 4: z = NODATA_value # column
     724                if i == 2 and j == 2: z = NODATA_value # random
     725                if i == 5 and j == 1: z = NODATA_value
     726                if i == 1: z = NODATA_value # row
     727                if i == 3 and j == 1: z = NODATA_value # two pts/row
     728                if i == 3 and j == 3: z = NODATA_value
     729
     730                if z <> NODATA_value:
     731                    new_ref_elev1.append(z)
     732                    new_ref_pts1.append( [x,y] )
     733
     734                ref_points.append( [x,y] )
     735                ref_elevation.append(z)
     736                fid.write('%f ' %z)
     737            fid.write('\n')
     738
     739        fid.close()
     740
     741        #Write prj file with metadata
     742        metafilename = root+'.prj'
     743        fid = open(metafilename, 'w')
     744
     745
     746        fid.write("""Projection UTM
     747Zone 56
     748Datum WGS84
     749Zunits NO
     750Units METERS
     751Spheroid WGS84
     752Xshift 0.0000000000
     753Yshift 10000000.0000000000
     754Parameters
     755""")
     756        fid.close()
     757
     758        #Convert to NetCDF pts
     759        convert_dem_from_ascii2netcdf(root)
     760        dem2pts(root, easting_min=2010.0, easting_max=2110.0,
     761                northing_min=3035.0, northing_max=3125.5)
     762
     763        #Check contents
     764        #Get NetCDF
     765        fid = NetCDFFile(root+'.pts', 'r')
     766
     767        # Get the variables
     768        #print fid.variables.keys()
     769        points = fid.variables['points']
     770        elevation = fid.variables['elevation']
     771
     772        #Check values
     773        assert fid.xllcorner[0] == 2010.0
     774        assert fid.yllcorner[0] == 3035.0
     775
     776        #create new reference points
     777        ref_points = []
     778        ref_elevation = []
     779        new_ref_pts2 = []
     780        new_ref_elev2 = []
     781        for i in range(4):
     782            y = (4-i)*25.0 + 25.0
     783            y_new = y + 3000.5 - 3035.0
     784            for j in range(4):
     785                x = j*25.0 + 25.0
     786                x_new = x + 2000.5 - 2010.0
     787                z = x+2*y
     788
     789                if j == 3: z = NODATA_value # column
     790                if i == 1 and j == 1: z = NODATA_value # random
     791                if i == 4 and j == 0: z = NODATA_value
     792                if i == 0: z = NODATA_value # row
     793                if i == 2 and j == 0: z = NODATA_value # two pts/row
     794                if i == 2 and j == 2: z = NODATA_value
     795
     796                if z <> NODATA_value:
     797                    new_ref_elev2.append(z)
     798                    new_ref_pts2.append( [x_new,y_new] )
     799
     800
     801                ref_points.append( [x_new,y_new] )
     802                ref_elevation.append(z)
     803
     804        #print points[:]
     805        #print ref_points
     806        #assert allclose(points, ref_points)
     807
     808        #print attributes[:]
     809        #print ref_elevation
     810        #assert allclose(elevation, ref_elevation)
     811
     812
     813        assert len(points) == len(new_ref_pts2), 'length of returned points not correct'
     814        assert allclose(points, new_ref_pts2), 'points do not align'
     815
     816        #print 'elevation', elevation[:]
     817        assert len(elevation) == len(new_ref_elev2), 'length of returned elevation not correct'
     818        assert allclose(elevation, new_ref_elev2), 'elevations do not align'
     819        #Cleanup
     820        fid.close()
     821
     822
     823        os.remove(root + '.pts')
     824        os.remove(root + '.dem')
     825        os.remove(root + '.asc')
     826        os.remove(root + '.prj')
     827
    393828#********************
    394829#*** END OF OBSOLETE FUNCTIONS
Note: See TracChangeset for help on using the changeset viewer.