GopherProxy

	tpism_python.py - pism-exp-gsw - ice stream and sediment transport experiments
	git clone git://src.adamsgaard.dk/pism-exp-gsw
	Log
	Files
	Refs
	README
	LICENSE
	---
	tpism_python.py (2829B)
	---
	1 #!/usr/bin/env python3
	2
	3 # Copyright (C) 2009-2015, 2018 the PISM Authors
	4
	5 # @package pism_python
	6 # \author the PISM authors
	7 # \brief Creates "from scratch" a boring dataset with the right format
	8 # to use as a PISM bootstrapping file.
	9 # \details Example use of Python for this purpose.
	10 #
	11 # Usage, including a minimal PISM call to bootstrap from this file:
	12 #
	13 # \verbatim $ pism_python.py # creates foo.nc \endverbatim
	14 # \verbatim $ pismr -i foo.nc -bootstrap -Mx 41 -My 41 -Mz 21 -Lz 4000 -…
	15
	16 import sys
	17 import time
	18 import numpy as np
	19
	20 # try different netCDF modules
	21 try:
	22 from netCDF4 import Dataset as CDF
	23 except:
	24 print("netCDF4 is not installed!")
	25 sys.exit(1)
	26
	27 # set up the grid:
	28 Lx = 1e6
	29 Ly = 1e6
	30 Mx = 51
	31 My = 71
	32 x = np.linspace(-Lx, Lx, Mx)
	33 y = np.linspace(-Ly, Ly, My)
	34
	35 # create dummy fields
	36 [xx, yy] = np.meshgrid(x, y) # if there were "ndgrid" in numpy we would…
	37 acab = np.zeros((Mx, My))
	38 artm = np.zeros((Mx, My)) + 273.15 + 10.0 # 10 degrees Celsius
	39 topg = 1000.0 + 200.0 * (xx + yy) / max(Lx, Ly) # change "1000.0" to "0…
	40 # flotation criterion, etc.
	41 thk = 3000.0 * (1.0 - 3.0 * (xx 2 + yy 2) / Lx ** 2)
	42 thk[thk < 0.0] = 0.0
	43
	44 # Output filename
	45 ncfile = 'foo.nc'
	46
	47 # Write the data:
	48 nc = CDF(ncfile, "w", format='NETCDF3_CLASSIC') # for netCDF4 module
	49
	50 # Create dimensions x and y
	51 nc.createDimension("x", size=Mx)
	52 nc.createDimension("y", size=My)
	53
	54 x_var = nc.createVariable("x", 'f4', dimensions=("x",))
	55 x_var.units = "m"
	56 x_var.long_name = "easting"
	57 x_var.standard_name = "projection_x_coordinate"
	58 x_var[:] = x
	59
	60 y_var = nc.createVariable("y", 'f4', dimensions=("y",))
	61 y_var.units = "m"
	62 y_var.long_name = "northing"
	63 y_var.standard_name = "projection_y_coordinate"
	64 y_var[:] = y
	65
	66 fill_value = np.nan
	67
	68
	69 def def_var(nc, name, units, fillvalue):
	70 # dimension transpose is standard: "float thk(y, x)" in NetCDF file
	71 var = nc.createVariable(name, 'f', dimensions=("y", "x"), fill_value…
	72 var.units = units
	73 return var
	74
	75
	76 bed_var = def_var(nc, "topg", "m", fill_value)
	77 bed_var.standard_name = "bedrock_altitude"
	78 bed_var[:] = topg
	79
	80 thk_var = def_var(nc, "thk", "m", fill_value)
	81 thk_var.standard_name = "land_ice_thickness"
	82 thk_var[:] = thk
	83
	84 acab_var = def_var(nc, "climatic_mass_balance", "m year-1", fill_value)
	85 acab_var.standard_name = "land_ice_surface_specific_mass_balance"
	86 acab_var[:] = acab
	87
	88 artm_var = def_var(nc, "ice_surface_temp", "K", fill_value)
	89 artm_var[:] = artm
	90
	91 # set global attributes
	92 nc.Conventions = "CF-1.4"
	93 historysep = ' '
	94 historystr = time.asctime() + ': ' + historysep.join(sys.argv) + '\n'
	95 setattr(nc, 'history', historystr)
	96
	97 nc.close()
	98 print(' PISM-bootable NetCDF file %s written' % ncfile)
	99 print(' for example, run:')
	100 print(' $ pismr -i foo.nc -bootstrap -Mx 41 -My 41 -Mz 21 -Lz 4000 -M…