- Timestamp:
- Jun 22, 2010, 12:04:32 PM (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/anuga_core/source/anuga/file_conversion/test_sww2dem.py
r7841 r7866 1520 1520 os.remove(ascfile) 1521 1521 os.remove(swwfile2) 1522 1523 1524 def test_export_grid(self): 1525 """ 1526 test_export_grid(self): 1527 Test that sww information can be converted correctly to asc/prj 1528 format readable by e.g. ArcView 1529 """ 1530 1531 import time, os 1532 from Scientific.IO.NetCDF import NetCDFFile 1533 1534 try: 1535 os.remove('teg*.sww') 1536 except: 1537 pass 1538 1539 #Setup 1540 self.domain.set_name('teg') 1541 1542 prjfile = self.domain.get_name() + '_elevation.prj' 1543 ascfile = self.domain.get_name() + '_elevation.asc' 1544 swwfile = self.domain.get_name() + '.sww' 1545 1546 self.domain.set_datadir('.') 1547 self.domain.smooth = True 1548 self.domain.set_quantity('elevation', lambda x,y: -x-y) 1549 self.domain.set_quantity('stage', 1.0) 1550 1551 self.domain.geo_reference = Geo_reference(56,308500,6189000) 1552 1553 sww = SWW_file(self.domain) 1554 sww.store_connectivity() 1555 sww.store_timestep() 1556 self.domain.evolve_to_end(finaltime = 0.01) 1557 sww.store_timestep() 1558 1559 cellsize = 0.25 1560 #Check contents 1561 #Get NetCDF 1562 1563 fid = NetCDFFile(sww.filename, netcdf_mode_r) 1564 1565 # Get the variables 1566 x = fid.variables['x'][:] 1567 y = fid.variables['y'][:] 1568 z = fid.variables['elevation'][:] 1569 time = fid.variables['time'][:] 1570 stage = fid.variables['stage'][:] 1571 1572 fid.close() 1573 1574 #Export to ascii/prj files 1575 sww2dem_batch(self.domain.get_name(), 1576 quantities = 'elevation', 1577 cellsize = cellsize, 1578 verbose = self.verbose, 1579 format = 'asc') 1580 1581 #Check asc file 1582 ascid = open(ascfile) 1583 lines = ascid.readlines() 1584 ascid.close() 1585 1586 L = lines[2].strip().split() 1587 assert L[0].strip().lower() == 'xllcorner' 1588 assert num.allclose(float(L[1].strip().lower()), 308500) 1589 1590 L = lines[3].strip().split() 1591 assert L[0].strip().lower() == 'yllcorner' 1592 assert num.allclose(float(L[1].strip().lower()), 6189000) 1593 1594 #Check grid values 1595 for j in range(5): 1596 L = lines[6+j].strip().split() 1597 y = (4-j) * cellsize 1598 for i in range(5): 1599 assert num.allclose(float(L[i]), -i*cellsize - y) 1600 1601 #Cleanup 1602 os.remove(prjfile) 1603 os.remove(ascfile) 1604 os.remove(swwfile) 1605 1606 def test_export_gridII(self): 1607 """ 1608 test_export_gridII(self): 1609 Test that sww information can be converted correctly to asc/prj 1610 format readable by e.g. ArcView 1611 """ 1612 1613 import time, os 1614 from Scientific.IO.NetCDF import NetCDFFile 1615 1616 try: 1617 os.remove('teg*.sww') 1618 except: 1619 pass 1620 1621 #Setup 1622 self.domain.set_name('tegII') 1623 1624 swwfile = self.domain.get_name() + '.sww' 1625 1626 self.domain.set_datadir('.') 1627 self.domain.smooth = True 1628 self.domain.set_quantity('elevation', lambda x,y: -x-y) 1629 self.domain.set_quantity('stage', 1.0) 1630 1631 self.domain.geo_reference = Geo_reference(56,308500,6189000) 1632 1633 sww = SWW_file(self.domain) 1634 sww.store_connectivity() 1635 sww.store_timestep() 1636 self.domain.evolve_to_end(finaltime = 0.01) 1637 sww.store_timestep() 1638 1639 cellsize = 0.25 1640 #Check contents 1641 #Get NetCDF 1642 1643 fid = NetCDFFile(sww.filename, netcdf_mode_r) 1644 1645 # Get the variables 1646 x = fid.variables['x'][:] 1647 y = fid.variables['y'][:] 1648 z = fid.variables['elevation'][:] 1649 time = fid.variables['time'][:] 1650 stage = fid.variables['stage'][:] 1651 xmomentum = fid.variables['xmomentum'][:] 1652 ymomentum = fid.variables['ymomentum'][:] 1653 1654 #print 'stage', stage 1655 #print 'xmom', xmomentum 1656 #print 'ymom', ymomentum 1657 1658 fid.close() 1659 1660 #Export to ascii/prj files 1661 if True: 1662 sww2dem_batch(self.domain.get_name(), 1663 quantities = ['elevation', 'depth'], 1664 cellsize = cellsize, 1665 verbose = self.verbose, 1666 format = 'asc') 1667 1668 else: 1669 sww2dem_batch(self.domain.get_name(), 1670 quantities = ['depth'], 1671 cellsize = cellsize, 1672 verbose = self.verbose, 1673 format = 'asc') 1674 1675 1676 export_grid(self.domain.get_name(), 1677 quantities = ['elevation'], 1678 cellsize = cellsize, 1679 verbose = self.verbose, 1680 format = 'asc') 1681 1682 prjfile = self.domain.get_name() + '_elevation.prj' 1683 ascfile = self.domain.get_name() + '_elevation.asc' 1684 1685 #Check asc file 1686 ascid = open(ascfile) 1687 lines = ascid.readlines() 1688 ascid.close() 1689 1690 L = lines[2].strip().split() 1691 assert L[0].strip().lower() == 'xllcorner' 1692 assert num.allclose(float(L[1].strip().lower()), 308500) 1693 1694 L = lines[3].strip().split() 1695 assert L[0].strip().lower() == 'yllcorner' 1696 assert num.allclose(float(L[1].strip().lower()), 6189000) 1697 1698 #print "ascfile", ascfile 1699 #Check grid values 1700 for j in range(5): 1701 L = lines[6+j].strip().split() 1702 y = (4-j) * cellsize 1703 for i in range(5): 1704 #print " -i*cellsize - y", -i*cellsize - y 1705 #print "float(L[i])", float(L[i]) 1706 assert num.allclose(float(L[i]), -i*cellsize - y) 1707 1708 #Cleanup 1709 os.remove(prjfile) 1710 os.remove(ascfile) 1711 1712 #Check asc file 1713 ascfile = self.domain.get_name() + '_depth.asc' 1714 prjfile = self.domain.get_name() + '_depth.prj' 1715 ascid = open(ascfile) 1716 lines = ascid.readlines() 1717 ascid.close() 1718 1719 L = lines[2].strip().split() 1720 assert L[0].strip().lower() == 'xllcorner' 1721 assert num.allclose(float(L[1].strip().lower()), 308500) 1722 1723 L = lines[3].strip().split() 1724 assert L[0].strip().lower() == 'yllcorner' 1725 assert num.allclose(float(L[1].strip().lower()), 6189000) 1726 1727 #Check grid values 1728 for j in range(5): 1729 L = lines[6+j].strip().split() 1730 y = (4-j) * cellsize 1731 for i in range(5): 1732 #print " -i*cellsize - y", -i*cellsize - y 1733 #print "float(L[i])", float(L[i]) 1734 assert num.allclose(float(L[i]), 1 - (-i*cellsize - y)) 1735 1736 #Cleanup 1737 os.remove(prjfile) 1738 os.remove(ascfile) 1739 os.remove(swwfile) 1740 1741 1742 def test_export_gridIII(self): 1743 """ 1744 test_export_gridIII 1745 Test that sww information can be converted correctly to asc/prj 1746 format readable by e.g. ArcView 1747 """ 1748 1749 import time, os 1750 from Scientific.IO.NetCDF import NetCDFFile 1751 1752 try: 1753 os.remove('teg*.sww') 1754 except: 1755 pass 1756 1757 #Setup 1758 1759 self.domain.set_name('tegIII') 1760 1761 swwfile = self.domain.get_name() + '.sww' 1762 1763 self.domain.set_datadir('.') 1764 self.domain.format = 'sww' 1765 self.domain.smooth = True 1766 self.domain.set_quantity('elevation', lambda x,y: -x-y) 1767 self.domain.set_quantity('stage', 1.0) 1768 1769 self.domain.geo_reference = Geo_reference(56,308500,6189000) 1770 1771 sww = SWW_file(self.domain) 1772 sww.store_connectivity() 1773 sww.store_timestep() #'stage') 1774 self.domain.evolve_to_end(finaltime = 0.01) 1775 sww.store_timestep() #'stage') 1776 1777 cellsize = 0.25 1778 #Check contents 1779 #Get NetCDF 1780 1781 fid = NetCDFFile(sww.filename, netcdf_mode_r) 1782 1783 # Get the variables 1784 x = fid.variables['x'][:] 1785 y = fid.variables['y'][:] 1786 z = fid.variables['elevation'][:] 1787 time = fid.variables['time'][:] 1788 stage = fid.variables['stage'][:] 1789 1790 fid.close() 1791 1792 #Export to ascii/prj files 1793 extra_name_out = 'yeah' 1794 if True: 1795 sww2dem_batch(self.domain.get_name(), 1796 quantities = ['elevation', 'depth'], 1797 extra_name_out = extra_name_out, 1798 cellsize = cellsize, 1799 verbose = self.verbose, 1800 format = 'asc') 1801 1802 else: 1803 sww2dem_batch(self.domain.get_name(), 1804 quantities = ['depth'], 1805 cellsize = cellsize, 1806 verbose = self.verbose, 1807 format = 'asc') 1808 1809 1810 sww2dem_batch(self.domain.get_name(), 1811 quantities = ['elevation'], 1812 cellsize = cellsize, 1813 verbose = self.verbose, 1814 format = 'asc') 1815 1816 prjfile = self.domain.get_name() + '_elevation_yeah.prj' 1817 ascfile = self.domain.get_name() + '_elevation_yeah.asc' 1818 1819 #Check asc file 1820 ascid = open(ascfile) 1821 lines = ascid.readlines() 1822 ascid.close() 1823 1824 L = lines[2].strip().split() 1825 assert L[0].strip().lower() == 'xllcorner' 1826 assert num.allclose(float(L[1].strip().lower()), 308500) 1827 1828 L = lines[3].strip().split() 1829 assert L[0].strip().lower() == 'yllcorner' 1830 assert num.allclose(float(L[1].strip().lower()), 6189000) 1831 1832 #print "ascfile", ascfile 1833 #Check grid values 1834 for j in range(5): 1835 L = lines[6+j].strip().split() 1836 y = (4-j) * cellsize 1837 for i in range(5): 1838 #print " -i*cellsize - y", -i*cellsize - y 1839 #print "float(L[i])", float(L[i]) 1840 assert num.allclose(float(L[i]), -i*cellsize - y) 1841 1842 #Cleanup 1843 os.remove(prjfile) 1844 os.remove(ascfile) 1845 1846 #Check asc file 1847 ascfile = self.domain.get_name() + '_depth_yeah.asc' 1848 prjfile = self.domain.get_name() + '_depth_yeah.prj' 1849 ascid = open(ascfile) 1850 lines = ascid.readlines() 1851 ascid.close() 1852 1853 L = lines[2].strip().split() 1854 assert L[0].strip().lower() == 'xllcorner' 1855 assert num.allclose(float(L[1].strip().lower()), 308500) 1856 1857 L = lines[3].strip().split() 1858 assert L[0].strip().lower() == 'yllcorner' 1859 assert num.allclose(float(L[1].strip().lower()), 6189000) 1860 1861 #Check grid values 1862 for j in range(5): 1863 L = lines[6+j].strip().split() 1864 y = (4-j) * cellsize 1865 for i in range(5): 1866 assert num.allclose(float(L[i]), 1 - (-i*cellsize - y)) 1867 1868 #Cleanup 1869 os.remove(prjfile) 1870 os.remove(ascfile) 1871 os.remove(swwfile) 1872 1873 def test_export_grid_bad(self): 1874 """Test that sww information can be converted correctly to asc/prj 1875 format readable by e.g. ArcView 1876 """ 1877 1878 try: 1879 sww2dem_batch('a_small_round-egg', 1880 quantities = ['elevation', 'depth'], 1881 cellsize = 99, 1882 verbose = self.verbose, 1883 format = 'asc') 1884 except IOError: 1885 pass 1886 else: 1887 self.failUnless(0 ==1, 'Bad input did not throw exception error!') 1522 1888 1523 1889
Note: See TracChangeset
for help on using the changeset viewer.