Skip to content

Enable Copernicus DEM #374

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 8 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
60 changes: 43 additions & 17 deletions run_scripts/extpar_mch.aster.sh
Original file line number Diff line number Diff line change
Expand Up @@ -3,20 +3,20 @@
#SBATCH --job-name="extpar"
#SBATCH --nodes=1
#SBATCH --output="job.out"
#SBATCH --time=01:00:00
#SBATCH --time=00:23:00
#SBATCH --partition=postproc
#SBATCH --cpus-per-task=12
#SBATCH --cpus-per-task=64


#--------------------------------------------------------------------------------
# variables to define by user
#--------------------------------------------------------------------------------

# define model for which Extpar should run: c1, c2, i1, i2, i1_dev, i2_dev
# define model for which Extpar should run: c1, c2, i05, i1, i2, i1_dev, i2_dev
model="i1"

# Sandbox (make sure you have enough disk place at that location)!
sandboxdir=$SCRATCH/output_extpar/i1
sandboxdir=$SCRATCH/ExtPar/output/${model}


#--------------------------------------------------------------------------------
Expand Down Expand Up @@ -210,10 +210,40 @@ elif [[ $model == "i2_dev" ]]; then
model_grid_type=1
name_model_grid='INPUT_ICON_GRID'

elif [[ $model == "i05" ]]; then

#output file names
netcdf_output_filename="extpar_icon_grid_00005_R19B09_DOM02.nc"

# grid definition
grid_dir="/store_new/mch/msopr/glori/glori-ch500-nested/grid/"
grid_nc="icon_grid_00005_R19B09_DOM02.nc"

lsso_param='.TRUE.'
lsubtract_mean_slope='.FALSE.'

# orography raw data
ntiles_column=2
ntiles_row=4
topo_files="'ASTER_orig_T006.nc' 'ASTER_orig_T007.nc' 'ASTER_orig_T018.nc' 'ASTER_orig_T019.nc' 'ASTER_orig_T030.nc' 'ASTER_orig_T031.nc' 'ASTER_orig_T042.nc' 'ASTER_orig_T043.nc'"

#orography smoothing
lsmooth_oro='.FALSE.'
ilow_pass_oro=1
numfilt_oro=2
eps_filter=1.7

# soil: tiles
itile_mode=1

# model grid type
model_grid_type=1
name_model_grid='INPUT_ICON_GRID'

elif [[ $model == "i1" ]]; then

#output file names
netcdf_output_filename="external_parameter_icon_grid_0001_R19B08_mch.nc"
netcdf_output_filename="extpar_icon_grid_0001_R19B08_mch.nc"

# grid definition
grid_dir="/oprusers/osm/opr/data/grid_descriptions/"
Expand Down Expand Up @@ -292,11 +322,11 @@ binary_era=extpar_era_to_buffer.py
binary_ahf=extpar_ahf_to_buffer.py
binary_isa=extpar_isa_to_buffer.py
binary_tclim=extpar_cru_to_buffer.py
binary_aot=extpar_aot_to_buffer.py

# fortran executables
binary_lu=extpar_landuse_to_buffer.exe
binary_topo=extpar_topo_to_buffer.exe
binary_aot=extpar_aot_to_buffer.exe
binary_soil=extpar_soil_to_buffer.exe
binary_flake=extpar_flake_to_buffer.exe
binary_consistency_check=extpar_consistency_check.exe
Expand Down Expand Up @@ -402,6 +432,13 @@ input_alb = {
'alb_buffer_file': '${buffer_alb}',
}

input_aot = {
'iaot_type': 1,
'raw_data_aot_path': '',
'raw_data_aot_filename': '${raw_data_aot}',
'aot_buffer_file': '${buffer_aot}',
}

input_tclim = {
'raw_data_t_clim_path': '',
'raw_data_tclim_coarse': '',
Expand Down Expand Up @@ -477,16 +514,6 @@ cat > INPUT_ICON_GRID << EOF_grid_icon
/
EOF_grid_icon

cat > INPUT_AOT << EOF_aot
&aerosol_raw_data
raw_data_aot_path='./',
raw_data_aot_filename='${raw_data_aot}'
/
&aerosol_io_extpar
aot_buffer_file='${buffer_aot}',
/
EOF_aot

cat > INPUT_LU << EOF_lu
&lu_raw_data
raw_data_lu_path='',
Expand Down Expand Up @@ -564,7 +591,6 @@ isoil_data = 3,
ldeep_soil = .FALSE.,
raw_data_soil_path='',
raw_data_soil_filename='${raw_data_soil_HWSD}',
raw_data_deep_soil_filename='${raw_data_deep_soil}'
/
&soil_io_extpar
soil_buffer_file='${buffer_soil}',
Expand Down
28 changes: 24 additions & 4 deletions src/extpar_topo_to_buffer.f90
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
! V1_0 2010/12/21 Hermann Asensio
! Initial release
! V1_1 2011/01/20 Hermann Asensio
! small bug fixes accroding to Fortran compiler warnings
! small bug fixes according to Fortran compiler warnings
! V1_2 2011/03/25 Hermann Asensio
! update to support ICON refinement grids
! V1_4 2011/04/21 Anne Roches
Expand All @@ -29,6 +29,8 @@
! Add namelist switch lfilter_topo
! V2_10 2018-02-19 Juergen Helmert
! lsubtract_mean_slope, ERA-I surface temp for land points
! V2_11 2024-12-10 Christian R. Steger
! introduction of the COPERNICUS topography raw data set for external parameters
!
! Code Description:
! Language: Fortran 2003.
Expand All @@ -37,7 +39,8 @@
!>
!! @par extpar_topo_to_buffer
!!
!! This program reads in the GLOBE/ASTER/MERIT orography data set and aggregates the orographic height to the target grid
!! This program reads in the GLOBE/ASTER/MERIT/COPERNICUS orography data set
!! and aggregates the orographic height to the target grid
!! and computes the subgrid-scale orography parameters (SSO) required by the SSO-parameterization.
!!
!> Purpose: read in GLOBE/ASTER orography data and aggregate to COSMO grid
Expand Down Expand Up @@ -86,6 +89,7 @@ PROGRAM extpar_topo_to_buffer

USE mo_topo_data, ONLY: topo_aster, &
& topo_merit, &
& topo_copernicus, &
& itopo_type, &
& topo_tiles_grid, &
& topo_grid, &
Expand All @@ -106,6 +110,10 @@ PROGRAM extpar_topo_to_buffer
& merit_lat_max, &
& merit_lon_min, &
& merit_lon_max, &
& copernicus_lat_min, &
& copernicus_lat_max, &
& copernicus_lon_min, &
& copernicus_lon_max, &
& num_tiles, &
& allocate_topo_data,&
& fill_topo_data, &
Expand Down Expand Up @@ -361,8 +369,20 @@ PROGRAM extpar_topo_to_buffer
CALL logging%warning(message_text)
CALL logging%error('The chosen latitude edges are not within the MERIT domain.',__FILE__,__LINE__)
END IF


CASE(topo_copernicus)
WRITE(message_text,*)'Edges of domain: ', copernicus_lon_min,' ', copernicus_lon_max,' ', copernicus_lat_min,' ' &
& ,copernicus_lat_max
CALL logging%info(message_text)
IF (lon_geo (tg%ie,tg%je,tg%ke) > copernicus_lon_max .OR. lon_geo(1,1,1) < copernicus_lon_min) THEN
WRITE(message_text,*) 'COPERNICUS min lon is: ', copernicus_lon_min, ' and COPERNICUS max lon is: ', copernicus_lon_max
CALL logging%warning(message_text)
CALL logging%error('The chosen longitude edges are not within the COPERNICUS domain.',__FILE__,__LINE__)
END IF
IF (lat_geo(tg%ie,tg%je,tg%ke) > copernicus_lat_max .OR. lat_geo(1,1,1) < copernicus_lat_min) THEN
WRITE(message_text,*) 'COPERNICUS min lat is: ', copernicus_lat_min, ' and COPERNICUS max lat is: ', copernicus_lat_max
CALL logging%warning(message_text)
CALL logging%error('The chosen latitude edges are not within the COPERNICUS domain.',__FILE__,__LINE__)
END IF
END SELECT

CALL det_topo_tiles_grid(topo_tiles_grid)
Expand Down
6 changes: 6 additions & 0 deletions src/mo_agg_sgsl.f90
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,7 @@ MODULE mo_agg_sgsl
& topo_gl, &
& topo_aster, &
& topo_merit, &
& topo_copernicus, &
& get_fill_value_sgsl, &
& nr_tot !< total number of rows in GLOBE/ASTER data

Expand Down Expand Up @@ -401,6 +402,11 @@ SUBROUTINE agg_sgsl_data_to_target_grid(sgsl_tiles_grid, &
ndata(ie,je,ke) = ndata(ie,je,ke) + 1
sgsl(ie,je,ke) = sgsl(ie,je,ke) + sl(i)
ENDIF
CASE(topo_copernicus)
IF (sl(i) /= undef_sgsl) THEN
ndata(ie,je,ke) = ndata(ie,je,ke) + 1
sgsl(ie,je,ke) = sgsl(ie,je,ke) + sl(i)
ENDIF
END SELECT
!$OMP END CRITICAL
ENDIF
Expand Down
23 changes: 23 additions & 0 deletions src/mo_agg_topo_cosmo.f90
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,7 @@ MODULE mo_agg_topo_cosmo
& itopo_type, &
& topo_gl, &
& topo_merit, &
& topo_copernicus, &
& topo_aster

USE mo_topo_sso, ONLY: auxiliary_sso_parameter_cosmo, &
Expand Down Expand Up @@ -271,6 +272,9 @@ SUBROUTINE agg_topo_data_to_target_grid_cosmo(topo_tiles_grid, &
CASE(topo_merit)
hh = undef_topo
h_3rows = undef_topo
CASE(topo_copernicus)
hh = undef_topo
h_3rows = undef_topo
END SELECT

! initialize some variables
Expand Down Expand Up @@ -623,6 +627,25 @@ SUBROUTINE agg_topo_data_to_target_grid_cosmo(topo_tiles_grid, &
h22(ie,je,ke) = h22(ie,je,ke) + dhdydy(i)
ENDIF
ENDIF
CASE(topo_copernicus)
IF (h_3rows(i,j_c) /= undef_topo) THEN
ndata(ie,je,ke) = ndata(ie,je,ke) + 1
hh_target(ie,je,ke) = hh_target(ie,je,ke) + h_3rows(i,j_c)
hh2_target(ie,je,ke) = hh2_target(ie,je,ke) + (h_3rows(i,j_c) * h_3rows(i,j_c))
IF (lscale_separation) THEN
hh_target_scale(ie,je,ke) = hh_target_scale(ie,je,ke) + &
& h_3rows_scale(i,j_c)
hh2_target_scale(ie,je,ke) = hh2_target_scale(ie,je,ke) + &
& (h_3rows_scale(i,j_c) * h_3rows_scale(i,j_c))
hh_sqr_diff(ie,je,ke) = hh_sqr_diff(ie,je,ke)+ &
& (h_3rows(i,j_c) - h_3rows_scale(i,j_c))**2
ENDIF
IF(lsso_param) THEN
h11(ie,je,ke) = h11(ie,je,ke) + dhdxdx(i)
h12(ie,je,ke) = h12(ie,je,ke) + dhdxdy(i)
h22(ie,je,ke) = h22(ie,je,ke) + dhdydy(i)
ENDIF
ENDIF
END SELECT
!$OMP END CRITICAL
ENDIF
Expand Down
8 changes: 6 additions & 2 deletions src/mo_agg_topo_icon.f90
Original file line number Diff line number Diff line change
Expand Up @@ -68,15 +68,16 @@ MODULE mo_agg_topo_icon

USE mo_icon_domain, ONLY: icon_domain, grid_cells

USE mo_topo_data, ONLY: ntiles, & !< there are 16/240 GLOBE/ASTER/MERIT tiles
USE mo_topo_data, ONLY: ntiles, & !< there are 16/240 GLOBE/ASTER/MERIT/COPERNICUS tiles
& max_tiles, &
& nc_tot, & !< number of total GLOBE/ASTER columns un a latitude circle
& nr_tot, & !< total number of rows in GLOBE/ASTER/MERIT data
& nr_tot, & !< total number of rows in GLOBE/ASTER/MERIT/COPERNICUS data
& get_fill_value, &
& itopo_type, &
& topo_gl, &
& topo_aster, &
& topo_merit, &
& topo_copernicus, &
& get_varname

USE mo_topo_sso, ONLY: auxiliary_sso_parameter_icon,&
Expand Down Expand Up @@ -272,6 +273,9 @@ SUBROUTINE agg_topo_data_to_target_grid_icon(topo_tiles_grid, &
CASE(topo_merit)
hh = default_topo
hh_red = default_topo
CASE(topo_copernicus)
hh = default_topo
hh_red = default_topo
END SELECT

! initialize some variables
Expand Down
11 changes: 10 additions & 1 deletion src/mo_extpar_output_nc.f90
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,7 @@ MODULE mo_extpar_output_nc

USE mo_soil_data, ONLY: HWSD_data

USE mo_topo_data, ONLY: itopo_type, topo_aster, topo_gl, topo_merit
USE mo_topo_data, ONLY: itopo_type, topo_aster, topo_gl, topo_merit, topo_copernicus

USE mo_cosmo_grid, ONLY: lon_rot, lat_rot

Expand Down Expand Up @@ -1717,6 +1717,9 @@ SUBROUTINE set_cdi_global_att_icon(global_attributes,itopo_type,name_lookup_tabl
CASE(3) ! topo_merit
global_attributes(3)%attributetext=TRIM(lu_dataset)//', FAO DSMW, MERIT, Lake Database'
IF(isoil_data == HWSD_data) global_attributes(3)%attributetext=TRIM(lu_dataset)//', HWSD, MERIT, Lake Database'
CASE(4) ! topo_copernicus
global_attributes(3)%attributetext=TRIM(lu_dataset)//', FAO DSMW, COPERNICUS, Lake Database'
IF(isoil_data == HWSD_data) global_attributes(3)%attributetext=TRIM(lu_dataset)//', HWSD, COPERNICUS, Lake Database'
END SELECT


Expand Down Expand Up @@ -1787,6 +1790,12 @@ SUBROUTINE set_global_att_extpar(global_attributes,name_lookup_table_lu,lu_datas
ELSE
global_attributes(3)%attributetext=TRIM(lu_dataset)//', FAO DSMW, MERIT, Lake Database'
ENDIF
CASE(topo_copernicus)
IF (isoil_data >= HWSD_data) THEN
global_attributes(3)%attributetext=TRIM(lu_dataset)//', HWSD, COPERNICUS, Lake Database'
ELSE
global_attributes(3)%attributetext=TRIM(lu_dataset)//', FAO DSMW, COPERNICUS, Lake Database'
ENDIF
END SELECT
global_attributes(4)%attname = 'note'
global_attributes(4)%attributetext='Landuse data look-up table: '//TRIM(name_lookup_table_lu)
Expand Down
17 changes: 13 additions & 4 deletions src/mo_preproc_for_sgsl.f90
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ MODULE mo_preproc_for_sgsl
CONTAINS

! wrapper function for the preprocessing of raw orography data
! and calls the right subroutine for itopo_type (GLOBE, ASTER or MERIT)
! and calls the right subroutine for itopo_type (GLOBE, ASTER, MERIT or COPERNICUS)
SUBROUTINE preproc_orography (raw_data_orography_path, &
& topo_files, &
& sgsl_files, &
Expand Down Expand Up @@ -108,6 +108,8 @@ SUBROUTINE topo_subgrid_slope( infile, outfile, itopo_type)
& dlon_aster = 1./3600. , &! resolution
& dlat_merit = 3./3600. , &! resolution
& dlon_merit = 3./3600. , &! resolution
& dlat_copernicus = 1./3600. , &! resolution
& dlon_copernicus = 1./3600. , &! resolution
& eps = 1.E-9, &
& add_offset = 0., &
& scale_factor = 0.001, &
Expand Down Expand Up @@ -160,8 +162,10 @@ SUBROUTINE topo_subgrid_slope( infile, outfile, itopo_type)
status = nf90_inq_varid(ncid,"altitude", varid)
ELSE IF (itopo_type == 2) THEN
status = nf90_inq_varid(ncid,"Z", varid)
ELSE
ELSE IF (itopo_type == 3) THEN
status = nf90_inq_varid(ncid,"Elevation", varid)
ELSE
status = nf90_inq_varid(ncid,"elevation", varid)
ENDIF
CALL check_err(status, __FILE__, __LINE__)

Expand All @@ -186,9 +190,12 @@ SUBROUTINE topo_subgrid_slope( infile, outfile, itopo_type)
ELSE IF (itopo_type == 2) THEN !ASTER
dlat = dlat_aster
dlon = dlon_aster
ELSE !MERIT
ELSE IF (itopo_type == 3) THEN !MERIT
dlat = dlat_merit
dlon = dlon_merit
ELSE !COPERNICUS
dlat = dlat_copernicus
dlon = dlon_copernicus
ENDIF !itopo_type

! compute values for inner domain
Expand Down Expand Up @@ -467,8 +474,10 @@ SUBROUTINE topo_subgrid_slope( infile, outfile, itopo_type)
status = nf90_put_att(ncido, outid,'comment',trim(comment))
ELSE IF (itopo_type == 2) THEN
status = nf90_put_att(ncido, outid,'comment','ASTER tile: '//infile)
ELSE IF (itopo_type == 3) THEN
status = nf90_put_att(ncido, outid,'comment','MERIT tile: '//infile)
ELSE
status = nf90_put_att(ncido, outid,'comment','MERIT tile: '//infile)
status = nf90_put_att(ncido, outid,'comment','COPERNICUS tile: '//infile)
ENDIF
CALL check_err(status, __FILE__, __LINE__)

Expand Down
Loading