[7776] | 1 | """ This module is responsible for loading and saving NetCDF NC files |
---|
| 2 | """ |
---|
| 3 | |
---|
| 4 | import os, sys |
---|
| 5 | import csv |
---|
| 6 | import string |
---|
| 7 | import shutil |
---|
| 8 | from struct import unpack |
---|
| 9 | import array as p_array |
---|
| 10 | from os import sep, path, remove, mkdir, access, F_OK, W_OK, getcwd |
---|
| 11 | |
---|
| 12 | import numpy as num |
---|
| 13 | |
---|
| 14 | from Scientific.IO.NetCDF import NetCDFFile |
---|
| 15 | from os.path import exists, basename, join |
---|
| 16 | from os import getcwd |
---|
| 17 | |
---|
| 18 | from anuga.coordinate_transforms.redfearn import redfearn, \ |
---|
| 19 | convert_from_latlon_to_utm |
---|
| 20 | from anuga.coordinate_transforms.geo_reference import Geo_reference, \ |
---|
| 21 | write_NetCDF_georeference, ensure_geo_reference |
---|
| 22 | from anuga.geospatial_data.geospatial_data import Geospatial_data,\ |
---|
| 23 | ensure_absolute |
---|
| 24 | from anuga.config import minimum_storable_height as \ |
---|
| 25 | default_minimum_storable_height |
---|
| 26 | from anuga.config import netcdf_mode_r, netcdf_mode_w, netcdf_mode_a |
---|
| 27 | from anuga.config import netcdf_float, netcdf_float32, netcdf_int |
---|
| 28 | from anuga.config import max_float |
---|
| 29 | from anuga.utilities.numerical_tools import ensure_numeric, mean |
---|
| 30 | from anuga.caching.caching import myhash |
---|
| 31 | from anuga.shallow_water.shallow_water_domain import Domain |
---|
| 32 | from anuga.abstract_2d_finite_volumes.pmesh2domain import \ |
---|
| 33 | pmesh_to_domain_instance |
---|
| 34 | from anuga.abstract_2d_finite_volumes.util import get_revision_number, \ |
---|
| 35 | remove_lone_verts, sww2timeseries, get_centroid_values |
---|
| 36 | |
---|
| 37 | from anuga.abstract_2d_finite_volumes.neighbour_mesh import segment_midpoints |
---|
| 38 | from anuga.load_mesh.loadASCII import export_mesh_file |
---|
| 39 | from anuga.geometry.polygon import intersection |
---|
| 40 | from anuga.file_conversion.sww2dem import sww2dem |
---|
| 41 | |
---|
| 42 | from anuga.utilities.system_tools import get_vars_in_expression |
---|
| 43 | import anuga.utilities.log as log |
---|
| 44 | |
---|
| 45 | from anuga.utilities.file_utils import create_filename,\ |
---|
| 46 | get_all_swwfiles |
---|
| 47 | from anuga.file.csv_file import load_csv_as_dict |
---|
| 48 | from sww_file import Read_sww, Write_sww |
---|
| 49 | |
---|
| 50 | from anuga.anuga_exceptions import DataMissingValuesError, \ |
---|
| 51 | DataFileNotOpenError, DataTimeError, DataDomainError, \ |
---|
| 52 | NewQuantity |
---|
| 53 | |
---|
| 54 | |
---|
| 55 | # Definitions of various NetCDF dimension names, etc. |
---|
| 56 | lon_name = 'LON' |
---|
| 57 | lat_name = 'LAT' |
---|
| 58 | time_name = 'TIME' |
---|
| 59 | precision = netcdf_float # So if we want to change the precision its done here |
---|
| 60 | |
---|
| 61 | ## |
---|
| 62 | # @brief Clas for a NetCDF data file writer. |
---|
| 63 | class Write_nc: |
---|
| 64 | """Write an nc file. |
---|
| 65 | |
---|
| 66 | Note, this should be checked to meet cdc netcdf conventions for gridded |
---|
| 67 | data. http://www.cdc.noaa.gov/cdc/conventions/cdc_netcdf_standard.shtml |
---|
| 68 | """ |
---|
| 69 | |
---|
| 70 | ## |
---|
| 71 | # @brief Instantiate a Write_nc instance. |
---|
| 72 | # @param quantity_name |
---|
| 73 | # @param file_name |
---|
| 74 | # @param time_step_count The number of time steps. |
---|
| 75 | # @param time_step The time_step size. |
---|
| 76 | # @param lon |
---|
| 77 | # @param lat |
---|
| 78 | def __init__(self, |
---|
| 79 | quantity_name, |
---|
| 80 | file_name, |
---|
| 81 | time_step_count, |
---|
| 82 | time_step, |
---|
| 83 | lon, |
---|
| 84 | lat): |
---|
| 85 | """Instantiate a Write_nc instance (NetCDF file writer). |
---|
| 86 | |
---|
| 87 | time_step_count is the number of time steps. |
---|
| 88 | time_step is the time step size |
---|
| 89 | |
---|
| 90 | pre-condition: quantity_name must be 'HA', 'UA'or 'VA'. |
---|
| 91 | """ |
---|
| 92 | |
---|
| 93 | self.quantity_name = quantity_name |
---|
| 94 | quantity_units = {'HA':'CENTIMETERS', |
---|
| 95 | 'UA':'CENTIMETERS/SECOND', |
---|
| 96 | 'VA':'CENTIMETERS/SECOND'} |
---|
| 97 | |
---|
| 98 | multiplier_dic = {'HA':100.0, # To convert from m to cm |
---|
| 99 | 'UA':100.0, # and m/s to cm/sec |
---|
| 100 | 'VA':-100.0} # MUX files have positive x in the |
---|
| 101 | # Southern direction. This corrects |
---|
| 102 | # for it, when writing nc files. |
---|
| 103 | |
---|
| 104 | self.quantity_multiplier = multiplier_dic[self.quantity_name] |
---|
| 105 | |
---|
| 106 | #self.file_name = file_name |
---|
| 107 | self.time_step_count = time_step_count |
---|
| 108 | self.time_step = time_step |
---|
| 109 | |
---|
| 110 | # NetCDF file definition |
---|
| 111 | self.outfile = NetCDFFile(file_name, netcdf_mode_w) |
---|
| 112 | outfile = self.outfile |
---|
| 113 | |
---|
| 114 | #Create new file |
---|
| 115 | nc_lon_lat_header(outfile, lon, lat) |
---|
| 116 | |
---|
| 117 | # TIME |
---|
| 118 | outfile.createDimension(time_name, None) |
---|
| 119 | outfile.createVariable(time_name, precision, (time_name,)) |
---|
| 120 | |
---|
| 121 | #QUANTITY |
---|
| 122 | outfile.createVariable(self.quantity_name, precision, |
---|
| 123 | (time_name, lat_name, lon_name)) |
---|
| 124 | outfile.variables[self.quantity_name].missing_value = -1.e+034 |
---|
| 125 | outfile.variables[self.quantity_name].units = \ |
---|
| 126 | quantity_units[self.quantity_name] |
---|
| 127 | outfile.variables[lon_name][:]= ensure_numeric(lon) |
---|
| 128 | outfile.variables[lat_name][:]= ensure_numeric(lat) |
---|
| 129 | |
---|
| 130 | #Assume no one will be wanting to read this, while we are writing |
---|
| 131 | #outfile.close() |
---|
| 132 | |
---|
| 133 | ## |
---|
| 134 | # @brief Write a time-step of quantity data. |
---|
| 135 | # @param quantity_slice The data to be stored for this time-step. |
---|
| 136 | def store_timestep(self, quantity_slice): |
---|
| 137 | """Write a time slice of quantity info |
---|
| 138 | |
---|
| 139 | quantity_slice is the data to be stored at this time step |
---|
| 140 | """ |
---|
| 141 | |
---|
| 142 | # Get the variables |
---|
| 143 | time = self.outfile.variables[time_name] |
---|
| 144 | quantity = self.outfile.variables[self.quantity_name] |
---|
| 145 | |
---|
| 146 | # get index oflice to write |
---|
| 147 | i = len(time) |
---|
| 148 | |
---|
| 149 | #Store time |
---|
| 150 | time[i] = i * self.time_step #self.domain.time |
---|
| 151 | quantity[i,:] = quantity_slice * self.quantity_multiplier |
---|
| 152 | |
---|
| 153 | ## |
---|
| 154 | # @brief Close file underlying the class instance. |
---|
| 155 | def close(self): |
---|
| 156 | self.outfile.close() |
---|
| 157 | |
---|
| 158 | |
---|
| 159 | |
---|
| 160 | ## |
---|
| 161 | # @brief Write an NC elevation file. |
---|
| 162 | # @param file_out Path to the output file. |
---|
| 163 | # @param lon ?? |
---|
| 164 | # @param lat ?? |
---|
| 165 | # @param depth_vector The elevation data to write. |
---|
| 166 | def write_elevation_nc(file_out, lon, lat, depth_vector): |
---|
| 167 | """Write an nc elevation file.""" |
---|
| 168 | |
---|
| 169 | # NetCDF file definition |
---|
| 170 | outfile = NetCDFFile(file_out, netcdf_mode_w) |
---|
| 171 | |
---|
| 172 | #Create new file |
---|
| 173 | nc_lon_lat_header(outfile, lon, lat) |
---|
| 174 | |
---|
| 175 | # ELEVATION |
---|
| 176 | zname = 'ELEVATION' |
---|
| 177 | outfile.createVariable(zname, precision, (lat_name, lon_name)) |
---|
| 178 | outfile.variables[zname].units = 'CENTIMETERS' |
---|
| 179 | outfile.variables[zname].missing_value = -1.e+034 |
---|
| 180 | |
---|
| 181 | outfile.variables[lon_name][:] = ensure_numeric(lon) |
---|
| 182 | outfile.variables[lat_name][:] = ensure_numeric(lat) |
---|
| 183 | |
---|
| 184 | depth = num.reshape(depth_vector, (len(lat), len(lon))) |
---|
| 185 | outfile.variables[zname][:] = depth |
---|
| 186 | |
---|
| 187 | outfile.close() |
---|
| 188 | |
---|
| 189 | |
---|
| 190 | ## |
---|
| 191 | # @brief Write lat/lon headers to a NetCDF file. |
---|
| 192 | # @param outfile Handle to open file to write to. |
---|
| 193 | # @param lon An iterable of the longitudes. |
---|
| 194 | # @param lat An iterable of the latitudes. |
---|
| 195 | # @note Defines lat/long dimensions and variables. Sets various attributes: |
---|
| 196 | # .point_spacing and .units |
---|
| 197 | # and writes lat/lon data. |
---|
| 198 | |
---|
| 199 | def nc_lon_lat_header(outfile, lon, lat): |
---|
| 200 | """Write lat/lon headers to a NetCDF file. |
---|
| 201 | |
---|
| 202 | outfile is the netcdf file handle. |
---|
| 203 | lon - a list/array of the longitudes |
---|
| 204 | lat - a list/array of the latitudes |
---|
| 205 | """ |
---|
| 206 | |
---|
| 207 | outfile.institution = 'Geoscience Australia' |
---|
| 208 | outfile.description = 'Converted from URS binary C' |
---|
| 209 | |
---|
| 210 | # Longitude |
---|
| 211 | outfile.createDimension(lon_name, len(lon)) |
---|
| 212 | outfile.createVariable(lon_name, precision, (lon_name,)) |
---|
| 213 | outfile.variables[lon_name].point_spacing = 'uneven' |
---|
| 214 | outfile.variables[lon_name].units = 'degrees_east' |
---|
| 215 | outfile.variables[lon_name].assignValue(lon) |
---|
| 216 | |
---|
| 217 | # Latitude |
---|
| 218 | outfile.createDimension(lat_name, len(lat)) |
---|
| 219 | outfile.createVariable(lat_name, precision, (lat_name,)) |
---|
| 220 | outfile.variables[lat_name].point_spacing = 'uneven' |
---|
| 221 | outfile.variables[lat_name].units = 'degrees_north' |
---|
| 222 | outfile.variables[lat_name].assignValue(lat) |
---|
| 223 | |
---|
| 224 | |
---|
| 225 | |
---|
| 226 | |
---|
| 227 | ## |
---|
| 228 | # @brief Filter data file, selecting timesteps first:step:last. |
---|
| 229 | # @param filename1 Data file to filter. |
---|
| 230 | # @param filename2 File to write filtered timesteps to. |
---|
| 231 | # @param first First timestep. |
---|
| 232 | # @param last Last timestep. |
---|
| 233 | # @param step Timestep stride. |
---|
| 234 | def filter_netcdf(filename1, filename2, first=0, last=None, step=1): |
---|
| 235 | """Filter data file, selecting timesteps first:step:last. |
---|
| 236 | |
---|
| 237 | Read netcdf filename1, pick timesteps first:step:last and save to |
---|
| 238 | nettcdf file filename2 |
---|
| 239 | """ |
---|
| 240 | |
---|
| 241 | from Scientific.IO.NetCDF import NetCDFFile |
---|
| 242 | |
---|
| 243 | # Get NetCDF |
---|
| 244 | infile = NetCDFFile(filename1, netcdf_mode_r) #Open existing file for read |
---|
| 245 | outfile = NetCDFFile(filename2, netcdf_mode_w) #Open new file |
---|
| 246 | |
---|
| 247 | # Copy dimensions |
---|
| 248 | for d in infile.dimensions: |
---|
| 249 | outfile.createDimension(d, infile.dimensions[d]) |
---|
| 250 | |
---|
| 251 | # Copy variable definitions |
---|
| 252 | for name in infile.variables: |
---|
| 253 | var = infile.variables[name] |
---|
| 254 | outfile.createVariable(name, var.dtype.char, var.dimensions) |
---|
| 255 | |
---|
| 256 | # Copy the static variables |
---|
| 257 | for name in infile.variables: |
---|
| 258 | if name == 'time' or name == 'stage': |
---|
| 259 | pass |
---|
| 260 | else: |
---|
| 261 | outfile.variables[name][:] = infile.variables[name][:] |
---|
| 262 | |
---|
| 263 | # Copy selected timesteps |
---|
| 264 | time = infile.variables['time'] |
---|
| 265 | stage = infile.variables['stage'] |
---|
| 266 | |
---|
| 267 | newtime = outfile.variables['time'] |
---|
| 268 | newstage = outfile.variables['stage'] |
---|
| 269 | |
---|
| 270 | if last is None: |
---|
| 271 | last = len(time) |
---|
| 272 | |
---|
| 273 | selection = range(first, last, step) |
---|
| 274 | for i, j in enumerate(selection): |
---|
| 275 | log.critical('Copying timestep %d of %d (%f)' |
---|
| 276 | % (j, last-first, time[j])) |
---|
| 277 | newtime[i] = time[j] |
---|
| 278 | newstage[i,:] = stage[j,:] |
---|
| 279 | |
---|
| 280 | # Close |
---|
| 281 | infile.close() |
---|
| 282 | outfile.close() |
---|
| 283 | |
---|