diff --git a/cmake/dependencies.cmake b/cmake/dependencies.cmake index df84a956..81d930e2 100644 --- a/cmake/dependencies.cmake +++ b/cmake/dependencies.cmake @@ -75,7 +75,7 @@ if (MUSICA_ENABLE_TUVX AND MUSICA_BUILD_C_CXX_INTERFACE) set(TUVX_INSTALL_INCLUDE_DIR ${MUSICA_INSTALL_INCLUDE_DIR} CACHE STRING "" FORCE) set_git_default(TUVX_GIT_REPOSITORY https://github.com/NCAR/tuv-x.git) - set_git_default(TUVX_GIT_TAG 80f896a0fb591987c2a79209377bd6f599b4fb6f) + set_git_default(TUVX_GIT_TAG 674ee1e72853bb44d23c36602fa73c955b2f021d) FetchContent_Declare(tuvx GIT_REPOSITORY ${TUVX_GIT_REPOSITORY} diff --git a/docker/Dockerfile b/docker/Dockerfile index 6fa6c930..31f319f2 100644 --- a/docker/Dockerfile +++ b/docker/Dockerfile @@ -9,7 +9,6 @@ RUN dnf -y update \ gfortran \ gdb \ git \ - lapack-devel \ lcov \ make \ netcdf-fortran-devel \ diff --git a/docker/Dockerfile.fortran-gcc b/docker/Dockerfile.fortran-gcc index d21acb7e..0a181191 100644 --- a/docker/Dockerfile.fortran-gcc +++ b/docker/Dockerfile.fortran-gcc @@ -13,7 +13,6 @@ RUN dnf -y update \ git \ hdf5-devel \ json-devel \ - lapack-devel \ lcov \ libcurl-devel \ m4 \ @@ -38,6 +37,7 @@ RUN cd musica \ && cmake -S . \ -B build \ -D CMAKE_BUILD_TYPE=${BUILD_TYPE} \ + -D MUSICA_GIT_TAG=${MUSICA_GIT_TAG} \ -D MUSICA_BUILD_FORTRAN_INTERFACE=ON \ -D MUSICA_ENABLE_MEMCHECK=ON \ && cd build \ diff --git a/docker/Dockerfile.fortran-gcc.integration b/docker/Dockerfile.fortran-gcc.integration index 5e7204ca..610eebc8 100644 --- a/docker/Dockerfile.fortran-gcc.integration +++ b/docker/Dockerfile.fortran-gcc.integration @@ -14,7 +14,6 @@ RUN dnf -y update \ git \ hdf5-devel \ json-devel \ - lapack-devel \ lcov \ libcurl-devel \ m4 \ diff --git a/docker/Dockerfile.memcheck b/docker/Dockerfile.memcheck index cd943615..8b330b51 100644 --- a/docker/Dockerfile.memcheck +++ b/docker/Dockerfile.memcheck @@ -9,7 +9,6 @@ RUN dnf -y update \ gfortran \ gdb \ git \ - lapack-devel \ lcov \ make \ netcdf-fortran-devel \ diff --git a/docker/Dockerfile.mpi b/docker/Dockerfile.mpi index fa85b48f..8935ca62 100644 --- a/docker/Dockerfile.mpi +++ b/docker/Dockerfile.mpi @@ -16,7 +16,6 @@ RUN sudo dnf -y install \ gcc-c++ \ gfortran \ git \ - lapack-devel \ lcov \ make \ netcdf-fortran-devel \ diff --git a/docker/Dockerfile.mpi_openmp b/docker/Dockerfile.mpi_openmp index f48a4605..7a5eb9b8 100644 --- a/docker/Dockerfile.mpi_openmp +++ b/docker/Dockerfile.mpi_openmp @@ -16,7 +16,6 @@ RUN sudo dnf -y install \ gcc-c++ \ gfortran \ git \ - lapack-devel \ lcov \ make \ netcdf-fortran-devel \ diff --git a/docker/Dockerfile.openmp b/docker/Dockerfile.openmp index efc433b8..cb612509 100644 --- a/docker/Dockerfile.openmp +++ b/docker/Dockerfile.openmp @@ -16,7 +16,6 @@ RUN sudo dnf -y install \ gcc-c++ \ gfortran \ git \ - lapack-devel \ lcov \ make \ netcdf-fortran-devel \ diff --git a/docker/Dockerfile.python b/docker/Dockerfile.python index 02788a0f..8484ddc7 100644 --- a/docker/Dockerfile.python +++ b/docker/Dockerfile.python @@ -9,7 +9,6 @@ RUN dnf -y update \ gcc-fortran \ gdb \ git \ - lapack-devel \ make \ netcdf-fortran-devel \ pip \ diff --git a/fortran/test/fetch_content_integration/CMakeLists.txt b/fortran/test/fetch_content_integration/CMakeLists.txt index 3b978644..d0707f7b 100644 --- a/fortran/test/fetch_content_integration/CMakeLists.txt +++ b/fortran/test/fetch_content_integration/CMakeLists.txt @@ -116,4 +116,4 @@ if (MUSICA_ENABLE_TUVX) copy_tuvx_data_dir ALL ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/../../../build/_deps/tuvx-src/data ${CMAKE_BINARY_DIR}/data ) -endif() +endif() \ No newline at end of file diff --git a/fortran/test/fetch_content_integration/test_tuvx_api.F90 b/fortran/test/fetch_content_integration/test_tuvx_api.F90 index 54b6e349..7dfdddff 100644 --- a/fortran/test/fetch_content_integration/test_tuvx_api.F90 +++ b/fortran/test/fetch_content_integration/test_tuvx_api.F90 @@ -3,7 +3,8 @@ ! program combined_tuvx_tests use iso_c_binding - use musica_tuvx, only: tuvx_t, grid_map_t, grid_t, profile_map_t, profile_t + use musica_tuvx, only: tuvx_t, grid_map_t, grid_t, profile_map_t, profile_t, & + radiator_map_t, radiator_t use musica_util, only: assert, error_t implicit none @@ -55,22 +56,36 @@ end subroutine test_tuvx_api_invalid_config subroutine test_tuvx_solve() - type(tuvx_t), pointer :: tuvx - type(error_t) :: error - type(grid_map_t), pointer :: grids - character(len=256) :: config_path - type(grid_t), pointer :: grid - type(profile_map_t), pointer :: profiles - type(profile_t), pointer :: profile, profile_copy - real*8, dimension(5), target :: edges, edge_values, temp_edge - real*8, dimension(4), target :: midpoints, midpoint_values, layer_densities, temp_midpoint - real*8 :: temp_real + type(tuvx_t), pointer :: tuvx + type(error_t) :: error + character(len=256) :: config_path + type(grid_map_t), pointer :: grids + type(grid_t), pointer :: grid, height_grid, wavelength_grid + type(profile_map_t), pointer :: profiles + type(profile_t), pointer :: profile, profile_copy + type(radiator_map_t), pointer :: radiators + type(radiator_t), pointer :: radiator, radiator_copy + real*8, dimension(5), target :: edges, edge_values, temp_edge + real*8, dimension(4), target :: midpoints, midpoint_values, layer_densities, temp_midpoint + real*8 :: temp_real + integer :: num_vertical_layers, num_wavelength_bins + real*8, dimension(3,2), target :: optical_depths, temp_od + real*8, dimension(3,2), target :: single_scattering_albedos, temp_ssa + real*8, dimension(3,2,1), target :: asymmetry_factors, temp_asym edges = (/ 1.0, 2.0, 3.0, 4.0, 5.0 /) midpoints = (/ 15.0, 25.0, 35.0, 45.0 /) edge_values = (/ 10.0, 20.0, 30.0, 40.0, 50.0 /) midpoint_values = (/ 15.0, 25.0, 35.0, 45.0 /) layer_densities = (/ 2.0, 4.0, 1.0, 7.0 /) + num_vertical_layers = 3 + num_wavelength_bins = 2 + optical_depths(:,1) = (/ 30.0, 20.0, 10.0 /) + optical_depths(:,2) = (/ 70.0, 80.0, 90.0 /) + single_scattering_albedos(:,1) = (/ 300.0, 200.0, 100.0 /) + single_scattering_albedos(:,2) = (/ 700.0, 800.0, 900.0 /) + asymmetry_factors(:,1,1) = (/ 3.0, 2.0, 1.0 /) + asymmetry_factors(:,2,1) = (/ 7.0, 8.0, 9.0 /) config_path = "examples/ts1_tsmlt.json" @@ -312,12 +327,146 @@ subroutine test_tuvx_solve() ASSERT_EQ( temp_edge(3), 36.0 ) ASSERT_EQ( temp_edge(4), 38.0 ) ASSERT_EQ( temp_edge(5), 40.0 ) - + + radiators => tuvx%get_radiators( error ) + ASSERT( error%is_success() ) + + radiator => radiators%get( "foo_radiator", error ) + ASSERT( .not. error%is_success() ) + deallocate( radiator ) + deallocate( radiators ) + + radiators =>radiator_map_t( error ) + ASSERT( error%is_success() ) + + height_grid => grid_t( "height", "km", num_vertical_layers, error ) + wavelength_grid => grid_t( "wavelength", "nm", num_wavelength_bins, error ) + radiator => radiator_t( "foo_radiator", height_grid, wavelength_grid, error ) + ASSERT( error%is_success() ) + + call radiator%set_optical_depths( optical_depths, error ) + ASSERT( error%is_success() ) + + call radiator%get_optical_depths( temp_od, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( temp_od(1,1), 30.0 ) + ASSERT_EQ( temp_od(2,1), 20.0 ) + ASSERT_EQ( temp_od(3,1), 10.0 ) + ASSERT_EQ( temp_od(1,2), 70.0 ) + ASSERT_EQ( temp_od(2,2), 80.0 ) + ASSERT_EQ( temp_od(3,2), 90.0 ) + + call radiator%set_single_scattering_albedos( single_scattering_albedos, error ) + ASSERT( error%is_success() ) + + call radiator%get_single_scattering_albedos( temp_ssa, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( temp_ssa(1,1), 300.0 ) + ASSERT_EQ( temp_ssa(2,1), 200.0 ) + ASSERT_EQ( temp_ssa(3,1), 100.0 ) + ASSERT_EQ( temp_ssa(1,2), 700.0 ) + ASSERT_EQ( temp_ssa(2,2), 800.0 ) + ASSERT_EQ( temp_ssa(3,2), 900.0 ) + + call radiator%set_asymmetry_factors( asymmetry_factors, error ) + ASSERT( error%is_success() ) + + call radiator%get_asymmetry_factors( temp_asym, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( temp_asym(1,1,1), 3.0 ) + ASSERT_EQ( temp_asym(2,1,1), 2.0 ) + ASSERT_EQ( temp_asym(3,1,1), 1.0 ) + ASSERT_EQ( temp_asym(1,2,1), 7.0 ) + ASSERT_EQ( temp_asym(2,2,1), 8.0 ) + ASSERT_EQ( temp_asym(3,2,1), 9.0 ) +! + call radiators%add( radiator, error ) + radiator_copy => radiators%get( "foo_radiator", error ) + + optical_depths(:,:) = 0.0 + single_scattering_albedos(:,:) = 0.0 + asymmetry_factors(:,:,:) = 0.0 + + call radiator_copy%get_optical_depths( optical_depths, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( optical_depths(1,1), 30.0 ) + ASSERT_EQ( optical_depths(2,1), 20.0 ) + ASSERT_EQ( optical_depths(3,1), 10.0 ) + ASSERT_EQ( optical_depths(1,2), 70.0 ) + ASSERT_EQ( optical_depths(2,2), 80.0 ) + ASSERT_EQ( optical_depths(3,2), 90.0 ) + + call radiator_copy%get_single_scattering_albedos( single_scattering_albedos, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( single_scattering_albedos(1,1), 300.0 ) + ASSERT_EQ( single_scattering_albedos(2,1), 200.0 ) + ASSERT_EQ( single_scattering_albedos(3,1), 100.0 ) + ASSERT_EQ( single_scattering_albedos(1,2), 700.0 ) + ASSERT_EQ( single_scattering_albedos(2,2), 800.0 ) + ASSERT_EQ( single_scattering_albedos(3,2), 900.0 ) + + call radiator_copy%get_asymmetry_factors( asymmetry_factors, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( asymmetry_factors(1,1,1), 3.0 ) + ASSERT_EQ( asymmetry_factors(2,1,1), 2.0 ) + ASSERT_EQ( asymmetry_factors(3,1,1), 1.0 ) + ASSERT_EQ( asymmetry_factors(1,2,1), 7.0 ) + ASSERT_EQ( asymmetry_factors(2,2,1), 8.0 ) + ASSERT_EQ( asymmetry_factors(3,2,1), 9.0 ) + + optical_depths(:,1) = (/ 90.0, 80.0, 70.0 /) + optical_depths(:,2) = (/ 75.0, 85.0, 95.0 /) + single_scattering_albedos(:,1) = (/ 900.0, 800.0, 700.0 /) + single_scattering_albedos(:,2) = (/ 750.0, 850.0, 950.0 /) + asymmetry_factors(:,1,1) = (/ 9.0, 8.0, 7.0 /) + asymmetry_factors(:,2,1) = (/ 5.0, 4.0, 3.0 /) + + call radiator_copy%set_optical_depths( optical_depths, error ) + call radiator_copy%set_single_scattering_albedos( single_scattering_albedos, error ) + call radiator_copy%set_asymmetry_factors( asymmetry_factors, error ) + ASSERT( error%is_success() ) + + optical_depths(:,:) = 0.0 + single_scattering_albedos(:,:) = 0.0 + asymmetry_factors(:,:,:) = 0.0 + + call radiator%get_optical_depths( optical_depths, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( optical_depths(1,1), 90.0 ) + ASSERT_EQ( optical_depths(2,1), 80.0 ) + ASSERT_EQ( optical_depths(3,1), 70.0 ) + ASSERT_EQ( optical_depths(1,2), 75.0 ) + ASSERT_EQ( optical_depths(2,2), 85.0 ) + ASSERT_EQ( optical_depths(3,2), 95.0 ) + + call radiator%get_single_scattering_albedos( single_scattering_albedos, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( single_scattering_albedos(1,1), 900.0 ) + ASSERT_EQ( single_scattering_albedos(2,1), 800.0 ) + ASSERT_EQ( single_scattering_albedos(3,1), 700.0 ) + ASSERT_EQ( single_scattering_albedos(1,2), 750.0 ) + ASSERT_EQ( single_scattering_albedos(2,2), 850.0 ) + ASSERT_EQ( single_scattering_albedos(3,2), 950.0 ) + + call radiator%get_asymmetry_factors( asymmetry_factors, error ) + ASSERT( error%is_success() ) + ASSERT_EQ( asymmetry_factors(1,1,1), 9.0 ) + ASSERT_EQ( asymmetry_factors(2,1,1), 8.0 ) + ASSERT_EQ( asymmetry_factors(3,1,1), 7.0 ) + ASSERT_EQ( asymmetry_factors(1,2,1), 5.0 ) + ASSERT_EQ( asymmetry_factors(2,2,1), 4.0 ) + ASSERT_EQ( asymmetry_factors(3,2,1), 3.0 ) + deallocate( grid ) deallocate( grids ) deallocate( profile ) deallocate( profile_copy ) deallocate( profiles ) + deallocate( radiator_copy ) + deallocate( radiator ) + deallocate( radiators ) + deallocate( height_grid ) + deallocate( wavelength_grid ) deallocate( tuvx ) end subroutine test_tuvx_solve diff --git a/fortran/tuvx/CMakeLists.txt b/fortran/tuvx/CMakeLists.txt index cb3eeb6c..68e400c9 100644 --- a/fortran/tuvx/CMakeLists.txt +++ b/fortran/tuvx/CMakeLists.txt @@ -4,5 +4,7 @@ target_sources(musica-fortran grid_map.F90 profile.F90 profile_map.F90 + radiator.F90 + radiator_map.F90 tuvx.F90 ) \ No newline at end of file diff --git a/fortran/tuvx/grid.F90 b/fortran/tuvx/grid.F90 index b0a99a04..202b1b94 100644 --- a/fortran/tuvx/grid.F90 +++ b/fortran/tuvx/grid.F90 @@ -82,7 +82,7 @@ end subroutine get_grid_midpoints_c procedure :: set_edges ! Get grid edges procedure :: get_edges - ! Set the grid edges and midpoints + ! Set the grid midpoints procedure :: set_midpoints ! Get the grid midpoints procedure :: get_midpoints diff --git a/fortran/tuvx/profile.F90 b/fortran/tuvx/profile.F90 index 68c8bd6a..757a08fc 100644 --- a/fortran/tuvx/profile.F90 +++ b/fortran/tuvx/profile.F90 @@ -16,7 +16,7 @@ module musica_tuvx_profile interface function create_profile_c(profile_name, profile_units, grid, error) & bind(C, name="CreateProfile") - use iso_c_binding, only: c_ptr, c_char, c_size_t + use iso_c_binding, only: c_ptr, c_char use musica_util, only: error_t_c character(len=1, kind=c_char), intent(in) :: profile_name(*) character(len=1, kind=c_char), intent(in) :: profile_units(*) diff --git a/fortran/tuvx/radiator.F90 b/fortran/tuvx/radiator.F90 new file mode 100644 index 00000000..60a79ee6 --- /dev/null +++ b/fortran/tuvx/radiator.F90 @@ -0,0 +1,361 @@ +! Copyright (C) 2023-2024 National Center for Atmospheric Research +! SPDX-License-Identifier: Apache-2.0 +! +module musica_tuvx_radiator + use iso_c_binding, only: c_ptr, c_null_ptr + + implicit none + +#define ASSERT( expr ) call assert( expr, __FILE__, __LINE__ ) + + private + public :: radiator_t + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + interface + function create_radiator_c(radiator_name, height_grid, wavelength_grid, error) & + bind(C, name="CreateRadiator") + use iso_c_binding, only : c_ptr, c_char + use musica_util, only: error_t_c + character(len=1, kind=c_char), intent(in) :: radiator_name(*) + type(c_ptr), value, intent(in) :: height_grid + type(c_ptr), value, intent(in) :: wavelength_grid + type(error_t_c), intent(inout) :: error + type(c_ptr) :: create_radiator_c + end function create_radiator_c + + subroutine delete_radiator_c(radiator, error) bind(C, name="DeleteRadiator") + use iso_c_binding, only : c_ptr + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(error_t_c), intent(inout) :: error + end subroutine delete_radiator_c + + subroutine set_optical_depths_c(radiator, optical_depths, num_vertical_layers, & + num_wavelength_bins, error) bind(C, name="SetRadiatorOpticalDepths") + use iso_c_binding, only : c_ptr, c_size_t + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(c_ptr), value, intent(in) :: optical_depths + integer(c_size_t), value, intent(in) :: num_vertical_layers + integer(c_size_t), value, intent(in) :: num_wavelength_bins + type(error_t_c), intent(inout) :: error + end subroutine set_optical_depths_c + + subroutine get_optical_depths_c(radiator, optical_depths, num_vertical_layers, & + num_wavelength_bins, error) bind(C, name="GetRadiatorOpticalDepths") + use iso_c_binding, only : c_ptr, c_size_t + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(c_ptr), value, intent(in) :: optical_depths + integer(c_size_t), value, intent(in) :: num_vertical_layers + integer(c_size_t), value, intent(in) :: num_wavelength_bins + type(error_t_c), intent(inout) :: error + end subroutine get_optical_depths_c + + subroutine set_single_scattering_albedos_c(radiator, single_scattering_albedos, & + num_vertical_layers, num_wavelength_bins, error) & + bind(C, name="SetRadiatorSingleScatteringAlbedos") + use iso_c_binding, only : c_ptr, c_size_t + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(c_ptr), value, intent(in) :: single_scattering_albedos + integer(c_size_t), value, intent(in) :: num_vertical_layers + integer(c_size_t), value, intent(in) :: num_wavelength_bins + type(error_t_c), intent(inout) :: error + end subroutine set_single_scattering_albedos_c + + subroutine get_single_scattering_albedos_c(radiator, single_scattering_albedos, & + num_vertical_layers, num_wavelength_bins, error) & + bind(C, name="GetRadiatorSingleScatteringAlbedos") + use iso_c_binding, only : c_ptr, c_size_t + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(c_ptr), value, intent(in) :: single_scattering_albedos + integer(c_size_t), value, intent(in) :: num_vertical_layers + integer(c_size_t), value, intent(in) :: num_wavelength_bins + type(error_t_c), intent(inout) :: error + end subroutine get_single_scattering_albedos_c + + subroutine set_asymmetry_factors_c(radiator, asymmetry_factors, num_vertical_layers, & + num_wavelength_bins, num_streams, error) bind(C, name="SetRadiatorAsymmetryFactors") + use iso_c_binding, only : c_ptr, c_size_t + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(c_ptr), value, intent(in) :: asymmetry_factors + integer(c_size_t), value, intent(in) :: num_vertical_layers + integer(c_size_t), value, intent(in) :: num_wavelength_bins + integer(c_size_t), value, intent(in) :: num_streams + type(error_t_c), intent(inout) :: error + end subroutine set_asymmetry_factors_c + + subroutine get_asymmetry_factors_c(radiator, asymmetry_factors, num_vertical_layers, & + num_wavelength_bins, num_streams, error) bind(C, name="GetRadiatorAsymmetryFactors") + use iso_c_binding, only : c_ptr, c_size_t + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator + type(c_ptr), value, intent(in) :: asymmetry_factors + integer(c_size_t), value, intent(in) :: num_vertical_layers + integer(c_size_t), value, intent(in) :: num_wavelength_bins + integer(c_size_t), value, intent(in) :: num_streams + type(error_t_c), intent(inout) :: error + end subroutine get_asymmetry_factors_c + end interface + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + type :: radiator_t + type(c_ptr) :: ptr_ = c_null_ptr + contains + ! Set radiator optical depths + procedure :: set_optical_depths + ! Get radiator optical depths + procedure :: get_optical_depths + ! Set radiator single scattering albedos + procedure :: set_single_scattering_albedos + ! Get radiator single scattering albedos + procedure :: get_single_scattering_albedos + ! Set radiator asymmetry_factors + procedure :: set_asymmetry_factors + ! Get radiator asymmetry factors + procedure :: get_asymmetry_factors + ! Deallocate radiator instance + final :: finalize_radiator_t + end type radiator_t + + interface radiator_t + procedure radiator_t_ptr_constructor + procedure radiator_t_constructor + end interface radiator_t + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +contains + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Constructs a radiator instance that wraps an existing TUV-x radiator + function radiator_t_ptr_constructor(radiator_c_ptr) result(this) + ! Arguments + type(c_ptr), intent(in) :: radiator_c_ptr + + ! Return value + type(radiator_t), pointer :: this + + allocate( this ) + this%ptr_ = radiator_c_ptr + + end function radiator_t_ptr_constructor + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Constructs a radiator instance that allocates a new TUV-x radiator + function radiator_t_constructor(radiator_name, height_grid, wavelength_grid, & + error) result(this) + use musica_tuvx_grid, only: grid_t + use musica_util, only: error_t, error_t_c, to_c_string + + ! Arguments + character(len=*), intent(in) :: radiator_name + type(grid_t), intent(in) :: height_grid + type(grid_t), intent(in) :: wavelength_grid + type(error_t), intent(inout) :: error + + ! Return value + type(radiator_t), pointer :: this + + ! Local variables + type(error_t_c) :: error_c + + allocate( this ) + this%ptr_ = create_radiator_c(to_c_string(radiator_name), height_grid%ptr_, & + wavelength_grid%ptr_, error_c) + error = error_t(error_c) + + end function radiator_t_constructor + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine set_optical_depths(this, optical_depths, error) + use iso_c_binding, only: c_size_t, c_loc + use musica_util, only: error_t, error_t_c, dk => musica_dk + + ! Arguments + class(radiator_t), intent(inout) :: this + real(dk), target, dimension(:,:), intent(in) :: optical_depths + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + integer(kind=c_size_t) :: num_vertical_layers + integer(kind=c_size_t) :: num_wavelength_bins + + num_vertical_layers = size(optical_depths, 1) + num_wavelength_bins = size(optical_depths, 2) + + call set_optical_depths_c(this%ptr_, c_loc(optical_depths), & + num_vertical_layers, num_wavelength_bins, error_c) + error = error_t(error_c) + + end subroutine set_optical_depths + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine get_optical_depths(this, optical_depths, error) + use iso_c_binding, only: c_size_t, c_loc + use musica_util, only: error_t, error_t_c, dk => musica_dk + + ! Arguments + class(radiator_t), intent(inout) :: this + real(dk), target, dimension(:,:), intent(in) :: optical_depths + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + integer(kind=c_size_t) :: num_vertical_layers + integer(kind=c_size_t) :: num_wavelength_bins + + num_vertical_layers = size(optical_depths, 1) + num_wavelength_bins = size(optical_depths, 2) + + call get_optical_depths_c(this%ptr_, c_loc(optical_depths), & + num_vertical_layers, num_wavelength_bins, error_c) + error = error_t(error_c) + + end subroutine get_optical_depths + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine set_single_scattering_albedos(this, single_scattering_albedos, & + error) + use iso_c_binding, only: c_size_t, c_loc + use musica_util, only: error_t, error_t_c, dk => musica_dk + + ! Arguments + class(radiator_t), intent(inout) :: this + real(dk), target, dimension(:,:), intent(in) :: single_scattering_albedos + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + integer(kind=c_size_t) :: num_vertical_layers + integer(kind=c_size_t) :: num_wavelength_bins + + num_vertical_layers = size(single_scattering_albedos, 1) + num_wavelength_bins = size(single_scattering_albedos, 2) + + call set_single_scattering_albedos_c(this%ptr_, & + c_loc(single_scattering_albedos), num_vertical_layers, & + num_wavelength_bins, error_c) + error = error_t(error_c) + + end subroutine set_single_scattering_albedos + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine get_single_scattering_albedos(this, single_scattering_albedos, & + error) + use iso_c_binding, only: c_size_t, c_loc + use musica_util, only: error_t, error_t_c, dk => musica_dk + + ! Arguments + class(radiator_t), intent(inout) :: this + real(dk), target, dimension(:,:), intent(in) :: single_scattering_albedos + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + integer(kind=c_size_t) :: num_vertical_layers + integer(kind=c_size_t) :: num_wavelength_bins + + num_vertical_layers = size(single_scattering_albedos, 1) + num_wavelength_bins = size(single_scattering_albedos, 2) + + call get_single_scattering_albedos_c(this%ptr_, & + c_loc(single_scattering_albedos), num_vertical_layers, & + num_wavelength_bins, error_c) + error = error_t(error_c) + + end subroutine get_single_scattering_albedos + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine set_asymmetry_factors(this, asymmetry_factors, error) + use iso_c_binding, only: c_size_t, c_loc + use musica_util, only: error_t, error_t_c, dk => musica_dk + + ! Arguments + class(radiator_t), intent(inout) :: this + real(dk), target, dimension(:,:,:), intent(in) :: asymmetry_factors + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + integer(kind=c_size_t) :: num_vertical_layers + integer(kind=c_size_t) :: num_wavelength_bins + integer(kind=c_size_t) :: num_streams + + num_vertical_layers = size(asymmetry_factors, 1) + num_wavelength_bins = size(asymmetry_factors, 2) + num_streams = size(asymmetry_factors, 3) + + call set_asymmetry_factors_c(this%ptr_, c_loc(asymmetry_factors), & + num_vertical_layers, num_wavelength_bins, num_streams, error_c) + error = error_t(error_c) + +end subroutine set_asymmetry_factors + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine get_asymmetry_factors(this, asymmetry_factors, error) + use iso_c_binding, only: c_size_t, c_loc + use musica_util, only: error_t, error_t_c, dk => musica_dk + + ! Arguments + class(radiator_t), intent(inout) :: this + real(dk), target, dimension(:,:,:), intent(in) :: asymmetry_factors + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + integer(kind=c_size_t) :: num_vertical_layers + integer(kind=c_size_t) :: num_wavelength_bins + integer(kind=c_size_t) :: num_streams + + num_vertical_layers = size(asymmetry_factors, 1) + num_wavelength_bins = size(asymmetry_factors, 2) + num_streams = size(asymmetry_factors, 3) + + call get_asymmetry_factors_c(this%ptr_, c_loc(asymmetry_factors), & + num_vertical_layers, num_wavelength_bins, num_streams, error_c) + error = error_t(error_c) + + end subroutine get_asymmetry_factors + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Deallocate the radiator instance + subroutine finalize_radiator_t(this) + use iso_c_binding, only: c_associated + use musica_util, only: error_t, error_t_c, assert + + ! Arguments + type(radiator_t), intent(inout) :: this + + ! Local variables + type(error_t_c) :: error_c + type(error_t) :: error + + if (c_associated(this%ptr_)) then + call delete_radiator_c(this%ptr_, error_c) + this%ptr_ = c_null_ptr + error = error_t(error_c) + ASSERT(error%is_success()) + end if + + end subroutine finalize_radiator_t + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +end module musica_tuvx_radiator \ No newline at end of file diff --git a/fortran/tuvx/radiator_map.F90 b/fortran/tuvx/radiator_map.F90 new file mode 100644 index 00000000..f9cbdd61 --- /dev/null +++ b/fortran/tuvx/radiator_map.F90 @@ -0,0 +1,180 @@ +! Copyright (C) 2023-2024 National Center for Atmospheric Research +! SPDX-License-Identifier: Apache-2.0 +! +module musica_tuvx_radiator_map + use iso_c_binding, only: c_ptr, c_null_ptr + + implicit none + +#define ASSERT( expr ) call assert( expr, __FILE__, __LINE__ ) + + private + public :: radiator_map_t + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + interface + function create_radiator_map_c(error) bind(C, name="CreateRadiatorMap") + use iso_c_binding, only: c_ptr + use musica_util, only: error_t_c + type(error_t_c), intent(inout) :: error + type(c_ptr) :: create_radiator_map_c + end function create_radiator_map_c + + subroutine delete_radiator_map_c(radiator_map, error) & + bind(C, name="DeleteRadiatorMap") + use iso_c_binding, only: c_ptr + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator_map + type(error_t_c), intent(inout) :: error + end subroutine delete_radiator_map_c + + subroutine add_radiator_c(radiator_map, radiator, error) & + bind(C, name="AddRadiator") + use iso_c_binding, only: c_ptr + use musica_util, only: error_t_c + type(c_ptr), value, intent(in) :: radiator_map + type(c_ptr), value, intent(in) :: radiator + type(error_t_c), intent(inout) :: error + end subroutine add_radiator_c + + function get_radiator_c(radiator_map, radiator_name, error) & + bind(C, name="GetRadiator") + use musica_util, only: error_t_c + use iso_c_binding, only: c_ptr, c_char + type(c_ptr), value, intent(in) :: radiator_map + character(len=1, kind=c_char), intent(in) :: radiator_name(*) + type(error_t_c), intent(inout) :: error + type(c_ptr) :: get_radiator_c + end function get_radiator_c + end interface + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + type :: radiator_map_t + type(c_ptr) :: ptr_ = c_null_ptr + contains + ! Adds a radiator to the radiator map + procedure :: add => add_radiator + ! Get a radiator given its name + procedure :: get => get_radiator + ! Deallocate the radiator map instance + final :: finalize_radiator_map_t + end type radiator_map_t + + interface radiator_map_t + procedure radiator_map_t_ptr_constructor + procedure radiator_map_t_constructor + end interface radiator_map_t + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +contains + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Wraps an existing radiator map + function radiator_map_t_ptr_constructor(radiator_map_c_ptr) result(this) + ! Arguments + type(c_ptr), intent(in) :: radiator_map_c_ptr + ! Return value + type(radiator_map_t), pointer :: this + + allocate( this ) + this%ptr_ = radiator_map_c_ptr + + end function radiator_map_t_ptr_constructor + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Creates a new radiator map + function radiator_map_t_constructor(error) result(this) + use musica_util, only: error_t, error_t_c, assert + + ! Arguments + type(error_t), intent(inout) :: error + + ! Return value + type(radiator_map_t), pointer :: this + + ! Local variables + type(error_t_c) :: error_c + + allocate( this ) + this%ptr_ = create_radiator_map_c(error_c) + error = error_t(error_c) + + end function radiator_map_t_constructor + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Adds a radiator to a radiator map + subroutine add_radiator(this, radiator, error) + use musica_tuvx_radiator, only: radiator_t + use musica_util, only: error_t, error_t_c, assert + + ! Arguments + class(radiator_map_t), intent(inout) :: this + type(radiator_t), intent(in) :: radiator + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + + call add_radiator_c(this%ptr_, radiator%ptr_, error_c) + error = error_t(error_c) + + end subroutine add_radiator + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Gets a radiator given its name + function get_radiator(this, radiator_name, error) result(radiator) + use iso_c_binding, only: c_char + use musica_tuvx_radiator, only : radiator_t + use musica_util, only: error_t, error_t_c, to_c_string + + ! Arguments + class(radiator_map_t), intent(in) :: this + character(len=*), intent(in) :: radiator_name + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + + ! Return value + type(radiator_t), pointer :: radiator + + radiator => radiator_t(get_radiator_c(this%ptr_, to_c_string(radiator_name), & + error_c)) + + error = error_t(error_c) + + end function get_radiator + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + !> Deallocates the radiator map instance + subroutine finalize_radiator_map_t(this) + use iso_c_binding, only: c_associated + use musica_util, only: error_t, error_t_c, assert + + ! Arguments + type(radiator_map_t), intent(inout) :: this + + ! Local variables + type(error_t_c) :: error_c + type(error_t) :: error + + if (c_associated(this%ptr_)) then + call delete_radiator_map_c(this%ptr_, error_c) + this%ptr_ = c_null_ptr + error = error_t(error_c) + ASSERT(error%is_success()) + end if + + end subroutine finalize_radiator_map_t + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +end module musica_tuvx_radiator_map diff --git a/fortran/tuvx/tuvx.F90 b/fortran/tuvx/tuvx.F90 index 94e4edd0..36a58d79 100644 --- a/fortran/tuvx/tuvx.F90 +++ b/fortran/tuvx/tuvx.F90 @@ -2,18 +2,21 @@ ! SPDX-License-Identifier: Apache-2.0 ! module musica_tuvx - use iso_c_binding, only: c_ptr, c_null_ptr - use musica_tuvx_grid, only : grid_t - use musica_tuvx_grid_map, only : grid_map_t - use musica_tuvx_profile, only : profile_t - use musica_tuvx_profile_map, only : profile_map_t + use iso_c_binding, only: c_ptr, c_null_ptr + use musica_tuvx_grid, only : grid_t + use musica_tuvx_grid_map, only : grid_map_t + use musica_tuvx_profile, only : profile_t + use musica_tuvx_profile_map, only : profile_map_t + use musica_tuvx_radiator, only : radiator_t + use musica_tuvx_radiator_map, only : radiator_map_t implicit none #define ASSERT( expr ) call assert( expr, __FILE__, __LINE__ ) private - public :: tuvx_t, grid_map_t, grid_t, profile_map_t, profile_t + public :: tuvx_t, grid_map_t, grid_t, profile_map_t, profile_t, & + radiator_map_t, radiator_t !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -45,9 +48,17 @@ function get_profile_map_c(tuvx, error) bind(C, name="GetProfileMap") use musica_util, only: error_t_c use iso_c_binding, only: c_ptr type(c_ptr), value, intent(in) :: tuvx - type(error_t_c), intent(inout) :: error - type(c_ptr) :: get_profile_map_c + type(error_t_c), intent(inout) :: error + type(c_ptr) :: get_profile_map_c end function get_profile_map_c + + function get_radiator_map_c(tuvx, error) bind(C, name="GetRadiatorMap") + use musica_util, only: error_t_c + use iso_c_binding, only: c_ptr + type(c_ptr), value, intent(in) :: tuvx + type(error_t_c), intent(inout) :: error + type(c_ptr) :: get_radiator_map_c + end function get_radiator_map_c end interface !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -59,6 +70,8 @@ end function get_profile_map_c procedure :: get_grids ! Create a profile map procedure :: get_profiles + ! Create a radiator map + procedure :: get_radiators ! Deallocate the tuvx instance final :: finalize end type tuvx_t @@ -154,6 +167,28 @@ end function get_profiles !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !> Get the radiator map + function get_radiators(this, error) result(radiator_map) + use musica_util, only: error_t, error_t_c + + ! Arguments + class(tuvx_t), intent(inout) :: this + type(error_t), intent(inout) :: error + + ! Local variables + type(error_t_c) :: error_c + + ! Return value + type(radiator_map_t), pointer :: radiator_map + + radiator_map => radiator_map_t(get_radiator_map_c(this%ptr_, error_c)) + + error = error_t(error_c) + + end function get_radiators + + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !> Deallocate the tuvx instance subroutine finalize(this) use musica_util, only: error_t, error_t_c, assert diff --git a/include/musica/tuvx/grid.hpp b/include/musica/tuvx/grid.hpp index 3d52174e..bd7a6d3b 100644 --- a/include/musica/tuvx/grid.hpp +++ b/include/musica/tuvx/grid.hpp @@ -2,7 +2,6 @@ // SPDX-License-Identifier: Apache-2.0 #pragma once -#include #include #include @@ -13,10 +12,12 @@ namespace musica { class GridMap; class Profile; + class Radiator; - /// @brief A grid struct used to access grid information in tuvx - struct Grid + /// @brief A grid class used to access grid information in tuvx + class Grid { + public: /// @brief Creates a grid instance /// @param grid_name The name of the grid /// @param units The units of the grid @@ -56,6 +57,7 @@ namespace musica friend class GridMap; friend class Profile; + friend class Radiator; /// @brief Wraps an existing grid instance. Used by GridMap /// @param updater The updater for the grid diff --git a/include/musica/tuvx/grid_map.hpp b/include/musica/tuvx/grid_map.hpp index 48a02e01..50686593 100644 --- a/include/musica/tuvx/grid_map.hpp +++ b/include/musica/tuvx/grid_map.hpp @@ -3,7 +3,6 @@ #pragma once #include -#include #include #include @@ -13,9 +12,10 @@ namespace musica { - /// @brief A grid map struct used to access grid information in tuvx - struct GridMap + /// @brief A grid map class used to access grid information in tuvx + class GridMap { + public: GridMap(void *grid_map) : grid_map_(grid_map), owns_grid_map_(false) diff --git a/include/musica/tuvx/profile.hpp b/include/musica/tuvx/profile.hpp index 147d970a..309f9bbc 100644 --- a/include/musica/tuvx/profile.hpp +++ b/include/musica/tuvx/profile.hpp @@ -14,9 +14,10 @@ namespace musica { class ProfileMap; - /// @brief A struct used to interact with TUV-x profiles (properties with values on a grid) - struct Profile + /// @brief A class used to interact with TUV-x profiles (properties with values on a grid) + class Profile { + public: /// @brief Creates a profile instance /// @param profile_name The name of the profile /// @param units The units of the profile diff --git a/include/musica/tuvx/profile_map.hpp b/include/musica/tuvx/profile_map.hpp index 35fc8857..0b5ee6fe 100644 --- a/include/musica/tuvx/profile_map.hpp +++ b/include/musica/tuvx/profile_map.hpp @@ -13,9 +13,10 @@ namespace musica { - /// @brief A struct used to store a collection of profiles - struct ProfileMap + /// @brief A class used to store a collection of profiles + class ProfileMap { + public: ProfileMap(void *profile_map) : profile_map_(profile_map), owns_profile_map_(false) diff --git a/include/musica/tuvx/radiator.hpp b/include/musica/tuvx/radiator.hpp new file mode 100644 index 00000000..f08b9b33 --- /dev/null +++ b/include/musica/tuvx/radiator.hpp @@ -0,0 +1,263 @@ +// Copyright (C) 2023-2024 National Center for Atmospheric Research +// SPDX-License-Identifier: Apache-2.0 +#pragma once + +#include +#include + +#include + +namespace musica +{ + class RadiatorMap; + + /// @brief Radiator class used to access radiator information in tuvx + class Radiator + { + public: + /// @brief Creates radiator + /// @param radiator_name Radiator name + /// @param height_grid Height grid + /// @param wavelength_grid Wavelength grid + /// @param error Error to indicate success or failure + Radiator(const char *radiator_name, Grid *height_grid, Grid *wavelength_grid, Error *error); + + ~Radiator(); + + /// @brief Sets optical depth values + /// @param optical_depths 2D array of optical depth values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void + SetOpticalDepths(double *optical_depths, std::size_t num_vertical_layers, std::size_t num_wavelength_bins, Error *error); + + /// @brief Gets optical depth values + /// @param optical_depths 2D array of optical depth values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void + GetOpticalDepths(double *optical_depths, std::size_t num_vertical_layers, std::size_t num_wavelength_bins, Error *error); + + /// @brief Sets single scattering albedos values + /// @param single_scattering_albedos 2D array of single scattering albedos values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void SetSingleScatteringAlbedos( + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error); + + /// @brief Gets single scattering albedos values + /// @param single_scattering_albedos 2D array of single scattering albedos values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void GetSingleScatteringAlbedos( + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error); + + /// @brief Sets asymmetry factor values + /// @param asymmetry_factor 3D array of asymmetery factor values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param num_streams Number of streams + /// @param error Error to indicate success or failure + void SetAsymmetryFactors( + double *asymmetry_factor, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error); + + /// @brief Gets asymmetry factor values + /// @param asymmetry_factor 3D array of asymmetery factor values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param num_streams Number of streams + /// @param error Error to indicate success or failure + void GetAsymmetryFactors( + double *asymmetry_factor, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error); + + private: + void *radiator_; // A valid pointer to a radiator instance indicates ownership by this wrapper + void *updater_; + + friend class RadiatorMap; + + /// @brief Wraps an existing radiator instance. Used by RadiatorMap + /// @param updater The updater for the radiator + Radiator(void *updater) + : radiator_(nullptr), + updater_(updater) + { + } + }; + +#ifdef __cplusplus + extern "C" + { +#endif + + // The external C API for TUVX + // callable by wrappers in other languages + + /// @brief Creates radiator + /// @param radiator_name Radiator name + /// @param height_grid Height grid + /// @param wavelength_grid Wavelength grid + /// @param error Error to indicate success or failure + Radiator *CreateRadiator(const char *radiator_name, Grid *height_grid, Grid *wavelength_grid, Error *error); + + /// @brief Deletes radiator + /// @param radiator Radiator + /// @param error Error to indicate success or failure + void DeleteRadiator(Radiator *radiator, Error *error); + + /// @brief Sets optical depth values + /// @param radiator Radiator + /// @param optical_depths 2D array of optical depth values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void SetRadiatorOpticalDepths( + Radiator *radiator, + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error); + + /// @brief Gets optical depth values + /// @param radiator Radiator + /// @param optical_depths 2D array of optical depth values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void GetRadiatorOpticalDepths( + Radiator *radiator, + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error); + + /// @brief Sets single scattering albedos values + /// @param radiator Radiator + /// @param single_scattering_albedos 2D array of single scattering albedos values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void SetRadiatorSingleScatteringAlbedos( + Radiator *radiator, + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error); + + /// @brief Gets single scattering albedos values + /// @param radiator Radiator + /// @param single_scattering_albedos 2D array of single scattering albedos values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param error Error to indicate success or failure + void GetRadiatorSingleScatteringAlbedos( + Radiator *radiator, + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error); + + /// @brief Sets asymmetry factor values + /// @param radiator Radiator + /// @param asymmetry_factor 3D array of asymmetery factor values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param num_streams Number of streams + /// @param error Error to indicate success or failure + void SetRadiatorAsymmetryFactors( + Radiator *radiator, + double *asymmetry_factor, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error); + + /// @brief Gets asymmetry factor values + /// @param radiator Radiator + /// @param asymmetry_factor 3D array of asymmetery factor values + /// @param num_vertical_layers Number of vertical layers + /// @param num_wavelength_bins Number of wavelength bins + /// @param num_streams Number of streams + /// @param error Error to indicate success or failure + void GetRadiatorAsymmetryFactors( + Radiator *radiator, + double *asymmetry_factor, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error); + + // INTERNAL USE. If tuvx ever gets rewritten in C++, these functions will + // go away but the C API will remain the same and downstream projects (like CAM-SIMA) will + // not need to change + void *InternalCreateRadiator( + const char *radiator_name, + std::size_t radiator_name_length, + void *height_grid, + void *wavelength_grid, + int *error_code); + void InternalDeleteRadiator(void *radiator, int *error_code); + void *InternalGetRadiatorUpdater(void *radiator, int *error_code); + void InternalDeleteRadiatorUpdater(void *updater, int *error_code); + void InternalSetOpticalDepths( + void *radiator, + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + int *error_code); + void InternalGetOpticalDepths( + void *radiator, + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + int *error_code); + void InternalSetSingleScatteringAlbedos( + void *radiator, + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + int *error_code); + void InternalGetSingleScatteringAlbedos( + void *radiator, + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + int *error_code); + void InternalSetAsymmetryFactors( + void *radiator, + double *asymmetry_factors, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + int *error_code); + void InternalGetAsymmetryFactors( + void *radiator, + double *asymmetry_factors, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + int *error_code); + +#ifdef __cplusplus + } +#endif + +} // namespace musica \ No newline at end of file diff --git a/include/musica/tuvx/radiator_map.hpp b/include/musica/tuvx/radiator_map.hpp new file mode 100644 index 00000000..c998c9ac --- /dev/null +++ b/include/musica/tuvx/radiator_map.hpp @@ -0,0 +1,93 @@ +// Copyright (C) 2023-2024 National Center for Atmospheric Research +// SPDX-License-Identifier: Apache-2.0 +#pragma once + +#include +#include + +#include +#include +#include + +namespace musica +{ + + /// @brief Radiator map used to access radiator information in tuvx + class RadiatorMap + { + public: + RadiatorMap(void *radiator_map) + : radiator_map_(radiator_map), + owns_radiator_map_(false) + { + } + + /// @brief Creates radiator map + /// @param error Error to indicate success or failure + RadiatorMap(Error *error); + + ~RadiatorMap(); + + /// @brief Adds a radiator to the radiator map + /// @param radiator Radiator to add + /// @param error Error to indicate success or failure + void AddRadiator(Radiator *radiator, Error *error); + + /// @brief Returns a radiator. For now, this calls the interal tuvx fortran api, but will allow the change to c++ later + /// on to be transparent to downstream projects + /// @param radiator_name Radiator name + /// @param error Error to indicate success or failure + /// @return Radiator + Radiator *GetRadiator(const char *radiator_name, Error *error); + + private: + void *radiator_map_; + bool owns_radiator_map_; + }; + +#ifdef __cplusplus + extern "C" + { +#endif + + // The external C API for TUVX + // callable by wrappers in other languages + + /// @brief Creates radiator map + /// @param error Error to indicate success or failure + /// @return Radiator map + RadiatorMap *CreateRadiatorMap(Error *error); + + /// @brief Deletes radiator map + /// @param radiator_map Radiator map to delete + /// @param error Error to indicate success or failure + void DeleteRadiatorMap(RadiatorMap *radiator_map, Error *error); + + /// @brief Adds a radiator to the radiator map + /// @param radiator_map Radiator map to add the radiator to + /// @param radiator Radiator to add + /// @param error Error to indicate success or failure + void AddRadiator(RadiatorMap *radiator_map, Radiator *radiator, Error *error); + + /// @brief Returns a radiator from the radiator map + /// @param radiator_map Radiator map to get the radiator from + /// @param radiator_name Radiator name + /// @param error Error to indicate success or failure + /// @return The radiator pointer, or nullptr if the radiator is not found + Radiator *GetRadiator(RadiatorMap *radiator_map, const char *radiator_name, Error *error); + + // INTERNAL USE. If tuvx ever gets rewritten in C++, these functions will + // go away but the C API will remain the same and downstream projects (like CAM-SIMA) will + // not need to change + void *InternalCreateRadiatorMap(int *error_code); + void InternalDeleteRadiatorMap(void *radiator_map, int *error_code); + void InternalAddRadiator(void *radiator_map, void *radiator, int *error_code); + void * + InternalGetRadiator(void *radiator_map, const char *radiator_name, std::size_t radiator_name_length, int *error_code); + void *InternalGetRadiatorUpdaterFromMap(void *radiator_map, void *radiator, int *error_code); + +#ifdef __cplusplus + } +#endif + +} // namespace musica diff --git a/include/musica/tuvx/tuvx.hpp b/include/musica/tuvx/tuvx.hpp index e89d7685..bd862fde 100644 --- a/include/musica/tuvx/tuvx.hpp +++ b/include/musica/tuvx/tuvx.hpp @@ -7,6 +7,7 @@ #include #include +#include #include #include @@ -38,6 +39,12 @@ namespace musica /// @return a profile map pointer ProfileMap *CreateProfileMap(Error *error); + /// @brief Create a radiator map. For now, this calls the interal tuvx fortran api, but will allow the change to c++ + /// later on to be transparent to downstream projects + /// @param error The error struct to indicate success or failure + /// @return a radiator map pointer + RadiatorMap *CreateRadiatorMap(Error *error); + ~TUVX(); private: @@ -55,6 +62,7 @@ namespace musica void DeleteTuvx(const TUVX *tuvx, Error *error); GridMap *GetGridMap(TUVX *tuvx, Error *error); ProfileMap *GetProfileMap(TUVX *tuvx, Error *error); + RadiatorMap *GetRadiatorMap(TUVX *tuvx, Error *error); // for use by musica interanlly. If tuvx ever gets rewritten in C++, these functions will // go away but the C API will remain the same and downstream projects (like CAM-SIMA) will @@ -63,6 +71,7 @@ namespace musica void InternalDeleteTuvx(void *tuvx, int *error_code); void *InternalGetGridMap(void *tuvx, int *error_code); void *InternalGetProfileMap(void *tuvx, int *error_code); + void *InternalGetRadiatorMap(void *tuvx, int *error_code); #ifdef __cplusplus } diff --git a/src/packaging/CMakeLists.txt b/src/packaging/CMakeLists.txt index 17c27bb3..9d04933f 100644 --- a/src/packaging/CMakeLists.txt +++ b/src/packaging/CMakeLists.txt @@ -76,6 +76,8 @@ if (MUSICA_ENABLE_TUVX) ${MUSICA_FORTRAN_SRC_DIR}/tuvx/grid_map.F90 ${MUSICA_FORTRAN_SRC_DIR}/tuvx/profile.F90 ${MUSICA_FORTRAN_SRC_DIR}/tuvx/profile_map.F90 + ${MUSICA_FORTRAN_SRC_DIR}/tuvx/radiator.F90 + ${MUSICA_FORTRAN_SRC_DIR}/tuvx/radiator_map.F90 ${MUSICA_FORTRAN_SRC_DIR}/tuvx/tuvx.F90 DESTINATION ${MUSICA_INSTALL_INCLUDE_DIR}/musica/fortran diff --git a/src/test/unit/tuvx/tuvx_c_api.cpp b/src/test/unit/tuvx/tuvx_c_api.cpp index b91ae930..91dbbc93 100644 --- a/src/test/unit/tuvx/tuvx_c_api.cpp +++ b/src/test/unit/tuvx/tuvx_c_api.cpp @@ -385,3 +385,340 @@ TEST_F(TuvxCApiTest, CanCreateProfileMap) ASSERT_TRUE(IsSuccess(error)); DeleteError(&error); } + +TEST_F(TuvxCApiTest, CannotGetConfiguredRadiator) +{ + const char* yaml_config_path = "examples/ts1_tsmlt.yml"; + SetUp(yaml_config_path); + Error error; + RadiatorMap* radiator_map = GetRadiatorMap(tuvx, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(radiator_map, nullptr); + Radiator* radiator = GetRadiator(radiator_map, "foo", &error); + ASSERT_FALSE(IsSuccess(error)); // non-host grid + ASSERT_EQ(radiator, nullptr); + DeleteRadiatorMap(radiator_map, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteError(&error); +} + +TEST_F(TuvxCApiTest, CanCreateRadiator) +{ + Error error; + Grid* height = CreateGrid("height", "km", 3, &error); + Grid* wavelength = CreateGrid("wavelength", "nm", 2, &error); + Radiator* radiator = CreateRadiator("foo", height, wavelength, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(radiator, nullptr); + + // Test for optical depths + std::size_t num_vertical_layers = 3; + std::size_t num_wavelength_bins = 2; + // Allocate array as 1D + double* optical_depths_1D = new double[num_wavelength_bins * num_vertical_layers]; + // Allocate an array of pointers to each row + double** optical_depths = new double*[num_vertical_layers]; + // Fill in the pointers to the rows + for (int row = 0; row < num_vertical_layers; row++) + { + optical_depths[row] = &optical_depths_1D[row * num_wavelength_bins]; + } + int i = 1; + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + optical_depths[row][col] = 10 * i; + i++; + } + } + SetRadiatorOpticalDepths(radiator, optical_depths[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + optical_depths[row][col] = -999.0; + } + } + GetRadiatorOpticalDepths(radiator, optical_depths[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(optical_depths[0][0], 10.0); + ASSERT_EQ(optical_depths[0][1], 20.0); + ASSERT_EQ(optical_depths[1][0], 30.0); + ASSERT_EQ(optical_depths[1][1], 40.0); + ASSERT_EQ(optical_depths[2][0], 50.0); + ASSERT_EQ(optical_depths[2][1], 60.0); + + // Test for single scattering albedos + double* albedos_1D = new double[num_wavelength_bins * num_vertical_layers]; + double** albedos = new double*[num_vertical_layers]; + for (int row = 0; row < num_vertical_layers; row++) + { + albedos[row] = &albedos_1D[row * num_wavelength_bins]; + } + i = 1; + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + albedos[row][col] = 100 * i; + i++; + } + } + SetRadiatorSingleScatteringAlbedos(radiator, albedos[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + albedos[row][col] = -999.0; + } + } + GetRadiatorSingleScatteringAlbedos(radiator, albedos[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(albedos[0][0], 100.0); + ASSERT_EQ(albedos[0][1], 200.0); + ASSERT_EQ(albedos[1][0], 300.0); + ASSERT_EQ(albedos[1][1], 400.0); + ASSERT_EQ(albedos[2][0], 500.0); + ASSERT_EQ(albedos[2][1], 600.0); + + // Test for asymmetery factors + double* factors_1D = new double[num_wavelength_bins * num_vertical_layers]; + double** factors = new double*[num_vertical_layers]; + for (int row = 0; row < num_vertical_layers; row++) + { + factors[row] = &factors_1D[row * num_wavelength_bins]; + } + i = 1; + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + factors[row][col] = 1 * i; + i++; + } + } + std::size_t num_streams = 1; + SetRadiatorAsymmetryFactors(radiator, factors[0], num_vertical_layers, num_wavelength_bins, num_streams, &error); + ASSERT_TRUE(IsSuccess(error)); + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + factors[row][col] = -999.0; + } + } + GetRadiatorAsymmetryFactors(radiator, factors[0], num_vertical_layers, num_wavelength_bins, num_streams, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(factors[0][0], 1); + ASSERT_EQ(factors[0][1], 2); + ASSERT_EQ(factors[1][0], 3); + ASSERT_EQ(factors[1][1], 4); + ASSERT_EQ(factors[2][0], 5); + ASSERT_EQ(factors[2][1], 6); + + // Clean up + DeleteRadiator(radiator, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteGrid(height, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteGrid(wavelength, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteError(&error); + delete[] optical_depths; + delete[] optical_depths_1D; + delete[] albedos; + delete[] albedos_1D; + delete[] factors; + delete[] factors_1D; +} + +TEST_F(TuvxCApiTest, CanCreateRadiatorMap) +{ + Error error; + RadiatorMap* radiator_map = CreateRadiatorMap(&error); + ASSERT_TRUE(IsSuccess(error)); + Grid* height = CreateGrid("height", "km", 3, &error); + Grid* wavelength = CreateGrid("wavelength", "nm", 2, &error); + Radiator* foo_radiator = CreateRadiator("foo", height, wavelength, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(foo_radiator, nullptr); + AddRadiator(radiator_map, foo_radiator, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(radiator_map, nullptr); + Grid* bar_height = CreateGrid("bar_height", "km", 3, &error); + Grid* bar_wavelength = CreateGrid("bar_wavelength", "nm", 2, &error); + Radiator* bar_radiator = CreateRadiator("bar", bar_height, bar_wavelength, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(bar_radiator, nullptr); + AddRadiator(radiator_map, bar_radiator, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(radiator_map, nullptr); + + // Test for optical depths + std::size_t num_vertical_layers = 3; + std::size_t num_wavelength_bins = 2; + double* optical_depths_1D = new double[num_wavelength_bins * num_vertical_layers]; + double** optical_depths = new double*[num_vertical_layers]; + for (int row = 0; row < num_vertical_layers; row++) + { + optical_depths[row] = &optical_depths_1D[row * num_wavelength_bins]; + } + int i = 1; + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + optical_depths[row][col] = 10 * i; + i++; + } + } + SetRadiatorOpticalDepths(foo_radiator, optical_depths[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + + // Test for single scattering albedos + double* albedos_1D = new double[num_wavelength_bins * num_vertical_layers]; + double** albedos = new double*[num_vertical_layers]; + for (int row = 0; row < num_vertical_layers; row++) + { + albedos[row] = &albedos_1D[row * num_wavelength_bins]; + } + i = 1; + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + albedos[row][col] = 100 * i; + i++; + } + } + SetRadiatorSingleScatteringAlbedos(foo_radiator, albedos[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + + // Test for asymmetery factors + std::size_t num_streams = 1; + double* factors_1D = new double[num_wavelength_bins * num_vertical_layers]; + double** factors = new double*[num_vertical_layers]; + for (int row = 0; row < num_vertical_layers; row++) + { + factors[row] = &factors_1D[row * num_wavelength_bins]; + } + i = 1; + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + factors[row][col] = 1 * i; + i++; + } + } + SetRadiatorAsymmetryFactors(foo_radiator, factors[0], num_vertical_layers, num_wavelength_bins, num_streams, &error); + ASSERT_TRUE(IsSuccess(error)); + + // Test for optical depths + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + optical_depths[row][col] = -999.0; + } + } + GetRadiatorOpticalDepths(foo_radiator, optical_depths[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(optical_depths[0][0], 10.0); + ASSERT_EQ(optical_depths[0][1], 20.0); + ASSERT_EQ(optical_depths[1][0], 30.0); + ASSERT_EQ(optical_depths[1][1], 40.0); + ASSERT_EQ(optical_depths[2][0], 50.0); + ASSERT_EQ(optical_depths[2][1], 60.0); + + // Test for single scattering albedos + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + albedos[row][col] = -999.0; + } + } + GetRadiatorSingleScatteringAlbedos(foo_radiator, albedos[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(albedos[0][0], 100.0); + ASSERT_EQ(albedos[0][1], 200.0); + ASSERT_EQ(albedos[1][0], 300.0); + ASSERT_EQ(albedos[1][1], 400.0); + ASSERT_EQ(albedos[2][0], 500.0); + ASSERT_EQ(albedos[2][1], 600.0); + + // Test for asymmetry factors + for (int row = 0; row < num_vertical_layers; row++) + { + for (int col = 0; col < num_wavelength_bins; col++) + { + factors[row][col] = -999.0; + } + } + GetRadiatorAsymmetryFactors(foo_radiator, factors[0], num_vertical_layers, num_wavelength_bins, num_streams, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(factors[0][0], 1); + ASSERT_EQ(factors[0][1], 2); + ASSERT_EQ(factors[1][0], 3); + ASSERT_EQ(factors[1][1], 4); + ASSERT_EQ(factors[2][0], 5); + ASSERT_EQ(factors[2][1], 6); + + // Test copy for radiator map + Radiator* foo_copy = GetRadiator(radiator_map, "foo", &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_NE(foo_copy, nullptr); + GetRadiatorOpticalDepths(foo_copy, optical_depths[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(optical_depths[0][0], 10.0); + ASSERT_EQ(optical_depths[0][1], 20.0); + ASSERT_EQ(optical_depths[1][0], 30.0); + ASSERT_EQ(optical_depths[1][1], 40.0); + ASSERT_EQ(optical_depths[2][0], 50.0); + ASSERT_EQ(optical_depths[2][1], 60.0); + GetRadiatorSingleScatteringAlbedos(foo_copy, albedos[0], num_vertical_layers, num_wavelength_bins, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(albedos[0][0], 100.0); + ASSERT_EQ(albedos[0][1], 200.0); + ASSERT_EQ(albedos[1][0], 300.0); + ASSERT_EQ(albedos[1][1], 400.0); + ASSERT_EQ(albedos[2][0], 500.0); + ASSERT_EQ(albedos[2][1], 600.0); + GetRadiatorAsymmetryFactors(foo_copy, factors[0], num_vertical_layers, num_wavelength_bins, 1, &error); + ASSERT_TRUE(IsSuccess(error)); + ASSERT_EQ(factors[0][0], 1); + ASSERT_EQ(factors[0][1], 2); + ASSERT_EQ(factors[1][0], 3); + ASSERT_EQ(factors[1][1], 4); + ASSERT_EQ(factors[2][0], 5); + ASSERT_EQ(factors[2][1], 6); + + // Clean up + DeleteRadiator(foo_radiator, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteRadiator(bar_radiator, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteRadiator(foo_copy, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteRadiatorMap(radiator_map, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteGrid(height, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteGrid(wavelength, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteGrid(bar_height, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteGrid(bar_wavelength, &error); + ASSERT_TRUE(IsSuccess(error)); + DeleteError(&error); + delete[] optical_depths; + delete[] optical_depths_1D; + delete[] albedos; + delete[] albedos_1D; + delete[] factors; + delete[] factors_1D; +} \ No newline at end of file diff --git a/src/tuvx/CMakeLists.txt b/src/tuvx/CMakeLists.txt index 4d9acdb4..46fa14bc 100644 --- a/src/tuvx/CMakeLists.txt +++ b/src/tuvx/CMakeLists.txt @@ -5,10 +5,14 @@ target_sources(musica interface_grid_map.F90 interface_profile.F90 interface_profile_map.F90 + interface_radiator.F90 + interface_radiator_map.F90 grid.cpp grid_map.cpp profile.cpp profile_map.cpp + radiator.cpp + radiator_map.cpp tuvx.cpp tuvx_util.F90 ) \ No newline at end of file diff --git a/src/tuvx/interface.F90 b/src/tuvx/interface.F90 index 8c756bf5..48ca85db 100644 --- a/src/tuvx/interface.F90 +++ b/src/tuvx/interface.F90 @@ -3,107 +3,112 @@ ! module tuvx_interface - use iso_c_binding, only : c_ptr, c_loc, c_int, c_size_t, c_char - use tuvx_core, only : core_t - use tuvx_grid_warehouse, only : grid_warehouse_t - use tuvx_profile_warehouse, only : profile_warehouse_t - use musica_tuvx_util, only : to_f_string, string_t_c - use musica_string, only : string_t - - implicit none +use iso_c_binding, only : c_ptr, c_loc, c_int, c_size_t, c_char +use tuvx_core, only : core_t +use tuvx_grid_warehouse, only : grid_warehouse_t +use tuvx_profile_warehouse, only : profile_warehouse_t +use tuvx_radiator_warehouse, only : radiator_warehouse_t +use musica_tuvx_util, only : to_f_string, string_t_c +use musica_string, only : string_t + +implicit none - private +private - contains +contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - function internal_create_tuvx(c_config_path, config_path_length, error_code) bind(C, name="InternalCreateTuvx") - use iso_c_binding, only: c_ptr, c_f_pointer + function internal_create_tuvx(c_config_path, config_path_length, error_code) & + bind(C, name="InternalCreateTuvx") + use iso_c_binding, only: c_ptr, c_f_pointer - ! arguments - character(kind=c_char), dimension(*), intent(in) :: c_config_path - integer(kind=c_size_t), value :: config_path_length - integer(kind=c_int), intent(out) :: error_code + ! arguments + character(kind=c_char), dimension(*), intent(in) :: c_config_path + integer(kind=c_size_t), value :: config_path_length + integer(kind=c_int), intent(out) :: error_code - ! local variables - character(len=:), allocatable :: f_config_path - type(c_ptr) :: internal_create_tuvx - type(core_t), pointer :: core - type(string_t) :: musica_config_path - integer :: i + ! local variables + character(len=:), allocatable :: f_config_path + type(c_ptr) :: internal_create_tuvx + type(core_t), pointer :: core + type(string_t) :: musica_config_path + integer :: i - allocate(character(len=config_path_length) :: f_config_path) - do i = 1, config_path_length - f_config_path(i:i) = c_config_path(i) - end do + allocate(character(len=config_path_length) :: f_config_path) + do i = 1, config_path_length + f_config_path(i:i) = c_config_path(i) + end do - musica_config_path = string_t(f_config_path) + musica_config_path = string_t(f_config_path) - core => core_t(musica_config_path) + core => core_t(musica_config_path) - deallocate(f_config_path) - error_code = 0 + deallocate(f_config_path) + error_code = 0 - internal_create_tuvx = c_loc(core) + internal_create_tuvx = c_loc(core) - end function internal_create_tuvx + end function internal_create_tuvx !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine internal_delete_tuvx(tuvx, error_code) bind(C, name="InternalDeleteTuvx") - use iso_c_binding, only: c_ptr, c_f_pointer + subroutine internal_delete_tuvx(tuvx, error_code) & + bind(C, name="InternalDeleteTuvx") + use iso_c_binding, only: c_ptr, c_f_pointer - ! arguments - type(c_ptr), value, intent(in) :: tuvx - integer(kind=c_int), intent(out) :: error_code + ! arguments + type(c_ptr), value, intent(in) :: tuvx + integer(kind=c_int), intent(out) :: error_code - ! local variables - type(core_t), pointer :: core - - call c_f_pointer(tuvx, core) - if (associated(core)) then - deallocate(core) - end if - end subroutine internal_delete_tuvx + ! local variables + type(core_t), pointer :: core + + call c_f_pointer(tuvx, core) + if (associated(core)) then + deallocate(core) + end if + end subroutine internal_delete_tuvx !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - function internal_get_grid_map(tuvx, error_code) result(grid_map_ptr) bind(C, name="InternalGetGridMap") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int - - ! arguments - type(c_ptr), intent(in), value :: tuvx - integer(kind=c_int), intent(out) :: error_code - - ! result - type(c_ptr) :: grid_map_ptr - - ! variables - type(core_t), pointer :: core - type(grid_warehouse_t), pointer :: grid_warehouse - - call c_f_pointer(tuvx, core) - grid_warehouse => core%get_grid_warehouse() - - grid_map_ptr = c_loc(grid_warehouse) + function internal_get_grid_map(tuvx, error_code) result(grid_map_ptr) & + bind(C, name="InternalGetGridMap") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int + + ! arguments + type(c_ptr), value, intent(in) :: tuvx + integer(kind=c_int), intent(out) :: error_code + + ! result + type(c_ptr) :: grid_map_ptr + + ! variables + type(core_t), pointer :: core + type(grid_warehouse_t), pointer :: grid_warehouse + + call c_f_pointer(tuvx, core) + grid_warehouse => core%get_grid_warehouse() + + grid_map_ptr = c_loc(grid_warehouse) - end function internal_get_grid_map + end function internal_get_grid_map !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - function internal_get_profile_map(tuvx, error_code) result(profile_map_ptr) bind(C, name="InternalGetProfileMap") + function internal_get_profile_map(tuvx, error_code) result(profile_map_ptr) & + bind(C, name="InternalGetProfileMap") use iso_c_binding, only: c_ptr, c_f_pointer, c_int ! arguments - type(c_ptr), intent(in), value :: tuvx + type(c_ptr), value, intent(in) :: tuvx integer(kind=c_int), intent(out) :: error_code ! result type(c_ptr) :: profile_map_ptr ! variables - type(core_t), pointer :: core + type(core_t), pointer :: core type(profile_warehouse_t), pointer :: profile_warehouse call c_f_pointer(tuvx, core) @@ -115,4 +120,28 @@ end function internal_get_profile_map !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -end module tuvx_interface + function internal_get_radiator_map(tuvx, error_code) result(radiator_map_ptr) & + bind(C, name="InternalGetRadiatorMap") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int + + ! arguments + type(c_ptr), value, intent(in) :: tuvx + integer(kind=c_int), intent(out) :: error_code + + ! result + type(c_ptr) :: radiator_map_ptr + + ! variables + type(core_t), pointer :: core + type(radiator_warehouse_t), pointer :: radiator_warehouse + + call c_f_pointer(tuvx, core) + radiator_warehouse => core%get_radiator_warehouse() + + radiator_map_ptr = c_loc(radiator_warehouse) + + end function internal_get_radiator_map + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +end module tuvx_interface \ No newline at end of file diff --git a/src/tuvx/interface_grid.F90 b/src/tuvx/interface_grid.F90 index 62bb09a1..6e211695 100644 --- a/src/tuvx/interface_grid.F90 +++ b/src/tuvx/interface_grid.F90 @@ -16,12 +16,11 @@ module tuvx_interface_grid function internal_create_grid(grid_name, grid_name_length, units, & units_length, num_sections, error_code) & bind(C, name="InternalCreateGrid") result(grid) - use iso_c_binding, only: c_ptr, c_f_pointer, c_char, c_loc, c_size_t, & - c_int + use iso_c_binding, only: c_ptr, c_f_pointer, c_char, c_loc, c_size_t, c_int use musica_tuvx_util, only: to_f_string use musica_string, only: string_t use tuvx_grid_from_host, only: grid_from_host_t - + ! arguments type(c_ptr) :: grid character(kind=c_char, len=1), dimension(*), intent(in) :: grid_name @@ -30,12 +29,12 @@ function internal_create_grid(grid_name, grid_name_length, units, & integer(kind=c_size_t), intent(in), value :: units_length integer(kind=c_size_t), intent(in), value :: num_sections integer(kind=c_int), intent(out) :: error_code - + ! variables type(grid_from_host_t), pointer :: f_grid type(string_t) :: f_name, f_units integer :: i - + allocate(character(len=grid_name_length) :: f_name%val_) do i = 1, grid_name_length f_name%val_(i:i) = grid_name(i) @@ -48,190 +47,190 @@ function internal_create_grid(grid_name, grid_name_length, units, & f_grid => grid_from_host_t(f_name, f_units, int(num_sections)) grid = c_loc(f_grid) - + end function internal_create_grid +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_delete_grid(grid, error_code) & + bind(C, name="InternalDeleteGrid") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int + + ! arguments + type(c_ptr), value, intent(in) :: grid + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(grid_t), pointer :: f_grid + + call c_f_pointer(grid, f_grid) + if (associated(f_grid)) then + deallocate(f_grid) + end if + + end subroutine internal_delete_grid + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! function internal_get_grid_updater(grid, error_code) & - bind(C, name="InternalGetGridUpdater") result(updater) - use iso_c_binding, only: c_ptr, c_f_pointer, c_loc, c_int - use tuvx_grid_from_host, only: grid_from_host_t, grid_updater_t - - ! arguments - type(c_ptr), value, intent(in) :: grid - integer(kind=c_int), intent(out) :: error_code - - ! output - type(c_ptr) :: updater - - ! variables - type(grid_from_host_t), pointer :: f_grid - type(grid_updater_t), pointer :: f_updater - - call c_f_pointer(grid, f_grid) - allocate(f_updater, source = grid_updater_t(f_grid)) - updater = c_loc(f_updater) - + bind(C, name="InternalGetGridUpdater") result(updater) + use iso_c_binding, only: c_ptr, c_f_pointer, c_loc, c_int + use tuvx_grid_from_host, only: grid_from_host_t, grid_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: grid + integer(kind=c_int), intent(out) :: error_code + + ! output + type(c_ptr) :: updater + + ! variables + type(grid_from_host_t), pointer :: f_grid + type(grid_updater_t), pointer :: f_updater + + call c_f_pointer(grid, f_grid) + allocate(f_updater, source = grid_updater_t(f_grid)) + updater = c_loc(f_updater) + end function internal_get_grid_updater !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine internal_delete_grid(grid, error_code) & - bind(C, name="InternalDeleteGrid") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int - - ! arguments - type(c_ptr), value, intent(in) :: grid - integer(kind=c_int), intent(out) :: error_code - - ! variables - type(grid_t), pointer :: f_grid - - call c_f_pointer(grid, f_grid) - if (associated(f_grid)) then - deallocate(f_grid) - end if - - end subroutine internal_delete_grid + subroutine internal_delete_grid_updater(updater, error_code) & + bind(C, name="InternalDeleteGridUpdater") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int + use tuvx_grid_from_host, only: grid_updater_t -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! arguments + type(c_ptr), value, intent(in) :: updater + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(grid_updater_t), pointer :: f_updater + + call c_f_pointer(updater, f_updater) + if (associated(f_updater)) then + deallocate(f_updater) + end if - subroutine internal_delete_grid_updater(updater, error_code) & - bind(C, name="InternalDeleteGridUpdater") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int - use tuvx_grid_from_host, only: grid_updater_t - - ! arguments - type(c_ptr), value, intent(in) :: updater - integer(kind=c_int), intent(out) :: error_code - - ! variables - type(grid_updater_t), pointer :: f_updater - - call c_f_pointer(updater, f_updater) - if (associated(f_updater)) then - deallocate(f_updater) - end if - - end subroutine internal_delete_grid_updater + end subroutine internal_delete_grid_updater !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine internal_set_edges(grid_updater, edges, num_edges, error_code) & - bind(C, name="InternalSetEdges") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t - use musica_constants, only: dk => musica_dk - use tuvx_grid_from_host, only: grid_updater_t - - ! arguments - type(c_ptr), value, intent(in) :: grid_updater - type(c_ptr), value, intent(in) :: edges - integer(kind=c_size_t), intent(in), value :: num_edges - integer(kind=c_int), intent(out) :: error_code - - ! variables - type(grid_updater_t), pointer :: f_updater - real(kind=dk), pointer :: f_edges(:) - - call c_f_pointer(grid_updater, f_updater) - call c_f_pointer(edges, f_edges, [num_edges]) - - if (size(f_updater%grid_%edge_) /= num_edges) then - error_code = 1 - return - end if - f_updater%grid_%edge_(:) = f_edges(:) - - end subroutine internal_set_edges + subroutine internal_set_edges(grid_updater, edges, num_edges, error_code) & + bind(C, name="InternalSetEdges") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_grid_from_host, only: grid_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: grid_updater + type(c_ptr), value, intent(in) :: edges + integer(kind=c_size_t), intent(in), value :: num_edges + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(grid_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_edges(:) + + call c_f_pointer(grid_updater, f_updater) + call c_f_pointer(edges, f_edges, [num_edges]) + + if (size(f_updater%grid_%edge_) /= num_edges) then + error_code = 1 + return + end if + f_updater%grid_%edge_(:) = f_edges(:) + + end subroutine internal_set_edges !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine internal_get_edges(grid_updater, edges, num_edges, error_code) & - bind(C, name="InternalGetEdges") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t - use musica_constants, only: dk => musica_dk - use tuvx_grid_from_host, only: grid_updater_t - - ! arguments - type(c_ptr), value, intent(in) :: grid_updater - type(c_ptr), value, intent(in) :: edges - integer(kind=c_size_t), intent(in), value :: num_edges - integer(kind=c_int), intent(out) :: error_code - - ! variables - type(grid_updater_t), pointer :: f_updater - real(kind=dk), pointer :: f_edges(:) - - call c_f_pointer(grid_updater, f_updater) - call c_f_pointer(edges, f_edges, [num_edges]) - - if (size(f_updater%grid_%edge_) /= num_edges) then - error_code = 1 - return - end if - f_edges(:) = f_updater%grid_%edge_(:) - - end subroutine internal_get_edges + subroutine internal_get_edges(grid_updater, edges, num_edges, error_code) & + bind(C, name="InternalGetEdges") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_grid_from_host, only: grid_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: grid_updater + type(c_ptr), value, intent(in) :: edges + integer(kind=c_size_t), intent(in), value :: num_edges + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(grid_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_edges(:) + + call c_f_pointer(grid_updater, f_updater) + call c_f_pointer(edges, f_edges, [num_edges]) + + if (size(f_updater%grid_%edge_) /= num_edges) then + error_code = 1 + return + end if + f_edges(:) = f_updater%grid_%edge_(:) + + end subroutine internal_get_edges !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine internal_set_midpoints(grid_updater, midpoints, num_midpoints, & - error_code) bind(C, name="InternalSetMidpoints") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_double, c_size_t - use musica_constants, only: dk => musica_dk - use tuvx_grid_from_host, only: grid_updater_t - - ! arguments - type(c_ptr), value, intent(in) :: grid_updater - type(c_ptr), value, intent(in) :: midpoints - integer(kind=c_int), intent(in), value :: num_midpoints - integer(kind=c_int), intent(out) :: error_code - - ! variables - type(grid_updater_t), pointer :: f_updater - real(kind=dk), pointer :: f_midpoints(:) - - call c_f_pointer(grid_updater, f_updater) - call c_f_pointer(midpoints, f_midpoints, [num_midpoints]) - - if (size(f_updater%grid_%mid_) /= num_midpoints) then - error_code = 1 - return - end if - f_updater%grid_%mid_(:) = f_midpoints(:) - - end subroutine internal_set_midpoints + subroutine internal_set_midpoints(grid_updater, midpoints, num_midpoints, & + error_code) bind(C, name="InternalSetMidpoints") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_double, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_grid_from_host, only: grid_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: grid_updater + type(c_ptr), value, intent(in) :: midpoints + integer(kind=c_int), intent(in), value :: num_midpoints + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(grid_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_midpoints(:) + + call c_f_pointer(grid_updater, f_updater) + call c_f_pointer(midpoints, f_midpoints, [num_midpoints]) + + if (size(f_updater%grid_%mid_) /= num_midpoints) then + error_code = 1 + return + end if + f_updater%grid_%mid_(:) = f_midpoints(:) + + end subroutine internal_set_midpoints !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine internal_get_midpoints(grid_updater, midpoints, num_midpoints, & - error_code) bind(C, name="InternalGetMidpoints") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_double, c_size_t - use musica_constants, only: dk => musica_dk - use tuvx_grid_from_host, only: grid_updater_t - - ! arguments - type(c_ptr), value, intent(in) :: grid_updater - type(c_ptr), value, intent(in) :: midpoints - integer(kind=c_int), intent(in), value :: num_midpoints - integer(kind=c_int), intent(out) :: error_code - - ! variables - type(grid_updater_t), pointer :: f_updater - real(kind=dk), pointer :: f_midpoints(:) - - call c_f_pointer(grid_updater, f_updater) - call c_f_pointer(midpoints, f_midpoints, [num_midpoints]) - - if (size(f_updater%grid_%mid_) /= num_midpoints) then - error_code = 1 - return - end if - f_midpoints(:) = f_updater%grid_%mid_(:) - - end subroutine internal_get_midpoints + subroutine internal_get_midpoints(grid_updater, midpoints, num_midpoints, & + error_code) bind(C, name="InternalGetMidpoints") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_double, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_grid_from_host, only: grid_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: grid_updater + type(c_ptr), value, intent(in) :: midpoints + integer(kind=c_int), intent(in), value :: num_midpoints + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(grid_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_midpoints(:) + + call c_f_pointer(grid_updater, f_updater) + call c_f_pointer(midpoints, f_midpoints, [num_midpoints]) + + if (size(f_updater%grid_%mid_) /= num_midpoints) then + error_code = 1 + return + end if + f_midpoints(:) = f_updater%grid_%mid_(:) + + end subroutine internal_get_midpoints !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -end module tuvx_interface_grid +end module tuvx_interface_grid \ No newline at end of file diff --git a/src/tuvx/interface_grid_map.F90 b/src/tuvx/interface_grid_map.F90 index 977eca1a..c0eaf250 100644 --- a/src/tuvx/interface_grid_map.F90 +++ b/src/tuvx/interface_grid_map.F90 @@ -17,31 +17,31 @@ module tuvx_interface_grid_map !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - function internal_create_grid_map(error_code) result(grid_map) & - bind(C, name="InternalCreateGridMap") - use iso_c_binding, only: c_ptr, c_int, c_null_ptr - use tuvx_grid_warehouse, only: grid_warehouse_t - - ! arguments - integer(kind=c_int), intent(out) :: error_code - - ! result - type(c_ptr) :: grid_map - - ! variables - class(grid_warehouse_t), pointer :: f_grid_warehouse - - f_grid_warehouse => grid_warehouse_t() - select type(f_grid_warehouse) - type is(grid_warehouse_t) - grid_map = c_loc(f_grid_warehouse) - error_code = 0 - class default - error_code = 1 - grid_map = c_null_ptr - end select - - end function internal_create_grid_map + function internal_create_grid_map(error_code) result(grid_map) & + bind(C, name="InternalCreateGridMap") + use iso_c_binding, only: c_ptr, c_int, c_null_ptr + use tuvx_grid_warehouse, only: grid_warehouse_t + + ! arguments + integer(kind=c_int), intent(out) :: error_code + + ! result + type(c_ptr) :: grid_map + + ! variables + class(grid_warehouse_t), pointer :: f_grid_warehouse + + f_grid_warehouse => grid_warehouse_t() + select type(f_grid_warehouse) + type is(grid_warehouse_t) + grid_map = c_loc(f_grid_warehouse) + error_code = 0 + class default + error_code = 1 + grid_map = c_null_ptr + end select + + end function internal_create_grid_map !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -90,56 +90,56 @@ end subroutine internal_add_grid !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - function interal_get_grid(grid_map, c_grid_name, c_grid_name_length, & - c_grid_units, c_grid_units_length, error_code) & - result(grid_ptr) bind(C, name="InternalGetGrid") - use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_char, c_size_t, & - c_null_ptr, c_loc - use tuvx_grid_from_host, only: grid_from_host_t - - ! arguments - type(c_ptr), intent(in), value :: grid_map - character(len=1, kind=c_char), dimension(*), intent(in) :: c_grid_name - integer(kind=c_size_t), value :: c_grid_name_length - character(len=1, kind=c_char), dimension(*), intent(in) :: c_grid_units - integer(kind=c_size_t), value :: c_grid_units_length - integer(kind=c_int), intent(out) :: error_code - - ! variables - class(grid_t), pointer :: f_grid - type(grid_warehouse_t), pointer :: grid_warehouse - character(len=:), allocatable :: f_grid_name - character(len=:), allocatable :: f_grid_units - integer :: i - - ! result - type(c_ptr) :: grid_ptr - - allocate(character(len=c_grid_name_length) :: f_grid_name) - do i = 1, c_grid_name_length - f_grid_name(i:i) = c_grid_name(i) - end do - - allocate(character(len=c_grid_units_length) :: f_grid_units) - do i = 1, c_grid_units_length - f_grid_units(i:i) = c_grid_units(i) - end do - - call c_f_pointer(grid_map, grid_warehouse) - - f_grid => grid_warehouse%get_grid(f_grid_name, f_grid_units) - - select type(f_grid) - type is(grid_from_host_t) - error_code = 0 - grid_ptr = c_loc(f_grid) - class default - error_code = 1 - deallocate(f_grid) - grid_ptr = c_null_ptr - end select - - end function interal_get_grid + function interal_get_grid(grid_map, c_grid_name, c_grid_name_length, & + c_grid_units, c_grid_units_length, error_code) & + result(grid_ptr) bind(C, name="InternalGetGrid") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_char, c_size_t, & + c_null_ptr, c_loc + use tuvx_grid_from_host, only: grid_from_host_t + + ! arguments + type(c_ptr), intent(in), value :: grid_map + character(len=1, kind=c_char), dimension(*), intent(in) :: c_grid_name + integer(kind=c_size_t), value :: c_grid_name_length + character(len=1, kind=c_char), dimension(*), intent(in) :: c_grid_units + integer(kind=c_size_t), value :: c_grid_units_length + integer(kind=c_int), intent(out) :: error_code + + ! variables + class(grid_t), pointer :: f_grid + type(grid_warehouse_t), pointer :: grid_warehouse + character(len=:), allocatable :: f_grid_name + character(len=:), allocatable :: f_grid_units + integer :: i + + ! result + type(c_ptr) :: grid_ptr + + allocate(character(len=c_grid_name_length) :: f_grid_name) + do i = 1, c_grid_name_length + f_grid_name(i:i) = c_grid_name(i) + end do + + allocate(character(len=c_grid_units_length) :: f_grid_units) + do i = 1, c_grid_units_length + f_grid_units(i:i) = c_grid_units(i) + end do + + call c_f_pointer(grid_map, grid_warehouse) + + f_grid => grid_warehouse%get_grid(f_grid_name, f_grid_units) + + select type(f_grid) + type is(grid_from_host_t) + error_code = 0 + grid_ptr = c_loc(f_grid) + class default + error_code = 1 + deallocate(f_grid) + grid_ptr = c_null_ptr + end select + + end function interal_get_grid !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! diff --git a/src/tuvx/interface_profile_map.F90 b/src/tuvx/interface_profile_map.F90 index a6057095..47e9cba7 100644 --- a/src/tuvx/interface_profile_map.F90 +++ b/src/tuvx/interface_profile_map.F90 @@ -3,11 +3,11 @@ ! module tuvx_interface_profile_map - use iso_c_binding, only : c_ptr, c_loc, c_int, c_size_t, c_char - use tuvx_profile, only : profile_t + use iso_c_binding, only : c_ptr, c_loc, c_int, c_size_t, c_char + use tuvx_profile, only : profile_t use tuvx_profile_warehouse, only : profile_warehouse_t - use musica_tuvx_util, only : to_f_string, string_t_c - use musica_string, only : string_t + use musica_tuvx_util, only : to_f_string, string_t_c + use musica_string, only : string_t implicit none diff --git a/src/tuvx/interface_radiator.F90 b/src/tuvx/interface_radiator.F90 new file mode 100644 index 00000000..6bfd06a6 --- /dev/null +++ b/src/tuvx/interface_radiator.F90 @@ -0,0 +1,324 @@ +! Copyright (C) 2023-2024 National Center for Atmospheric Research +! SPDX-License-Identifier: Apache-2.0 +! +module tuvx_interface_radiator + use tuvx_radiator, only : radiator_t + + implicit none + + private + + contains + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + function internal_create_radiator(radiator_name, radiator_name_length, & + height_grid_updater_c, wavelength_grid_updater_c, error_code) & + result(radiator) bind(C, name="InternalCreateRadiator") + use iso_c_binding, only: c_ptr, c_f_pointer, c_char, c_loc, c_size_t, c_int + use musica_string, only: string_t + use tuvx_radiator_from_host, only: radiator_from_host_t + use tuvx_grid_from_host, only: grid_updater_t + + ! arguments + type(c_ptr) :: radiator + character(kind=c_char, len=1), dimension(*), intent(in) :: radiator_name + integer(kind=c_size_t), value, intent(in) :: radiator_name_length + type(c_ptr), value, intent(in) :: height_grid_updater_c + type(c_ptr), value, intent(in) :: wavelength_grid_updater_c + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_from_host_t), pointer :: f_radiator + type(string_t) :: f_name + type(grid_updater_t), pointer :: f_height_grid_updater + type(grid_updater_t), pointer :: f_wavelength_grid_updater + integer :: i + + allocate(character(len=radiator_name_length) :: f_name%val_) + do i = 1, radiator_name_length + f_name%val_(i:i) = radiator_name(i) + end do + + call c_f_pointer(height_grid_updater_c, f_height_grid_updater) + call c_f_pointer(wavelength_grid_updater_c, f_wavelength_grid_updater) + f_radiator => radiator_from_host_t(f_name, f_height_grid_updater%grid_, & + f_wavelength_grid_updater%grid_) + radiator = c_loc(f_radiator) + + end function internal_create_radiator + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_delete_radiator(radiator, error_code) & + bind(C, name="InternalDeleteRadiator") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int + + ! arguments + type(c_ptr), value, intent(in) :: radiator + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_t), pointer :: f_radiator + + call c_f_pointer(radiator, f_radiator) + if (associated(f_radiator)) then + deallocate(f_radiator) + end if + + end subroutine internal_delete_radiator + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + function internal_get_radiator_updater(radiator, error_code) & + bind(C, name="InternalGetRadiatorUpdater") result(updater) + use iso_c_binding, only: c_ptr, c_f_pointer, c_loc, c_int + use tuvx_radiator_from_host, only: radiator_from_host_t, radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator + integer(kind=c_int), intent(out) :: error_code + + ! output + type(c_ptr) :: updater + + ! variables + type(radiator_from_host_t), pointer :: f_radiator + type(radiator_updater_t), pointer :: f_updater + + call c_f_pointer(radiator, f_radiator) + allocate(f_updater, source = radiator_updater_t(f_radiator)) + updater = c_loc(f_updater) + + end function internal_get_radiator_updater + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_delete_radiator_updater(updater, error_code) & + bind(C, name="InternalDeleteRadiatorUpdater") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: updater + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + + call c_f_pointer(updater, f_updater) + if (associated(f_updater)) then + deallocate(f_updater) + end if + + end subroutine internal_delete_radiator_updater + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_set_optical_depths(radiator_updater, optical_depths, & + num_vertical_layers, num_wavelength_bins, error_code) & + bind(C, name="InternalSetOpticalDepths") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_updater + type(c_ptr), value, intent(in) :: optical_depths + integer(kind=c_size_t), value, intent(in) :: num_vertical_layers + integer(kind=c_size_t), value, intent(in) :: num_wavelength_bins + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_optical_depths(:,:) + + call c_f_pointer(radiator_updater, f_updater) + call c_f_pointer(optical_depths, f_optical_depths, & + [num_vertical_layers, num_wavelength_bins]) + + if ((size(f_updater%radiator_%state_%layer_OD_, 1) /= num_vertical_layers) & + .or. (size(f_updater%radiator_%state_%layer_OD_, 2) /= num_wavelength_bins)) & + then + error_code = 1 + return + end if + f_updater%radiator_%state_%layer_OD_(:,:) = f_optical_depths(:,:) + + end subroutine internal_set_optical_depths + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_get_optical_depths(radiator_updater, optical_depths, & + num_vertical_layers, num_wavelength_bins, error_code) & + bind(C, name="InternalGetOpticalDepths") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_updater + type(c_ptr), value, intent(in) :: optical_depths + integer(kind=c_size_t), value, intent(in) :: num_vertical_layers + integer(kind=c_size_t), value, intent(in) :: num_wavelength_bins + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_optical_depths(:,:) + + call c_f_pointer(radiator_updater, f_updater) + call c_f_pointer(optical_depths, f_optical_depths, & + [num_vertical_layers, num_wavelength_bins]) + + if ((size(f_updater%radiator_%state_%layer_OD_, 1) /= num_vertical_layers) & + .or. (size(f_updater%radiator_%state_%layer_OD_, 2) /= num_wavelength_bins)) & + then + error_code = 1 + return + end if + f_optical_depths(:,:) = f_updater%radiator_%state_%layer_OD_(:,:) + + end subroutine internal_get_optical_depths + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_set_single_scattering_albedos(radiator_updater, & + single_scattering_albedos, num_vertical_layers, num_wavelength_bins, & + error_code) bind(C, name="InternalSetSingleScatteringAlbedos") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_updater + type(c_ptr), value, intent(in) :: single_scattering_albedos + integer(kind=c_size_t), value, intent(in) :: num_vertical_layers + integer(kind=c_size_t), value, intent(in) :: num_wavelength_bins + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_single_scattering_albedos(:,:) + + call c_f_pointer(radiator_updater, f_updater) + call c_f_pointer(single_scattering_albedos, f_single_scattering_albedos, & + [num_vertical_layers, num_wavelength_bins]) + + if ((size(f_updater%radiator_%state_%layer_SSA_, 1) /= num_vertical_layers) & + .or. (size(f_updater%radiator_%state_%layer_SSA_, 2) /= num_wavelength_bins)) & + then + error_code = 1 + return + end if + f_updater%radiator_%state_%layer_SSA_(:,:) = f_single_scattering_albedos(:,:) + + end subroutine internal_set_single_scattering_albedos + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_get_single_scattering_albedos(radiator_updater, & + single_scattering_albedos, num_vertical_layers, num_wavelength_bins, & + error_code) bind(C, name="InternalGetSingleScatteringAlbedos") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_updater + type(c_ptr), value, intent(in) :: single_scattering_albedos + integer(kind=c_size_t), value, intent(in) :: num_vertical_layers + integer(kind=c_size_t), value, intent(in) :: num_wavelength_bins + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_single_scattering_albedos(:,:) + + call c_f_pointer(radiator_updater, f_updater) + call c_f_pointer(single_scattering_albedos, f_single_scattering_albedos, & + [num_vertical_layers, num_wavelength_bins]) + + if ((size(f_updater%radiator_%state_%layer_SSA_, 1) /= num_vertical_layers) & + .or. (size(f_updater%radiator_%state_%layer_SSA_, 2) /= num_wavelength_bins)) & + then + error_code = 1 + return + end if + f_single_scattering_albedos(:,:) = f_updater%radiator_%state_%layer_SSA_(:,:) + + end subroutine internal_get_single_scattering_albedos + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_set_asymmetry_factors(radiator_updater, & + asymmetry_factors, num_vertical_layers, num_wavelength_bins, num_streams, & + error_code) bind(C, name="InternalSetAsymmetryFactors") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_updater + type(c_ptr), value, intent(in) :: asymmetry_factors + integer(kind=c_size_t), value, intent(in) :: num_vertical_layers + integer(kind=c_size_t), value, intent(in) :: num_wavelength_bins + integer(kind=c_size_t), value, intent(in) :: num_streams + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_asymmetry_factors(:,:,:) + + call c_f_pointer(radiator_updater, f_updater) + call c_f_pointer(asymmetry_factors, f_asymmetry_factors, & + [num_vertical_layers, num_wavelength_bins, num_streams]) + + if ((size(f_updater%radiator_%state_%layer_G_, 1) /= num_vertical_layers) & + .or. (size(f_updater%radiator_%state_%layer_G_, 2) /= num_wavelength_bins) & + .or. (size(f_updater%radiator_%state_%layer_G_, 3) /= num_streams)) then + error_code = 1 + return + end if + f_updater%radiator_%state_%layer_G_(:,:,:) = f_asymmetry_factors(:,:,:) + + end subroutine internal_set_asymmetry_factors + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_get_asymmetry_factors(radiator_updater, & + asymmetry_factors, num_vertical_layers, num_wavelength_bins, num_streams, & + error_code) bind(C, name="InternalGetAsymmetryFactors") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_size_t + use musica_constants, only: dk => musica_dk + use tuvx_radiator_from_host, only: radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_updater + type(c_ptr), value, intent(in) :: asymmetry_factors + integer(kind=c_size_t), value, intent(in) :: num_vertical_layers + integer(kind=c_size_t), value, intent(in) :: num_wavelength_bins + integer(kind=c_size_t), value, intent(in) :: num_streams + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_updater_t), pointer :: f_updater + real(kind=dk), pointer :: f_asymmetry_factors(:,:,:) + + call c_f_pointer(radiator_updater, f_updater) + call c_f_pointer(asymmetry_factors, f_asymmetry_factors, & + [num_vertical_layers, num_wavelength_bins, num_streams]) + + if ((size(f_updater%radiator_%state_%layer_G_, 1) /= num_vertical_layers) & + .or. (size(f_updater%radiator_%state_%layer_G_, 2) /= num_wavelength_bins) & + .or. (size(f_updater%radiator_%state_%layer_G_, 3) /= num_streams)) then + error_code = 1 + return + end if + f_asymmetry_factors(:,:,:) = f_updater%radiator_%state_%layer_G_(:,:,:) + +end subroutine internal_get_asymmetry_factors + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +end module tuvx_interface_radiator \ No newline at end of file diff --git a/src/tuvx/interface_radiator_map.F90 b/src/tuvx/interface_radiator_map.F90 new file mode 100644 index 00000000..2ddea28e --- /dev/null +++ b/src/tuvx/interface_radiator_map.F90 @@ -0,0 +1,177 @@ +! Copyright (C) 2023-2024 National Center for Atmospheric Research +! SPDX-License-Identifier: Apache-2.0 +! +module tuvx_interface_radiator_map + + use iso_c_binding, only : c_ptr, c_loc, c_int, c_size_t, c_char + use tuvx_radiator, only : radiator_t + use tuvx_radiator_warehouse, only : radiator_warehouse_t + use musica_tuvx_util, only : to_f_string, string_t_c + use musica_string, only : string_t + + implicit none + + private + + contains + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + function internal_create_radiator_map(error_code) result(radiator_map) & + bind(C, name="InternalCreateRadiatorMap") + use iso_c_binding, only: c_ptr, c_int, c_null_ptr + use tuvx_radiator_warehouse, only: radiator_warehouse_t + + ! arguments + integer(kind=c_int), intent(out) :: error_code + + ! result + type(c_ptr) :: radiator_map + + ! variables + class(radiator_warehouse_t), pointer :: f_radiator_warehouse + + f_radiator_warehouse => radiator_warehouse_t() + select type(f_radiator_warehouse) + type is(radiator_warehouse_t) + radiator_map = c_loc(f_radiator_warehouse) + error_code = 0 + class default + error_code = 1 + radiator_map = c_null_ptr + end select + + end function internal_create_radiator_map + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_delete_radiator_map(radiator_map, error_code) & + bind(C, name="InternalDeleteRadiatorMap") + use iso_c_binding, only: c_ptr, c_int, c_f_pointer + use tuvx_radiator_warehouse, only: radiator_warehouse_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_map + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_warehouse_t), pointer :: f_radiator_warehouse + + call c_f_pointer(radiator_map, f_radiator_warehouse) + deallocate(f_radiator_warehouse) + error_code = 0 + +end subroutine internal_delete_radiator_map + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine internal_add_radiator(radiator_map, radiator, error_code) & + bind(C, name="InternalAddRadiator") + use iso_c_binding, only: c_ptr, c_int, c_f_pointer + use tuvx_radiator_warehouse, only: radiator_warehouse_t + use tuvx_radiator_from_host, only: radiator_from_host_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_map + type(c_ptr), value, intent(in) :: radiator + integer(kind=c_int), intent(out) :: error_code + + ! variables + type(radiator_warehouse_t), pointer :: f_radiator_warehouse + type(radiator_from_host_t), pointer :: f_radiator + + call c_f_pointer(radiator_map, f_radiator_warehouse) + call c_f_pointer(radiator, f_radiator) + + error_code = 0 + call f_radiator_warehouse%add(f_radiator) + + end subroutine internal_add_radiator + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + function internal_get_radiator(radiator_map, c_radiator_name, & + c_radiator_name_length, error_code) & + result(radiator_ptr) bind(C, name="InternalGetRadiator") + use iso_c_binding, only: c_ptr, c_f_pointer, c_int, c_char, c_size_t, & + c_null_ptr, c_loc + use tuvx_radiator_from_host, only: radiator_from_host_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_map + character(len=1, kind=c_char), dimension(*), intent(in) :: c_radiator_name + integer(kind=c_size_t), value, intent(in) :: c_radiator_name_length + integer(kind=c_int), intent(out) :: error_code + + ! variables + class(radiator_t), pointer :: f_radiator + class(radiator_t), pointer :: f_radiator_ptr + type(radiator_warehouse_t), pointer :: radiator_warehouse + character(len=:), allocatable :: f_radiator_name + integer :: i + + ! result + type(c_ptr) :: radiator_ptr + + allocate(character(len=c_radiator_name_length) :: f_radiator_name) + do i = 1, c_radiator_name_length + f_radiator_name(i:i) = c_radiator_name(i) + end do + + call c_f_pointer(radiator_map, radiator_warehouse) + + if (.not. radiator_warehouse%exists(f_radiator_name)) then + error_code = 1 + radiator_ptr = c_null_ptr + else + f_radiator_ptr => radiator_warehouse%get_radiator(f_radiator_name) + allocate(f_radiator, source = f_radiator_ptr) + nullify(f_radiator_ptr) + + select type(f_radiator) + type is(radiator_from_host_t) + error_code = 0 + radiator_ptr = c_loc(f_radiator) + class default + error_code = 1 + deallocate(f_radiator) + radiator_ptr = c_null_ptr + end select + end if + + end function internal_get_radiator + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + function internal_get_radiator_updater_from_map(radiator_map, radiator, error_code) & + result(updater) bind(C, name="InternalGetRadiatorUpdaterFromMap") + use iso_c_binding, only: c_ptr, c_f_pointer, c_loc + use tuvx_radiator_warehouse, only: radiator_warehouse_t + use tuvx_radiator_from_host, only: radiator_from_host_t, radiator_updater_t + + ! arguments + type(c_ptr), value, intent(in) :: radiator_map + type(c_ptr), value, intent(in) :: radiator + integer(kind=c_int), intent(out) :: error_code + + ! output + type(c_ptr) :: updater + + ! variables + type(radiator_warehouse_t), pointer :: f_radiator_warehouse + type(radiator_from_host_t), pointer :: f_radiator + type(radiator_updater_t), pointer :: f_updater + + call c_f_pointer(radiator_map, f_radiator_warehouse) + call c_f_pointer(radiator, f_radiator) + + error_code = 0 + allocate(f_updater) + f_updater = f_radiator_warehouse%get_updater(f_radiator) + updater = c_loc(f_updater) + + end function internal_get_radiator_updater_from_map + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + +end module tuvx_interface_radiator_map \ No newline at end of file diff --git a/src/tuvx/radiator.cpp b/src/tuvx/radiator.cpp new file mode 100644 index 00000000..57271029 --- /dev/null +++ b/src/tuvx/radiator.cpp @@ -0,0 +1,268 @@ +// Copyright (C) 2023-2024 National Center for Atmospheric Research +// SPDX-License-Identifier: Apache-2.0 +#include + +#include +#include +#include + +namespace musica +{ + + // Radiator external C API functions + + Radiator *CreateRadiator(const char *radiator_name, Grid *height_grid, Grid *wavelength_grid, Error *error) + { + DeleteError(error); + return new Radiator(radiator_name, height_grid, wavelength_grid, error); + } + + void DeleteRadiator(Radiator *radiator, Error *error) + { + DeleteError(error); + try + { + delete radiator; + } + catch (const std::system_error &e) + { + *error = ToError(e); + return; + } + *error = NoError(); + } + + void SetRadiatorOpticalDepths( + Radiator *radiator, + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + DeleteError(error); + radiator->SetOpticalDepths(optical_depths, num_vertical_layers, num_wavelength_bins, error); + } + + void GetRadiatorOpticalDepths( + Radiator *radiator, + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + DeleteError(error); + radiator->GetOpticalDepths(optical_depths, num_vertical_layers, num_wavelength_bins, error); + } + + void SetRadiatorSingleScatteringAlbedos( + Radiator *radiator, + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + DeleteError(error); + radiator->SetSingleScatteringAlbedos(single_scattering_albedos, num_vertical_layers, num_wavelength_bins, error); + } + + void GetRadiatorSingleScatteringAlbedos( + Radiator *radiator, + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + DeleteError(error); + radiator->GetSingleScatteringAlbedos(single_scattering_albedos, num_vertical_layers, num_wavelength_bins, error); + } + + void SetRadiatorAsymmetryFactors( + Radiator *radiator, + double *asymmetry_factors, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error) + { + DeleteError(error); + radiator->SetAsymmetryFactors(asymmetry_factors, num_vertical_layers, num_wavelength_bins, num_streams, error); + } + + void GetRadiatorAsymmetryFactors( + Radiator *radiator, + double *asymmetry_factors, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error) + { + DeleteError(error); + radiator->GetAsymmetryFactors(asymmetry_factors, num_vertical_layers, num_wavelength_bins, num_streams, error); + } + + // Radiation class functions + + Radiator::Radiator(const char *radiator_name, Grid *height_grid, Grid *wavelength_grid, Error *error) + { + int error_code = 0; + radiator_ = InternalCreateRadiator( + radiator_name, strlen(radiator_name), height_grid->updater_, wavelength_grid->updater_, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to create radiator") }; + return; + } + updater_ = InternalGetRadiatorUpdater(radiator_, &error_code); + if (error_code != 0) + { + InternalDeleteRadiator(radiator_, &error_code); + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to get radiator updater") }; + return; + } + *error = NoError(); + } + + Radiator::~Radiator() + { + int error_code = 0; + if (radiator_ != nullptr) + InternalDeleteRadiator(radiator_, &error_code); + if (updater_ != nullptr) + InternalDeleteRadiatorUpdater(updater_, &error_code); + radiator_ = nullptr; + updater_ = nullptr; + } + + void Radiator::SetOpticalDepths( + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + int error_code = 0; + if (updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator is not updatable") }; + return; + } + InternalSetOpticalDepths(updater_, optical_depths, num_vertical_layers, num_wavelength_bins, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to set optical depths") }; + return; + } + *error = NoError(); + } + + void Radiator::GetOpticalDepths( + double *optical_depths, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + int error_code = 0; + if (updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator is not updatable") }; + return; + } + InternalGetOpticalDepths(updater_, optical_depths, num_vertical_layers, num_wavelength_bins, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to get optical depths") }; + return; + } + *error = NoError(); + } + + void Radiator::SetSingleScatteringAlbedos( + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + int error_code = 0; + if (updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator is not updatable") }; + return; + } + InternalSetSingleScatteringAlbedos( + updater_, single_scattering_albedos, num_vertical_layers, num_wavelength_bins, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to set single scattering albedos") }; + return; + } + *error = NoError(); + } + + void Radiator::GetSingleScatteringAlbedos( + double *single_scattering_albedos, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + Error *error) + { + int error_code = 0; + if (updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator is not updatable") }; + return; + } + InternalGetSingleScatteringAlbedos( + updater_, single_scattering_albedos, num_vertical_layers, num_wavelength_bins, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to get single scattering albedos") }; + return; + } + *error = NoError(); + } + + void Radiator::SetAsymmetryFactors( + double *asymmetry_factors, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error) + { + int error_code = 0; + if (updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator is not updatable") }; + return; + } + InternalSetAsymmetryFactors( + updater_, asymmetry_factors, num_vertical_layers, num_wavelength_bins, num_streams, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to set asymmetry factors") }; + return; + } + *error = NoError(); + } + + void Radiator::GetAsymmetryFactors( + double *asymmetry_factors, + std::size_t num_vertical_layers, + std::size_t num_wavelength_bins, + std::size_t num_streams, + Error *error) + { + int error_code = 0; + if (updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator is not updatable") }; + return; + } + InternalGetAsymmetryFactors( + updater_, asymmetry_factors, num_vertical_layers, num_wavelength_bins, num_streams, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to get asymmetry factors") }; + return; + } + *error = NoError(); + } + +} // namespace musica \ No newline at end of file diff --git a/src/tuvx/radiator_map.cpp b/src/tuvx/radiator_map.cpp new file mode 100644 index 00000000..92263dda --- /dev/null +++ b/src/tuvx/radiator_map.cpp @@ -0,0 +1,180 @@ +// Copyright (C) 2023-2024 National Center for Atmospheric Research +// SPDX-License-Identifier: Apache-2.0 +#include + +#include +#include +#include + +namespace musica +{ + + // RadiatordMap external C API functions + + RadiatorMap *CreateRadiatorMap(Error *error) + { + DeleteError(error); + return new RadiatorMap(error); + } + + void DeleteRadiatorMap(RadiatorMap *radiator_map, Error *error) + { + DeleteError(error); + try + { + delete radiator_map; + } + catch (const std::system_error &e) + { + *error = ToError(e); + } + *error = NoError(); + } + + void AddRadiator(RadiatorMap *radiator_map, Radiator *radiator, Error *error) + { + DeleteError(error); + radiator_map->AddRadiator(radiator, error); + } + + Radiator *GetRadiator(RadiatorMap *radiator_map, const char *radiator_name, Error *error) + { + DeleteError(error); + return radiator_map->GetRadiator(radiator_name, error); + } + + // RadiatordMap class functions + + RadiatorMap::RadiatorMap(Error *error) + { + int error_code = 0; + radiator_map_ = InternalCreateRadiatorMap(&error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to create radiator map") }; + } + owns_radiator_map_ = true; + *error = NoError(); + } + + RadiatorMap::~RadiatorMap() + { + int error_code = 0; + if (radiator_map_ != nullptr && owns_radiator_map_) + { + InternalDeleteRadiatorMap(radiator_map_, &error_code); + } + radiator_map_ = nullptr; + owns_radiator_map_ = false; + } + + void RadiatorMap::AddRadiator(Radiator *radiator, Error *error) + { + if (radiator_map_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator map is null") }; + return; + } + if (radiator->radiator_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Cannot add unowned radiator to radiator map") }; + return; + } + if (radiator->updater_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Cannot add radiator in invalid state") }; + return; + } + + int error_code = 0; + + try + { + InternalAddRadiator(radiator_map_, radiator->radiator_, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to add radiator to radiator map") }; + } + InternalDeleteRadiatorUpdater(radiator->updater_, &error_code); + if (error_code != 0) + { + *error = Error{ 1, + CreateString(MUSICA_ERROR_CATEGORY), + CreateString("Failed to delete updater after transfer of ownership to radiator map") }; + } + radiator->updater_ = InternalGetRadiatorUpdaterFromMap(radiator_map_, radiator->radiator_, &error_code); + if (error_code != 0) + { + *error = Error{ 1, + CreateString(MUSICA_ERROR_CATEGORY), + CreateString("Failed to get updater after transfer of ownership to radiator map") }; + } + InternalDeleteRadiator(radiator->radiator_, &error_code); + if (error_code != 0) + { + *error = Error{ 1, + CreateString(MUSICA_ERROR_CATEGORY), + CreateString("Failed to delete radiator during transfer of ownership to radiator map") }; + } + radiator->radiator_ = nullptr; + } + catch (const std::system_error &e) + { + *error = ToError(e); + } + catch (...) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Internal error adding radiator") }; + } + *error = NoError(); + } + + Radiator *RadiatorMap::GetRadiator(const char *radiator_name, Error *error) + { + if (radiator_map_ == nullptr) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Radiator map is null") }; + return nullptr; + } + + Radiator *radiator = nullptr; + + try + { + int error_code = 0; + void *radiator_ptr = InternalGetRadiator(radiator_map_, radiator_name, strlen(radiator_name), &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to get radiator from radiator map") }; + return nullptr; + } + void *updater_ptr = InternalGetRadiatorUpdaterFromMap(radiator_map_, radiator_ptr, &error_code); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to get radiator updater") }; + InternalDeleteRadiator(radiator_ptr, &error_code); + return nullptr; + } + InternalDeleteRadiator(radiator_ptr, &error_code); + if (error_code != 0) + { + *error = + Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to delete radiator after getting updater") }; + InternalDeleteRadiatorUpdater(updater_ptr, &error_code); + return nullptr; + } + radiator = new Radiator(updater_ptr); + } + catch (const std::system_error &e) + { + *error = ToError(e); + } + catch (...) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Internal error getting radiator") }; + } + *error = NoError(); + return radiator; + } + +} // namespace musica \ No newline at end of file diff --git a/src/tuvx/tuvx.cpp b/src/tuvx/tuvx.cpp index 7c17cf02..9abd977d 100644 --- a/src/tuvx/tuvx.cpp +++ b/src/tuvx/tuvx.cpp @@ -58,6 +58,12 @@ namespace musica return tuvx->CreateProfileMap(error); } + RadiatorMap *GetRadiatorMap(TUVX *tuvx, Error *error) + { + DeleteError(error); + return tuvx->CreateRadiatorMap(error); + } + // TUVX class functions TUVX::TUVX() @@ -131,4 +137,17 @@ namespace musica return profile_map; } + RadiatorMap *TUVX::CreateRadiatorMap(Error *error) + { + *error = NoError(); + int error_code = 0; + RadiatorMap *radiator_map = new RadiatorMap(InternalGetRadiatorMap(tuvx_, &error_code)); + if (error_code != 0) + { + *error = Error{ 1, CreateString(MUSICA_ERROR_CATEGORY), CreateString("Failed to create radiator map") }; + return nullptr; + } + return radiator_map; + } + } // namespace musica