From dc6eb92b72af348b4a60346f801a53ad04ccadc8 Mon Sep 17 00:00:00 2001 From: Nick Wogan Date: Thu, 27 Oct 2022 15:58:21 -0700 Subject: [PATCH] added ability to specify the wavelength bins file --- CMakeLists.txt | 2 +- src/clima_types.f90 | 1 + src/clima_types_create.f90 | 8 ++++++++ src/radtran/clima_radtran.f90 | 6 +++--- src/radtran/clima_radtran_types.f90 | 4 +++- src/radtran/clima_radtran_types_create.f90 | 10 ++++++++-- 6 files changed, 24 insertions(+), 7 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 5ad28a5..baa2738 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,6 +1,6 @@ cmake_minimum_required(VERSION "3.14") -project(CLIMA LANGUAGES Fortran C VERSION "0.2.2") +project(CLIMA LANGUAGES Fortran C VERSION "0.2.3") set(CMAKE_Fortran_MODULE_DIRECTORY "${CMAKE_BINARY_DIR}/modules") diff --git a/src/clima_types.f90 b/src/clima_types.f90 index 9b4f751..0777fb8 100644 --- a/src/clima_types.f90 +++ b/src/clima_types.f90 @@ -50,6 +50,7 @@ module clima_types real(dp), allocatable :: solar_zenith ! optical-properties + character(:), allocatable :: wavelength_bins_file character(s_str_len), allocatable :: gases(:) character(s_str_len), allocatable :: particles(:) type(SettingsOpacity), allocatable :: sol diff --git a/src/clima_types_create.f90 b/src/clima_types_create.f90 index 5b0215b..d748ff9 100644 --- a/src/clima_types_create.f90 +++ b/src/clima_types_create.f90 @@ -671,9 +671,17 @@ subroutine unpack_settingsopticalproperties(op_prop, filename, s, err) type(type_dictionary), pointer :: tmp_dict type(type_list), pointer :: tmp + type(type_scalar), pointer :: scalar type (type_error), allocatable :: io_err integer :: ind + ! wavelength bins file + scalar => op_prop%get_scalar('wavelength-bins-file',required=.false.,error = io_err) + if (allocated(io_err)) then; err = trim(filename)//trim(io_err%message); return; endif + if (associated(scalar)) then + s%wavelength_bins_file = trim(scalar%string) + endif + ! species tmp_dict => op_prop%get_dictionary("species", required=.false., error=io_err) if (allocated(io_err)) then; err = trim(filename)//trim(io_err%message); return; endif diff --git a/src/radtran/clima_radtran.f90 b/src/radtran/clima_radtran.f90 index 640ff6c..7633fd3 100644 --- a/src/radtran/clima_radtran.f90 +++ b/src/radtran/clima_radtran.f90 @@ -157,7 +157,7 @@ function create_RadtranIR_2(datadir, species_names, particle_names, s, nz, err) err = '"'//s%filename//'/optical-properties/ir" does not exist.' return endif - rad%ir = OpticalProperties(datadir, IROpticalProperties, species_names, particle_names, s%ir, err) + rad%ir = OpticalProperties(datadir, IROpticalProperties, species_names, particle_names, s%ir, s%wavelength_bins_file, err) if (allocated(err)) return ! work arrays @@ -291,14 +291,14 @@ function create_Radtran_2(datadir, species_names, particle_names, s, star_f, sol err = '"'//s%filename//'/optical-properties/ir" does not exist.' return endif - rad%ir = OpticalProperties(datadir, IROpticalProperties, species_names, particle_names, s%ir, err) + rad%ir = OpticalProperties(datadir, IROpticalProperties, species_names, particle_names, s%ir, s%wavelength_bins_file, err) if (allocated(err)) return if (.not. allocated(s%sol)) then err = '"'//s%filename//'/optical-properties/solar" does not exist.' return endif - rad%sol = OpticalProperties(datadir, SolarOpticalProperties, species_names, particle_names, s%sol, err) + rad%sol = OpticalProperties(datadir, SolarOpticalProperties, species_names, particle_names, s%sol, s%wavelength_bins_file, err) if (allocated(err)) return ! photons hitting the planet diff --git a/src/radtran/clima_radtran_types.f90 b/src/radtran/clima_radtran_types.f90 index 2968e72..d62c544 100644 --- a/src/radtran/clima_radtran_types.f90 +++ b/src/radtran/clima_radtran_types.f90 @@ -124,13 +124,15 @@ module clima_radtran_types end type interface - module function create_OpticalProperties(datadir, optype, species_names, particle_names, sop, err) result(op) + module function create_OpticalProperties(datadir, optype, species_names, & + particle_names, sop, wavelength_bins_file, err) result(op) use clima_types, only: SettingsOpacity character(*), intent(in) :: datadir integer, intent(in) :: optype character(*), intent(in) :: species_names(:) character(*), intent(in) :: particle_names(:) type(SettingsOpacity), intent(in) :: sop + character(:), allocatable, intent(in) :: wavelength_bins_file character(:), allocatable, intent(out) :: err type(OpticalProperties) :: op end function diff --git a/src/radtran/clima_radtran_types_create.f90 b/src/radtran/clima_radtran_types_create.f90 index 95e2069..8c4591b 100644 --- a/src/radtran/clima_radtran_types_create.f90 +++ b/src/radtran/clima_radtran_types_create.f90 @@ -174,7 +174,8 @@ function create_Ksettings(sop) result(kset) end function - module function create_OpticalProperties(datadir, optype, species_names, particle_names, sop, err) result(op) + module function create_OpticalProperties(datadir, optype, species_names, & + particle_names, sop, wavelength_bins_file, err) result(op) use fortran_yaml_c, only: YamlFile use clima_const, only: c_light, s_str_len use clima_types, only: SettingsOpacity @@ -183,6 +184,7 @@ module function create_OpticalProperties(datadir, optype, species_names, particl character(*), intent(in) :: species_names(:) character(*), intent(in) :: particle_names(:) type(SettingsOpacity), intent(in) :: sop + character(:), allocatable, intent(in) :: wavelength_bins_file character(:), allocatable, intent(out) :: err type(OpticalProperties) :: op @@ -196,7 +198,11 @@ module function create_OpticalProperties(datadir, optype, species_names, particl op%op_type = optype ! get the bins - filename = datadir//"/kdistributions/bins.h5" + if (allocated(wavelength_bins_file)) then + filename = wavelength_bins_file ! custom bins file, specified in settings + else + filename = datadir//"/kdistributions/bins.h5" ! default + endif call read_wavl(filename, op%op_type, op%wavl, err) if (allocated(err)) return