GridWetData old tut

<source lang="python">

1. !/usr/bin/env python
2. Essential tutorial demonstrating the GridWetData API
4. Each demonstration section is independent
5. The examples assume data and grid definitions in a directory "DMI_data",
6. where the kode is executed

import os import sys

from GridWetData import * # import the GridWetData environment

3. 1) Single point interpolations of data somewhere at sometime
5. As example, we use ocean temperature
6. First argument of data_manager.DataManager is the directory with
7. data and grid descriptors; second argument (here DMIGrid_3D) is the grid class that
8. should be used to represent grids
9. TemperatureData_3DwithTime is a grid data sub class that can resolve
10. which specific data files should be loaded (when needed)
11. position is the position, where we assess the temperature, as
12. any indexable type as (lon[deg E], lat[deg N], meters below the sea surface)
13. Temperature is assessed using the call hook of the TemperatureData_3DwithTime
14. instance: tt(position, when), where when is a datetime instance
15. Necessary data for interpoaltion is first loaded, when the time argument is received,
16. to avoid loading excessive amount of data in the data base
17. The second interpolation call illustrates the caching property of
18. GridData_3DwithTime sub classes (like TemperatureData_3DwithTime).
19. GridData_3DwithTime determines that cache data is sufficient to
20. serve the requested second interpolation, and does not load more data from the data base

from GridWetData.HBM import DataManager, HBMGrid_3D, TemperatureData_3DwithTime dmg = DataManager("DMI_data", HBMGrid_3D) # interface to data base tt = TemperatureData_3DwithTime(dmg) # interface to temperature

pos = [4.2083, 56.025, 1.2] # lon[deg E], lat[deg N], meters below the sea surface print tt(pos, datetime(2015, 05, 20, 00, 16, 00)), " - should give 8.55091339007" # interpolation using call hook print tt(pos, datetime(2015, 05, 20, 00, 36, 00)), " - should give 8.54531238808" # interpolation using call hook

3. 2) Vertical scan down the water column of data somewhere at sometime
4. by manually selected data frames (i.e. not using the data manager)
6. This time, we'll not use the data manager, but do it manually to see
7. what goes on behind the scene. As example, we'll try to assess the
8. salinity near the sea bed
9. First line instantiates a DMIGrid_3D instance from a grid descriptor file
10. Second line updates grid topography corresponding to a specific sea level elevation data frame
11. Third line loads a salinity data set from a specific file (you are responsible
12. for the correspondance between data and dynamic sea level when circumventing
13. the data manager and *_3DwithTime classes)
14. Fifth line assesses full water depth at selected position (corresponding to
15. loaded sea level elevation)
16. Line 6-9 scans down the water column, interpolating salinity along the way

g3D = HBMGrid_3D( os.path.join("DMI_data", "ns_grid.nc")) # load this grid g3D.update_sea_level( os.path.join("DMI_data", "z_ns_2015_05_20_00_15_00.nc")) # manually select frame salt = GridData_3D(g3D, os.path.join("DMI_data", "s_ns_2015_05_20_00_15_00.nc")) # manually select frame

pos = [4.2083, 56.025, 0] # lon[deg E], lat[deg N], meters below the sea surface cwd = g3D.interpolate_wdepth(pos) # assess water full depth at position pos for depth in arange(0, cwd, 0.01):

   pos[2] = depth
   print depth, salt(pos) # only a position argument for GridData_3D interpolation call hook

3. 3) Surface layer / bottom layer / column average extraction of data on native grid
5. This time, we'll use the DataManager to create a GridData_3D
6. by interpolating two time frames. snapshot_at() interpolates
7. both data and sea level interpolation between time frames

dmg = DataManager("DMI_data", HBMGrid_3D) # interface to data base tt = TemperatureData_3DwithTime(dmg) # interface to temperature tframe = tt.snapshot_at(datetime(2015, 05, 20, 00, 16, 00)) # GridData_3D instance

tsurf = tframe.get_surface_layer() tbott = tframe.get_bottom_layer() tavg = tframe.get_vertical_average()

print "temperature in water column  :", tframe.data[100,200,:] print "cell thickness in water column  :", tframe.grid.cellw[100,200] print "surf cell temperature from layer  :", tsurf.data[100,200] print "bottom cell temperature from layer  :", tbott.data[100,200] print "vertical average temperature from layer :", tavg.data[100,200]

3. 4) Interpolation of data at custom provided mesh at sometime
5. GridData_3D and GridData_3DwithTime call hooks accepts a list of positions
6. and returns corresponding list of interplated values

dmg = DataManager("DMI_data", HBMGrid_3D) # interface to data base tt = TemperatureData_3DwithTime(dmg)

my_mesh = [] # create example mesh for x in arange(4, 6, 0.01):

   my_mesh.append([x, 56.025, 1.33])
   

data_at_my_mesh = tt(my_mesh, datetime(2015, 05, 20, 00, 16, 00)) print data_at_my_mesh

3. 5) Derived data: thermal front index layer on native grid
5. The module derived_properties contain functionality to
6. generate derived data from raw physical data sets

dmg = DataManager("DMI_data", HBMGrid_3D) # interface to data base temp = TemperatureData_3DwithTime(dmg).snapshot_at(datetime(2015, 05, 20, 00, 16, 00))

sindex = derived_layers.StratificationFront(temp) # GridData_2D instance sindex.write_data("sindex.nc") # save data to file - netcdf format in this case

1. 5b) the example above writes a single frame to a netcdf file.
2. It is also possible to write a series of frames - in this case you need to address the
3. specific writer (write_data_as_netCDF) as below, because suffix resolution can
4. not be applied

import netCDF4 as netcdf ncfile = netcdf.Dataset("test.nc", "w") for minute in range(16,20):

   temp    = TemperatureData_3DwithTime(dmg).snapshot_at(datetime(2015, 05, 20, 00, minute, 00))
   sindex  = derived_layers.StratificationFront(temp) # generate 2D data frame
   sindex.write_data_as_netCDF(ncfile, index=minute)       # add frame to file

ncfile.close()

<source lang="python">