offline_exf_seaice/input_ad.obcs/gendata.py

af61e5eb16 Mart*0001 #!/usr/bin/env python
                0002 # -*- coding: iso-8859-15 -*-
                0003 ######################## -*- coding: utf-8 -*-
                0004 """Usage: gendata.py
                0005 generates extra input files of input_ad.obcs test
                0006 """
                0007
                0008 import numpy as np
                0009 import sys
                0010
                0011 nx=80
                0012 ny=42
                0013 nr=1
                0014 nt=1
                0015
                0016 def writefield(fname,data):
                0017     """Call signatures::
                0018
                0019     writefield(filename, numpy.ndarray)
                0020
                0021     Write unblocked binary data.
                0022     """
                0023
                0024     if sys.byteorder == 'little': data.byteswap(True)
                0025
                0026     fid = open(fname,"wb")
                0027     data.tofile(fid)
                0028     fid.close()
                0029
                0030     # go back to original byte ordering
                0031     if sys.byteorder == 'little': data.byteswap(True)
                0032
                0033 def readfield(fname,dims,datatype):
                0034     """Call signatures::
                0035
                0036     readfield(filename, dims, numpy.datatype)
                0037
                0038     Read unblocked binary data with dimentions "dims".
                0039     """
                0040     fid = open(fname,"rb")
                0041     v   = np.fromfile(fid, datatype)
                0042     fid.close()
                0043
                0044     if sys.byteorder == 'little': v.byteswap(True)
                0045
                0046     if   len(v) == np.prod(dims):     v = v.reshape(dims)
                0047     elif len(v) == np.prod(dims[1:]): v = v.reshape(dims[1:])
                0048     else:
                0049         errstr = (  "dimensions do not match: \n len(data) = " + str(len(v))
                0050                   + ", but prod(dims) = " + str(np.prod(dims)) )
                0051         raise RuntimeError(errstr)
                0052
                0053     return v
                0054
                0055 # use default topography and modify it
                0056 fld = readfield('../input/bathy_3c.bin',[ny,nx],'float64')
                0057
                0058 fmin = fld.min()
                0059
                0060 bnew = np.copy(fld)
                0061 bnew[0,:17] = fmin
                0062 bnew[0,63:] = fmin
                0063
                0064 writefield('./bathy_3c.obcs',bnew)
                0065
                0066 # use initial conditions (of constant ocean) to derive obcs conditions
                0067 # for ocean
                0068 ufld = readfield('../input/uVel_3c0.bin',[ny,nx],'float64')
                0069 vfld = readfield('../input/vVel_3c0.bin',[ny,nx],'float64')
                0070
                0071 uS = ufld[1,:]
                0072 vS = vfld[2,:]
                0073 uN = ufld[-1,:]
                0074 vN = vfld[-1,:]
                0075
                0076 uW = ufld[:,1]
                0077 vW = vfld[:,2]
                0078 uE = ufld[:,-1]
                0079 vE = vfld[:,-1]
                0080
                0081 writefield('OBSu.bin',uS)
                0082 writefield('OBSv.bin',vS)
                0083 writefield('OBNu.bin',uN)
                0084 writefield('OBNv.bin',vN)
                0085 writefield('OBEu.bin',uE)
                0086 writefield('OBEv.bin',vE)
                0087 writefield('OBWu.bin',uW)
                0088 writefield('OBWv.bin',vW)
                0089
                0090 # initial conditions for sea ice variables
                0091 iceC0 = 1.
                0092 iceConc = np.ones((ny,nx))*iceC0
                0093 iceConc[bnew==0] = 0.
                0094 writefield('ice0_area.obcs',iceConc)
                0095
                0096 iceH0=0.2
                0097 iceVol=iceConc*iceH0
                0098 writefield('ice0_heff.obcs',iceVol)
                0099
                0100 wind0 = 10.
                0101 windy = np.ones((ny,nx))*wind0
                0102 windy[:,:nx//2] = -wind0
                0103 writefield('windy.bin',windy)