Commit fd52d60b authored by Modellers Operational's avatar Modellers Operational

Add h and land_mask time invariant calculation to regrid code

parent 9bcc4223
...@@ -15,7 +15,7 @@ if rank == 0: ...@@ -15,7 +15,7 @@ if rank == 0:
grid_res = float(sys.argv[4]) grid_res = float(sys.argv[4])
depth_layers = np.asarray(sys.argv[5].split(','), dtype=float) depth_layers = np.asarray(sys.argv[5].split(','), dtype=float)
var_list_raw = sys.argv[6].split(',') var_list_raw = sys.argv[6].split(',')
varmatch = {'temp':'nodes', 'salinity':'nodes', 'u':'elements', 'v':'elements', 'zeta':'surface', 'h':'static_nodes'} varmatch = {'temp':'nodes', 'salinity':'nodes', 'u':'elements', 'v':'elements', 'zeta':'surface'}
varlist = {} varlist = {}
for this_var in var_list_raw: for this_var in var_list_raw:
varlist[this_var] = varmatch[this_var] varlist[this_var] = varmatch[this_var]
...@@ -62,8 +62,23 @@ if rank == 0: ...@@ -62,8 +62,23 @@ if rank == 0:
output_list = [] output_list = []
for this_dict in all_interped_data: for this_dict in all_interped_data:
output_list.append(this_dict[this_var]) output_list.append(this_dict[this_var])
collected_interp_data[this_var] = np.ma.masked_invalid(np.vstack(output_list)) collected_interp_data[this_var] = np.ma.masked_invalid(np.vstack(output_list))
worker.fvcom_grid.select_ll = worker.fvcom_grid.ll
h = worker._Interpolater(worker.fvcom_grid.h)
grid_in_domain = pf.read.FileReader(worker.fvcom_file).in_domain(worker.regular_grid.mesh_lons, worker.regular_grid.mesh_lats)
h[~grid_in_domain] = np.nan
h = h.T
land_mask = np.zeros([h.shape[0], h.shape[1], len(depth_layers)])
for i, this_depth in enumerate(depth_layers):
if this_depth == 0:
land_mask[~np.isnan(h)] = 1
else:
land_mask[:,:,i][h > this_depth] = 1
# write to cmems format # write to cmems format
fvcom = pf.read.FileReader(fvcom_file) fvcom = pf.read.FileReader(fvcom_file)
...@@ -80,7 +95,8 @@ if rank == 0: ...@@ -80,7 +95,8 @@ if rank == 0:
'salinity':{'standard_name':'sea_water_salinity','units':'psu', 'missing_value':-32768, 'long_name':'Sea Water Salinity'}, 'salinity':{'standard_name':'sea_water_salinity','units':'psu', 'missing_value':-32768, 'long_name':'Sea Water Salinity'},
'u':{'standard_name':'eastward_sea_water_velocity','units':'m s-1', 'missing_value':-32768, 'long_name':'Eastward Current Velocity'}, 'u':{'standard_name':'eastward_sea_water_velocity','units':'m s-1', 'missing_value':-32768, 'long_name':'Eastward Current Velocity'},
'v':{'standard_name':'northward_sea_water_velocity','units':'m s-1', 'missing_value':-32768, 'long_name':'Northward Current Velocity'}, 'v':{'standard_name':'northward_sea_water_velocity','units':'m s-1', 'missing_value':-32768, 'long_name':'Northward Current Velocity'},
'zeta':{'standard_name':'sea_surface_height_above_geoid','units':'m', 'missing_value':-32768,'long_name':'Sea surface height above geoid'}} 'zeta':{'standard_name':'sea_surface_height_above_geoid','units':'m', 'missing_value':-32768,'long_name':'Sea surface height above geoid'},
'mask':{'standard_name':'sea_binary_mask','units':'1', 'long_name':'Land-sea mask: sea = 1 ; land = 0'}}
globals = {'type': 'OPERATIONAL MODEL REGULAR GRID OUTPUT', globals = {'type': 'OPERATIONAL MODEL REGULAR GRID OUTPUT',
'title': 'Regularly gridded data interpolated from FVCOM output', 'title': 'Regularly gridded data interpolated from FVCOM output',
...@@ -102,14 +118,13 @@ if rank == 0: ...@@ -102,14 +118,13 @@ if rank == 0:
outfile.add_variable('lon', worker.regular_grid.lons, ['lon'], attributes=all_attributes['lon']) outfile.add_variable('lon', worker.regular_grid.lons, ['lon'], attributes=all_attributes['lon'])
outfile.add_variable('lat', worker.regular_grid.lats, ['lat'], attributes=all_attributes['lat']) outfile.add_variable('lat', worker.regular_grid.lats, ['lat'], attributes=all_attributes['lat'])
outfile.add_variable('depth', worker.regular_grid.dep_lays, ['depth'], attributes=all_attributes['depth']) outfile.add_variable('depth', worker.regular_grid.dep_lays, ['depth'], attributes=all_attributes['depth'])
outfile.add_variable('mask', land_mask, ['lon', 'lat', 'depth'], attributes=all_attributes['mask'])
outfile.add_variable('h', h, ['lon', 'lat'], attributes=all_attributes['h'])
for this_var, this_mode in varlist.items(): for this_var, this_mode in varlist.items():
if this_mode == 'surface': if this_mode == 'surface':
outfile.add_variable(this_var, collected_interp_data[this_var], ['time', 'lon', 'lat'], outfile.add_variable(this_var, collected_interp_data[this_var], ['time', 'lon', 'lat'],
attributes=all_attributes[this_var], ncopts=ncopts) attributes=all_attributes[this_var], ncopts=ncopts)
elif this_mode == 'static_nodes':
outfile.add_variable(this_var, collected_interp_data[this_var], ['lon', 'lat'],
attributes=all_attributes[this_var], ncopts=ncopts)
else: else:
outfile.add_variable(this_var, collected_interp_data[this_var], ['time', 'lon', 'lat', 'depth'], outfile.add_variable(this_var, collected_interp_data[this_var], ['time', 'lon', 'lat', 'depth'],
attributes=all_attributes[this_var], ncopts=ncopts) attributes=all_attributes[this_var], ncopts=ncopts)
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment