Skip to content

GitLab

Commit 4755ecbd authored by Christoph Heim's avatar Christoph Heim
Browse files

MPAS new version. FV3 planned.

parent 3d5c4e37
Hide whitespace changes
Inline Side-by-side
01_nicam.py
View file @ 4755ecbd
......@@ -76,11 +76,11 @@ if __name__ == '__main__':
var_names = ['QV', 'QC', 'T', 'W',
'U10M', 'V10M', 'T2M', 'LWUTOA', 'SWDSFC',
'SLHFLX', 'SSHFLX', 'TQC']
#var_names = ['T2M']
var_names = ['T']
# model resolutions [km] of simulations
ress = [7, 3.5]
#ress = [7]
ress = [7]
#ress = [3.5]
# date range
......
03_icon.py
View file @ 4755ecbd
......@@ -147,6 +147,9 @@ if __name__ == '__main__':
'2.5km/icon_grid_0017_R02B10_G.nc'),
},
}
os.environ['GRIB_DEFINITION_PATH'] = ('/mnt/lustre01/sw/rhel6-x64/eccodes'+
'/definitions')
###########################################################################
## PREPARING STEPS
......
05_mpas.py
View file @ 4755ecbd
......@@ -4,10 +4,9 @@
description: Extract lat-lon box of data from model MPAS.
author: Christoph Heim
date created: 05.07.2019
date changed: 18.07.2019
date changed: 19.07.2019
usage: arguments:
1st: n jobs for multiprocessing pool
for 2D field: 8 jobs possible
python: 3.5.2
"""
###############################################################################
......@@ -51,66 +50,102 @@ def sellatlon_MPAS(inp_file, out_file, dt, box, options, var_dict,
TM = Timer()
file_code = '{}km_{}_{:%Y%m%d%H%M}'.format(res, var_name, dt)
if os.path.exists(out_file) and not options['recompute']:
TM.start('prepr')
TM.start('cdo')
print('\t\t'+file_code+' already computed')
TM.stop('cdo')
TM.stop('prepr')
else:
print('\t\t'+file_code)
split = os.path.split(out_file)
tmp_file = os.path.join(split[0],'tmp_'+split[1])
# first step
TM.start('prepr')
if not os.path.exists(tmp_file):
ofile = cdo.setattribute('*@axis="txz"',
input=(
'-selname,'+var_dict['file']+' '+inp_file
),
output=tmp_file, options='-f nc4')
TM.stop('prepr')
TM.start('grid')
if not os.path.exists(weights_file):
comp_weights_file(target_grid, weights_file,
tmp_file, grid_def_file,
res, box, options)
TM.stop('grid')
# cdo
TM.start('cdo')
#ofile = cdo.remap(target_grid, weights_file,
# input=(
# ' -sellevidx,'+
# str(box['vert0'])+'/'+str(box['vert1'])+
# ' -setgrid,mpas:'+grid_def_file+
# ' '+tmp_file),
# output=out_file, options='-f nc')
input = ("-sellonlatbox,{},{},{},{} -sellevidx,{}/{}"+
" -setgrid,mpas:{} {}").format(
print(file_code)
domain_str = "{},{},{},{}".format(
box['lon'].start, box['lon'].stop,
box['lat'].start, box['lat'].stop,
box['vert0'], box['vert1'],
grid_def_file,
tmp_file)
ofile = cdo.remap(target_grid, weights_file,
input=input, output=out_file, options='-f nc')
box['lat'].start, box['lat'].stop)
levels_str = "{}/{}".format(
box['vert0'], box['vert1'])
# delete tmp_file
TM.start('prepr')
if options['rm_tmp_files']:
os.remove(tmp_file)
TM.stop('prepr')
TM.stop('cdo')
TM.start('run')
subprocess.call(['./run_MPAS.sh', domain_str, levels_str,
str(res), os.path.split(out_file)[0],
os.path.split(out_file)[1][:-3],
var_dict['file'], inp_file, target_grid,
])
TM.stop('run')
#quit()
#quit()
#if os.path.exists(out_file) and not options['recompute']:
# TM.start('prepr')
# TM.start('cdo')
# print('\t\t'+file_code+' already computed')
# TM.stop('cdo')
# TM.stop('prepr')
#else:
# print('\t\t'+file_code)
# split = os.path.split(out_file)
# tmp_file = os.path.join(split[0],'tmp_'+split[1])
# # first step
# TM.start('prepr')
# #if not os.path.exists(tmp_file):
# # ofile = cdo.setattribute('*@axis="txz"',
# # input=(
# # '-selname,'+var_dict['file']+' '+inp_file
# # ),
# # output=tmp_file, options='-f nc4')
# TM.stop('prepr')
# TM.start('grid')
# if not os.path.exists(weights_file):
# comp_weights_file(target_grid, weights_file,
# tmp_file, grid_def_file,
# res, box, options)
# TM.stop('grid')
# # cdo
# TM.start('cdo')
# #ofile = cdo.remap(target_grid, weights_file,
# # input=(
# # ' -sellevidx,'+
# # str(box['vert0'])+'/'+str(box['vert1'])+
# # ' -setgrid,mpas:'+grid_def_file+
# # ' '+tmp_file),
# # output=out_file, options='-f nc')
# #input = ("-sellonlatbox,{},{},{},{} -sellevidx,{}/{}"+
# # " -setgrid,mpas:{} {}").format(
# # box['lon'].start, box['lon'].stop,
# # box['lat'].start, box['lat'].stop,
# # box['vert0'], box['vert1'],
# # grid_def_file,
# # tmp_file)
# #ofile = cdo.remap(target_grid, weights_file,
# # input=input, output=out_file, options='-f nc')
# zaxis_def_file = ('/work/ka1081/2019_06_Hackathon_Mainz/christoph_heim/'+
# 'newdata/MPAS_3.75/MPAS_3.75km_zgrid.nc')
# input = ("-sellonlatbox,{},{},{},{} -sellevidx,{}/{}"+
# " -setgridtype,regular -selname,{}"+
# " -setgrid,mpas:{} -selgrid,1 {}").format(
# box['lon'].start, box['lon'].stop,
# box['lat'].start, box['lat'].stop,
# box['vert0'], box['vert1'],
# var_dict['file'],
# grid_def_file,# zaxis_def_file,
# inp_file)
# print(input)
# ofile = cdo.remap(target_grid, weights_file,
# input=input, output=out_file, options='-f nc4 -O')
# quit()
# # delete tmp_file
# TM.start('prepr')
# if options['rm_tmp_files']:
# os.remove(tmp_file)
# TM.stop('prepr')
# TM.stop('cdo')
return(TM)
......@@ -137,11 +172,16 @@ if __name__ == '__main__':
var_names = ['QV', 'QC', 'T', 'W',
'U10M', 'V10M', 'T2M', 'LWUTOA', 'SWDSFC',
'TQC']
var_names = ['U10M', 'V10M', 'T2M', 'LWUTOA', 'SWDSFC',
'TQC']
#var_names = ['QV', 'QC', 'T', 'W']
var_names = ['QV']
# model resolutions [km] of simulations
ress = [7.5, 3.75]
ress = [3.75]
ress = [7.5]
#ress = [7.5]
# date range
first_date = datetime(2016,8,10)
......
06_ifs.py
View file @ 4755ecbd
......@@ -92,10 +92,6 @@ if __name__ == '__main__':
var_names = ['QV', 'QC', 'T', 'W',
'U10M', 'V10M', 'T2M', 'LWUTOA', 'SWDSFC',
'SLHFLX', 'SSHFLX', 'TQC']
var_names = ['QV', 'T', 'W',
'U10M', 'V10M', 'T2M', 'LWUTOA', 'SWDSFC',
'SLHFLX', 'SSHFLX', 'TQC']
#var_names = ['W']
# model resolutions [km] of simulations
......@@ -116,13 +112,38 @@ if __name__ == '__main__':
# IFS SPECIFIC SETTINGS
###########################################################################
#'QV' :{'file':'QV', 'dim':'3D',
# 'group':'mars_out_ml_moist',},
#'QC' :{'file':'clwc', 'dim':'3D',
# 'group':'mars_out_ml_moist',},
#'T' :{'file':'T', 'dim':'3D',
# 'group':'gg_mars_out_ml_upper_sh',},
#'W' :{'file':'param120.128.192','dim':'3D',
# 'group':'gg_mars_out_ml_upper_sh',},
#'U10M' :{'file':'10u', 'dim':'2D',
# 'group':'mars_out',},
#'V10M' :{'file':'10v', 'dim':'2D',
# 'group':'mars_out',},
#'T2M' :{'file':'2t', 'dim':'2D',
# 'group':'mars_out',},
#'LWUTOA':{'file':'ttr', 'dim':'2D',
# 'group':'mars_out',},
#'SWDSFC':{'file':'ssrd', 'dim':'2D',
# 'group':'mars_out',},
#'SLHFLX':{'file':'slhf', 'dim':'2D',
# 'group':'mars_out',},
#'SSHFLX':{'file':'sshf', 'dim':'2D',
# 'group':'mars_out',},
#'TQC' :{'file':'tclw', 'dim':'2D',
# 'group':'mars_out',},
var_dict = {
'QV' :{'file':'QV', 'dim':'3D',
'QV' :{'file':'q', 'dim':'3D',
'group':'mars_out_ml_moist',},
'QC' :{'file':'clwc', 'dim':'3D',
'group':'mars_out_ml_moist',},
'T' :{'file':'T', 'dim':'3D',
'T' :{'file':'t', 'dim':'3D',
'group':'gg_mars_out_ml_upper_sh',},
'W' :{'file':'param120.128.192','dim':'3D',
'group':'gg_mars_out_ml_upper_sh',},
......@@ -146,6 +167,9 @@ if __name__ == '__main__':
}
base_time = datetime(2016,8,1)
os.environ['GRIB_DEFINITION_PATH'] = ('/work/ka1081/programs/'+
'eccodes-2.9.0/share/eccodes/definitions')
###########################################################################
## PREPARING STEPS
......
07_geos.py
View file @ 4755ecbd
......@@ -72,6 +72,9 @@ if __name__ == '__main__':
var_names = ['QV', 'QC', 'T', 'H', 'W',
'U10M', 'V10M', 'T2M', 'LWUTOA', 'SWDSFC',
'SLHFLX', 'SSHFLX', 'TQC']
var_names = ['SWDSFC',
'SLHFLX', 'SSHFLX', 'TQC']
# model resolutions [km] of simulations
ress = [3]
......@@ -100,7 +103,7 @@ if __name__ == '__main__':
'V10M' :{'file':'geosgcm_surf','key':'V10M'},
'T2M' :{'file':'geosgcm_surf','key':'T2M'},
'LWUTOA':{'file':'geosgcm_surf','key':'OLR'},
'SWDSRF':{'file':'geosgcm_surf','key':'SWGDWN'},
'SWDSFC':{'file':'geosgcm_surf','key':'SWGDWN'},
'SLHFLX':{'file':'geosgcm_surf','key':'LHFX'},
'SSHFLX':{'file':'geosgcm_surf','key':'SHFX'},
'TQC' :{'file':'geosgcm_conv','key':'CWP'},
......
08_arpege.py
View file @ 4755ecbd
......@@ -257,8 +257,8 @@ if __name__ == '__main__':
# first variable that is read from grib file
main_vars = {'3D':'T','2D':'LWUTOA'}
run_var_type = '3D'
#run_var_type = '2D'
#run_var_type = '3D'
run_var_type = '2D'
var_names = var_namess[run_var_type]
main_var = main_vars[run_var_type]
......
09_fv3.py 0 → 100755
View file @ 4755ecbd
#!/usr/bin/python
# -*- coding: utf-8 -*-
"""
description: Extract lat-lon box of data from model FV3.
author: Christoph Heim
date created: 18.07.2019
date changed: 18.07.2019
usage: arguments:
1st: n jobs for multiprocessing pool
python: 3.5.2
"""
###############################################################################
import os, glob, subprocess, sys, time, shutil
import numpy as np
from datetime import datetime, timedelta
from multiprocessing import Pool
from pathlib import Path
from cdo import Cdo
from package.utilities import Timer, write_grid_file, cdo_mergetime
from namelist import domain, padding
from functions import paste_dir_names
###############################################################################
def comp_weights_file(target_grid, weights_file, inp_file, grid_def_file,
res, box, options):
"""
"""
#if (not os.path.exists(target_grid)) or (options['recompute']):
# write_grid_file(box, target_grid, res)
print('Compute weights file')
ofile = cdo.gennn(target_grid, input=(' -setgrid,'+grid_def_file+
' '+inp_file), output=weights_file,
options='-P 1')
def sellatlon_ICON(inp_file, out_file, grid_def_file, weights_file,
target_grid, dt, box, options, var_dict, var_name,
res):
"""
"""
TM = Timer()
file_code = '{}km_{}_{:%Y%m%d}'.format(res, var_name, dt)
if os.path.exists(out_file) and not options['recompute']:
TM.start('cdo')
print('\t\t'+file_code+' already computed')
TM.stop('cdo')
else:
# cdo
TM.start('cdo')
print('\t\t'+file_code)
if var_dict['dim'] == '3d':
ofile = cdo.remap(target_grid, weights_file,
input=(
' -sellevidx,'+
str(box['vert0'])+'/'+str(box['vert1'])+
' -setgrid,'+grid_def_file+
' '+inp_file),
output=out_file, options='-f nc')
elif var_dict['dim'] == '2d':
ofile = cdo.remap(target_grid, weights_file,
input=(
' -setgrid,'+grid_def_file+
' -selname,'+var_dict['key']+
' '+inp_file),
output=out_file, options='-f nc')
TM.stop('cdo')
return(TM)
if __name__ == '__main__':
# GENERAL SETTINGS
###########################################################################
# input and output directories
raw_data_dir = os.path.join('/work','ka1081','DYAMOND')
out_base_dir = os.path.join('/work','ka1081','2019_06_Hackathon_Mainz',
'christoph_heim','newdata')
# vert box to subselect
box = domain
box.update({'vert0':73-14,'vert1':91-14})
box['lon'] = slice(box['lon'].start - padding, box['lon'].stop + padding)
box['lat'] = slice(box['lat'].start - padding, box['lat'].stop + padding)
# name of model
model_name = 'ICON'
# variables to extract
var_names = ['QV', 'QC', 'T', 'W',
'U10M', 'V10M', 'T2M', 'LWUTOA', 'SWNSFC', 'SWDIFFUSFC',
'SLHFLX', 'SSHFLX', 'TQC']
#var_names = ['LWUTOA']
# model resolutions [km] of simulations
#ress = [10,5,2.5]
ress = [10,2.5]
ress = [2.5]
#ress = [10]
# date range
first_date = datetime(2016,8,10)
last_date = datetime(2016,8,19)
# options for computation
options = {}
options['recompute'] = 0
options['rm_tmp_folder'] = 0
###########################################################################
# ICON SPECIFIC SETTINGS
###########################################################################
grid_def_base_dir = os.path.join('/work','bk1040','experiments', 'input')
var_dict = {
'QV' :{'file':'qv', 'dim':'3d', },
'QC' :{'file':'tot_qc_dia', 'dim':'3d', },
'T' :{'file':'t', 'dim':'3d', },
'W' :{'file':'w', 'dim':'3d', },
'U10M' :{'file':'atm3_2d_ml', 'dim':'2d', 'key':'U_10M'},
'V10M' :{'file':'atm3_2d_ml', 'dim':'2d', 'key':'V_10M'},
'T2M' :{'file':'atm3_2d_ml', 'dim':'2d', 'key':'T_2M'},
'LWUTOA' :{'file':'atm_2d_avg_ml', 'dim':'2d', 'key':'ATHB_T'},
'SWNSFC' :{'file':'atm_2d_avg_ml', 'dim':'2d', 'key':'ASOB_S'},
'SWDIFFUSFC':{'file':'atm_2d_avg_ml', 'dim':'2d', 'key':'ASWDIFU_S'},
'SLHFLX' :{'file':'atm2_2d_ml', 'dim':'2d', 'key':'LHFL_S'},
'SSHFLX' :{'file':'atm2_2d_ml', 'dim':'2d', 'key':'SHFL_S'},
'TQC' :{'file':'atm1_2d_ml', 'dim':'2d', 'key':'TQC_DIA'},
}
grid_dict = {
10: {'grid_def_file':os.path.join(grid_def_base_dir,
'10km/icon_grid_0025_R02B08_G.nc'),
},
5: {'grid_def_file':os.path.join(grid_def_base_dir,
'5km_2/icon_grid_0015_R02B09_G.nc'),
},
2.5:{'grid_def_file':os.path.join(grid_def_base_dir,
'2.5km/icon_grid_0017_R02B10_G.nc'),
},
}
os.environ['GRIB_DEFINITION_PATH'] = ('/mnt/lustre01/sw/rhel6-x64/eccodes'+
'/definitions')
###########################################################################
## PREPARING STEPS
TM = Timer()
dt_range = np.arange(first_date, last_date + timedelta(days=1),
timedelta(days=1)).tolist()
cdo = Cdo()
if len(sys.argv) > 1:
n_tasks = int(sys.argv[1])
else:
n_tasks = 1
print('Using ' + str(n_tasks) + ' taks.')
## EXTRACT VARIABLES FROM SIMULATIONS
# args used as input to parallel executed function
args = []
for res in ress:
print('############## res ' + str(res) + ' ##################')
for var_name in var_names:
print('\t############## var ' + var_name + ' ##################')
sim_name = model_name + '-' + str(res) + 'km'
inp_dir = os.path.join(raw_data_dir, sim_name)
# out_dir = directory for final model output (after mergetime)
# out_tmp_dir = directory for output of files in time merge
# level of raw model output
out_dir, out_tmp_dir = paste_dir_names(out_base_dir,
model_name, res, domain)
Path(out_dir).mkdir(parents=True, exist_ok=True)
# remove temporary fils if desired
if options['rm_tmp_folder'] and os.path.exists(out_tmp_dir):
shutil.rmtree(out_tmp_dir)
Path(out_tmp_dir).mkdir(parents=True, exist_ok=True)
# prepare grid interpolation
# icon grid definition file
grid_def_file = grid_dict[res]['grid_def_file']
# weights file to compute once for a grid
weights_file = os.path.join(out_tmp_dir,
'weights_{}km_dom_{}'.format(res,
domain['code']))
# target grid on which to interpolate the model output
target_grid = os.path.join('grids','latlon_{}km_dom_{}'.format(
res, domain['code']))
write_grid_file(box, target_grid, res)
# find times and files that should be extracted
# and prepare arguments for function
for dt in dt_range:
inp_files_glob = glob.glob(os.path.join(inp_dir,
'*_{}_*{:%Y%m%d}*'.format(
var_dict[var_name]['file'], dt)))
inp_file = os.path.join(inp_files_glob[0])
out_file = os.path.join(out_tmp_dir,
var_name+'_{:%Y%m%d}'.format(dt)+'.nc')
args.append( (inp_file, out_file, grid_def_file,
weights_file, target_grid,
dt, box, options, var_dict[var_name],
var_name, res) )
TM.start('grid')
if ((not os.path.exists(weights_file)) or
(not os.path.exists(target_grid))):
comp_weights_file(target_grid, weights_file,
inp_file, grid_def_file,
res, box, options)
TM.stop('grid')
# run function serial or parallel
if n_tasks > 1:
with Pool(processes=n_tasks) as pool:
results = pool.starmap(sellatlon_ICON, args)
else:
results = []
for arg in args:
results.append(sellatlon_ICON(*arg))
# collect timings from subtasks
for task_TM in results:
TM.merge_timings(task_TM)
# merge all time step files to one
TM.start('merge')
for res in ress:
for var_name in var_names:
out_dir, out_tmp_dir = paste_dir_names(out_base_dir,
model_name, res, domain)
cdo_mergetime(out_tmp_dir, out_dir, var_name)
TM.stop('merge')
TM.print_report()
OLD.05_mpas.py 0 → 100755
View file @ 4755ecbd
#!/usr/bin/python
# -*- coding: utf-8 -*-
"""
description: Extract lat-lon box of data from model MPAS.
author: Christoph Heim
date created: 05.07.2019
date changed: 18.07.2019
usage: arguments:
1st: n jobs for multiprocessing pool
for 2D field: 8 jobs possible
python: 3.5.2
"""
###############################################################################
import os, glob, subprocess, sys, time, shutil
import numpy as np
from datetime import datetime, timedelta
from multiprocessing import Pool
from pathlib import Path
from cdo import Cdo
from package.utilities import Timer, write_grid_file, cdo_mergetime
from namelist import domain, padding
from functions import paste_dir_names
###############################################################################
def comp_weights_file(target_grid, weights_file, inp_file, grid_def_file,
res, box, options):
"""
"""
#if (not os.path.exists(target_grid)) or (options['recompute']):
# write_grid_file(box, target_grid, res)
print('Compute weights file')
#ofile = cdo.gennn(target_grid,
# input=(' -setgrid,mpas:'+grid_def_file+
# ' '+inp_file),
# output=weights_file)
input = "-sellonlatbox,{},{},{},{} -setgrid,mpas:{} {}".format(
box['lon'].start,box['lon'].stop,
box['lat'].start,box['lat'].stop,
grid_def_file,
inp_file)
ofile = cdo.gennn(target_grid,
input=input, output=weights_file)
def sellatlon_MPAS(inp_file, out_file, dt, box, options, var_dict,
weights_file, target_grid, var_name, res):
TM = Timer()
file_code = '{}km_{}_{:%Y%m%d%H%M}'.format(res, var_name, dt)
if os.path.exists(out_file) and not options['recompute']:
TM.start('prepr')
TM.start('cdo')
print('\t\t'+file_code+' already computed')
TM.stop('cdo')
TM.stop('prepr')
else:
print('\t\t'+file_code)
split = os.path.split(out_file)
tmp_file = os.path.join(split[0],'tmp_'+split[1])