Changeset 2718


Ignore:
Timestamp:
Apr 17, 2006, 10:47:48 PM (17 years ago)
Author:
steve
Message:
 
Location:
inundation/pyvolution
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • inundation/pyvolution/data_manager.py

    r2648 r2718  
    4848
    4949"""
     50
     51import exceptions
     52class DataMissingValuesError(exceptions.Exception): pass
     53class DataFileNotOpenError(exceptions.Exception): pass
     54class DataTimeError(exceptions.Exception): pass
     55class DataDomainError(exceptions.Exception): pass
     56
     57
     58
    5059import os
    5160
    … …  
    156165
    157166    if size is not None:
    158         FN += '_size%d' %size
     167        FN += '_size%d' %size
    159168
    160169    if time is not None:
    … …  
    199208
    200209        self.filename = create_filename(domain.get_datadir(),
    201                         domain.get_name(), extension)
     210            domain.get_name(), extension)
    202211
    203212        #print 'F', self.filename
    … …  
    259268            if hasattr(domain, 'texture'):
    260269                fid.texture = domain.texture
    261             #else:
     270        #else:
    262271            #    fid.texture = 'None'
    263272
    … …  
    337346        Q = domain.quantities['elevation']
    338347        X,Y,Z,V = Q.get_vertex_values(xy=True,
    339                                       precision = self.precision)
     348                      precision = self.precision)
    340349
    341350
    … …  
    387396        if not file_open:
    388397            msg = 'File %s could not be opened for append' %self.filename
    389             raise msg
     398            raise DataFileNotOpenError, msg
    390399
    391400
    … …  
    473482
    474483                #As in....
    475                 #eval( name + '[i,:] = A.astype(self.precision)' )
    476                 #FIXME: But we need a UNIT test for that before refactoring
     484            #eval( name + '[i,:] = A.astype(self.precision)' )
     485            #FIXME: But we need a UNIT test for that before refactoring
    477486
    478487
    … …  
    566575        Q = domain.quantities['elevation']
    567576        X,Y,Z,V = Q.get_vertex_values(xy=True,
    568                                       precision = self.precision)
     577                      precision = self.precision)
    569578
    570579
    … …  
    596605                #In that case, wait a while and try again
    597606                msg = 'Warning (store_timestep): File %s could not be opened'\
    598                           %self.filename
     607                  %self.filename
    599608                msg += ' - trying again'
    600609                print msg
    … …  
    604613                file_open = True
    605614
    606             if not file_open:
    607                 msg = 'File %s could not be opened for append' %self.filename
    608             raise msg
     615        if not file_open:
     616            msg = 'File %s could not be opened for append' %self.filename
     617            raise DataFileNotOPenError, msg
    609618
    610619
    … …  
    622631        Q = domain.quantities[name]
    623632        A,V = Q.get_vertex_values(xy=False,
    624                                   precision = self.precision)
     633                  precision = self.precision)
    625634
    626635        stage[i,:] = A.astype(self.precision)
    … …  
    866875            t = float(filename[i0:i1])
    867876        else:
    868             raise 'Hmmmm'
     877            raise DataTimeError, 'Hmmmm'
    869878
    870879
    … …  
    15241533            origin = None,
    15251534            datum = 'WGS84',
    1526             format = 'ers'):
     1535        format = 'ers'):
    15271536
    15281537    """Read SWW file and convert to Digitial Elevation model format (.asc or .ers)
    … …  
    17351744    for i in xrange(nrows):
    17361745        if format.lower() == 'asc':
    1737             yg = i*cellsize
    1738         else:
    1739             #this will flip the order of the y values for ers
     1746        yg = i*cellsize
     1747    else:
     1748        #this will flip the order of the y values for ers
    17401749            yg = (nrows-i)*cellsize
    17411750
    … …  
    17821791
    17831792    if format.lower() == 'ers':
    1784         # setup ERS header information
    1785         grid_values = reshape(grid_values,(nrows, ncols))
    1786         header = {}
    1787         header['datum'] = '"' + datum + '"'
    1788         # FIXME The use of hardwired UTM and zone number needs to be made optional
    1789         # FIXME Also need an automatic test for coordinate type (i.e. EN or LL)
    1790         header['projection'] = '"UTM-' + str(zone) + '"'
    1791         header['coordinatetype'] = 'EN'
    1792         if header['coordinatetype'] == 'LL':
    1793             header['longitude'] = str(newxllcorner)
    1794             header['latitude'] = str(newyllcorner)
    1795         elif header['coordinatetype'] == 'EN':
    1796             header['eastings'] = str(newxllcorner)
    1797             header['northings'] = str(newyllcorner)
    1798         header['nullcellvalue'] = str(NODATA_value)
    1799         header['xdimension'] = str(cellsize)
    1800         header['ydimension'] = str(cellsize)
    1801         header['value'] = '"' + quantity + '"'
    1802         #header['celltype'] = 'IEEE8ByteReal'  #FIXME: Breaks unit test
    1803 
    1804 
    1805         #Write
    1806         if verbose: print 'Writing %s' %demfile
    1807         import ermapper_grids
    1808         ermapper_grids.write_ermapper_grid(demfile, grid_values, header)
    1809 
    1810         fid.close()
     1793        # setup ERS header information
     1794        grid_values = reshape(grid_values,(nrows, ncols))
     1795        header = {}
     1796        header['datum'] = '"' + datum + '"'
     1797        # FIXME The use of hardwired UTM and zone number needs to be made optional
     1798        # FIXME Also need an automatic test for coordinate type (i.e. EN or LL)
     1799        header['projection'] = '"UTM-' + str(zone) + '"'
     1800        header['coordinatetype'] = 'EN'
     1801        if header['coordinatetype'] == 'LL':
     1802            header['longitude'] = str(newxllcorner)
     1803            header['latitude'] = str(newyllcorner)
     1804        elif header['coordinatetype'] == 'EN':
     1805            header['eastings'] = str(newxllcorner)
     1806            header['northings'] = str(newyllcorner)
     1807        header['nullcellvalue'] = str(NODATA_value)
     1808        header['xdimension'] = str(cellsize)
     1809        header['ydimension'] = str(cellsize)
     1810        header['value'] = '"' + quantity + '"'
     1811        #header['celltype'] = 'IEEE8ByteReal'  #FIXME: Breaks unit test
     1812
     1813
     1814        #Write
     1815        if verbose: print 'Writing %s' %demfile
     1816        import ermapper_grids
     1817        ermapper_grids.write_ermapper_grid(demfile, grid_values, header)
     1818
     1819        fid.close()
    18111820    else:
    18121821        #Write to Ascii format
    … …  
    18151824        prjfile = basename_out + '.prj'
    18161825
    1817         if verbose: print 'Writing %s' %prjfile
    1818         prjid = open(prjfile, 'w')
    1819         prjid.write('Projection    %s\n' %'UTM')
    1820         prjid.write('Zone          %d\n' %zone)
    1821         prjid.write('Datum         %s\n' %datum)
    1822         prjid.write('Zunits        NO\n')
    1823         prjid.write('Units         METERS\n')
    1824         prjid.write('Spheroid      %s\n' %datum)
    1825         prjid.write('Xshift        %d\n' %false_easting)
    1826         prjid.write('Yshift        %d\n' %false_northing)
    1827         prjid.write('Parameters\n')
    1828         prjid.close()
    1829 
    1830 
    1831 
    1832         if verbose: print 'Writing %s' %demfile
    1833 
    1834         ascid = open(demfile, 'w')
    1835 
    1836         ascid.write('ncols         %d\n' %ncols)
    1837         ascid.write('nrows         %d\n' %nrows)
    1838         ascid.write('xllcorner     %d\n' %newxllcorner)
    1839         ascid.write('yllcorner     %d\n' %newyllcorner)
    1840         ascid.write('cellsize      %f\n' %cellsize)
    1841         ascid.write('NODATA_value  %d\n' %NODATA_value)
    1842 
    1843 
    1844         #Get bounding polygon from mesh
    1845         #P = interp.mesh.get_boundary_polygon()
    1846         #inside_indices = inside_polygon(grid_points, P)
    1847 
    1848         for i in range(nrows):
    1849             if verbose and i%((nrows+10)/10)==0:
    1850                 print 'Doing row %d of %d' %(i, nrows)
     1826        if verbose: print 'Writing %s' %prjfile
     1827        prjid = open(prjfile, 'w')
     1828        prjid.write('Projection    %s\n' %'UTM')
     1829        prjid.write('Zone          %d\n' %zone)
     1830        prjid.write('Datum         %s\n' %datum)
     1831        prjid.write('Zunits        NO\n')
     1832        prjid.write('Units         METERS\n')
     1833        prjid.write('Spheroid      %s\n' %datum)
     1834        prjid.write('Xshift        %d\n' %false_easting)
     1835        prjid.write('Yshift        %d\n' %false_northing)
     1836        prjid.write('Parameters\n')
     1837        prjid.close()
     1838
     1839
     1840
     1841        if verbose: print 'Writing %s' %demfile
     1842
     1843        ascid = open(demfile, 'w')
     1844
     1845        ascid.write('ncols         %d\n' %ncols)
     1846        ascid.write('nrows         %d\n' %nrows)
     1847        ascid.write('xllcorner     %d\n' %newxllcorner)
     1848        ascid.write('yllcorner     %d\n' %newyllcorner)
     1849        ascid.write('cellsize      %f\n' %cellsize)
     1850        ascid.write('NODATA_value  %d\n' %NODATA_value)
     1851
     1852
     1853        #Get bounding polygon from mesh
     1854        #P = interp.mesh.get_boundary_polygon()
     1855        #inside_indices = inside_polygon(grid_points, P)
     1856
     1857        for i in range(nrows):
     1858            if verbose and i%((nrows+10)/10)==0:
     1859                print 'Doing row %d of %d' %(i, nrows)
    18511860
    18521861            base_index = (nrows-i-1)*ncols
    … …  
    18581867
    18591868            #print
    1860             #for j in range(ncols):
    1861             #    index = base_index+j#
     1869            #for j in range(ncols):
     1870            #    index = base_index+j#
    18621871            #    print grid_values[index],
    18631872            #    ascid.write('%f ' %grid_values[index])
    1864             #ascid.write('\n')
    1865 
    1866         #Close
    1867         ascid.close()
    1868         fid.close()
     1873            #ascid.write('\n')
     1874
     1875        #Close
     1876        ascid.close()
     1877        fid.close()
    18691878
    18701879#Backwards compatibility
    … …  
    18851894            verbose = verbose,
    18861895            origin = origin,
    1887             datum = 'WGS84',
    1888             format = 'asc')
     1896        datum = 'WGS84',
     1897        format = 'asc')
    18891898
    18901899def sww2ers(basename_in, basename_out = None,
    … …  
    19051914            verbose = verbose,
    19061915            origin = origin,
    1907             datum = datum,
    1908             format = 'ers')
     1916        datum = datum,
     1917        format = 'ers')
    19091918################################# END COMPATIBILITY ##############
    19101919
    … …  
    23102319            msg = 'NetCDFFile %s contains missing values'\
    23112320                  %(basename_in+'_ha.nc')
    2312             raise msg
     2321            raise DataMissingValuesError, msg
    23132322        else:
    23142323            amplitudes = amplitudes*(missing==0) + missing*NaN_filler
    … …  
    23192328            msg = 'NetCDFFile %s contains missing values'\
    23202329                  %(basename_in+'_ua.nc')
    2321             raise msg
     2330            raise DataMissingValuesError, msg
    23222331        else:
    23232332            uspeed = uspeed*(missing==0) + missing*NaN_filler
    … …  
    23282337            msg = 'NetCDFFile %s contains missing values'\
    23292338                  %(basename_in+'_va.nc')
    2330             raise msg
     2339            raise DataMissingValuesError, msg
    23312340        else:
    23322341            vspeed = vspeed*(missing==0) + missing*NaN_filler
    … …  
    23382347            msg = 'NetCDFFile %s contains missing values'\
    23392348                  %(basename_in+'_e.nc')
    2340             raise msg
     2349            raise DataMissingValuesError, msg
    23412350        else:
    23422351            elevations = elevations*(missing==0) + missing*NaN_filler
    … …  
    25162525        else:
    25172526            z = elevations
    2518         #FIXME: z should be obtained from MOST and passed in here
     2527    #FIXME: z should be obtained from MOST and passed in here
    25192528
    25202529    from Numeric import resize
    … …  
    26232632        msg += ' date-time with format %s.\n' %time_format
    26242633        msg += 'I got %s instead.' %fields[0]
    2625         raise msg
     2634        raise DataTimeError, msg
    26262635
    26272636
    … …  
    28282837        fid.close()
    28292838        msg = 'Domain could not be created: %s. Perhaps use "fail_if_NaN=False and NaN_filler = ..."' %e
    2830         raise msg
     2839        raise DataDomainError, msg
    28312840
    28322841    if not boundary is None:
    … …  
    28552864            if fail_if_NaN:
    28562865                msg = 'quantity "%s" contains no_data entry'%quantity
    2857                 raise msg
     2866                raise DataMissingValuesError, msg
    28582867            else:
    28592868                data = (X<>NaN)
    … …  
    28802889            if fail_if_NaN:
    28812890                msg = 'quantity "%s" contains no_data entry'%quantity
    2882                 raise msg
     2891                raise DataMissingValuesError, msg
    28832892            else:
    28842893                data = (X<>NaN)
    … …  
    29172926            %('FIXMEfilename', T[0], T[-1])
    29182927        msg += ' does not match model time: %s' %tau
    2919         if tau < time[0]: raise msg
    2920         if tau > time[-1]: raise msg
     2928        if tau < time[0]: raise DataTimeError, msg
     2929        if tau > time[-1]: raise DataTimeError, msg
    29212930        while tau > time[index]: index += 1
    29222931        while tau < time[index]: index -= 1
    … …  
    34433452                msg = 'File %s contains missing values'\
    34443453                      %(elevation_files[j])
    3445                 raise msg
     3454                raise DataMissingValuesError, msg
    34463455            else:
    34473456                elevation_grid = elevation_grid*(missing==0) + \

  • inundation/pyvolution/test_all.py

    r2516 r2718  
    9696
    9797    #print regressionTest()
    98     unittest.main(defaultTest='regressionTest')
     98    #unittest.main(defaultTest='regressionTest')
     99   
     100    suite = regressionTest()
     101    runner = unittest.TextTestRunner(verbosity=2)
     102    runner.run(suite)
Note: See TracChangeset for help on using the changeset viewer.