sfx_main.f90

!> @brief main run sfx on dataset subroutine
! ----------------------------------------------------------------------------
subroutine run_dataset(filename_out, dataset, model)

    use sfx_phys_const
    use sfx_common
    use sfx_config
    use sfx_io
    use sfx_data

    use sfx_esm, only: &
            get_surface_fluxes_vec_esm => get_surface_fluxes_vec, &
            numericsType_esm => numericsType
    use sfx_log, only: &
            get_surface_fluxes_vec_log => get_surface_fluxes_vec, &
            numericsType_log => numericsType
    use sfx_most, only: &
            get_surface_fluxes_vec_most => get_surface_fluxes_vec, &
            numericsType_most => numericsType
    use sfx_sheba, only: &
            get_surface_fluxes_vec_sheba => get_surface_fluxes_vec, &
            numericsType_sheba => numericsType
    ! --------------------------------------------------------------------------------

    ! directives list
    ! --------------------------------------------------------------------------------
    implicit none
    ! --------------------------------------------------------------------------------

    character(len=*), intent(in) :: filename_out
    type(sfxDatasetType), intent(in) :: dataset
    integer, intent(in) :: model


    ! input/output model data
    ! --------------------------------------------------------------------------------
    type(meteoDataVecType) :: meteo         !< meteorological data (input)
    type(meteoDataType) :: meteo_cell

    type(sfxDataVecType) :: sfx             !< surface fluxes (output)

    type(numericsType_esm) :: numerics_esm      !< surface flux module (ESM) numerics parameters
    type(numericsType_log) :: numerics_log      !< surface flux module (LOG) numerics parameters
    type(numericsType_most) :: numerics_most    !< surface flux module (MOST) numerics parameters
    type(numericsType_sheba) :: numerics_sheba  !< surface flux module (SHEBA) numerics parameters

    integer :: num                          !< number of 'cells' in input
    ! --------------------------------------------------------------------------------

    ! local variables
    ! --------------------------------------------------------------------------------
    integer :: i
    integer :: status
    ! --------------------------------------------------------------------------------


    write(*, *) ' Running SFX:'
    !write(*, '(a,a)') '    model = ', trim(get_model_tag(model))
    !write(*, '(a,a)') '    dataset = ', trim(get_dataset_tag(dataset%id))
    !write(*, '(a,a)') '    filename[IN] = ', trim(dataset%filename)
    !write(*, '(a,a)') '    filename[OUT] = ', trim(filename_out)
    !write(*, '(a, g0)') '    surface type = ', dataset%surface_type
    !write(*, '(a, g0)') '    h = ', dataset%h
    !write(*, '(a, g0)') '    z0(m) = ', dataset%z0_m
    !write(*, '(a, g0)') '    z0(h) = ', dataset%z0_h

    write(*, *) '    model = ', trim(get_model_tag(model))
    write(*, *) '    dataset = ', trim(get_dataset_tag(dataset%id))
    write(*, *) '    filename[IN] = ', trim(dataset%filename)
    write(*, *) '    filename[OUT] = ', trim(filename_out)
    write(*, *) '    surface type = ', dataset%surface_type
    write(*, *) '    h = ', dataset%h
    write(*, *) '    z0(m) = ', dataset%z0_m
    write(*, *) '    z0(h) = ', dataset%z0_h


    !< @brief define number of cells
    open(32, file = dataset%filename, iostat = status, status ='old')
    if (status /= 0) then
        write(*, *) ' FAILURE! > unable to open file: ', trim(dataset%filename)
        return
    end if

    num = 0
    status = 0
    do while (status.eq.0)
        read (32, *, iostat = status) meteo_cell%U, meteo_cell%dT, meteo_cell%Tsemi, meteo_cell%dQ
        num = num + 1
    enddo
    num = num - 1

    close(32)

    ! --- print number of elements in dataset
    write(*, '(a, g0)') '    size = ', num
    if (dataset%nmax > 0) then
        write(*, '(a, g0)') '    nmax = ', dataset%nmax
        num = min(num, dataset%nmax)
    end if


    !< @brief allocate input & output data
    call allocate_meteo_vec(meteo, num)
    call allocate_sfx_vec(sfx, num)


    !< @brief setting height & roughness
    meteo_cell%h = dataset%h
    meteo_cell%z0_m = dataset%z0_m

    !< @brief read input data
    open(32, file = dataset%filename, iostat = status, status = 'old')
    if (status /= 0) then
        write(*, *) ' FAILURE! > unable to open file: ', trim(dataset%filename)
        return
    end if
    do i = 1, num
        read(32, *) meteo_cell%U,  meteo_cell%dT,  meteo_cell%Tsemi,  meteo_cell%dQ

        meteo%h(i) = meteo_cell%h
        meteo%U(i) = meteo_cell%U
        meteo%dT(i) = meteo_cell%dT
        meteo%Tsemi(i) = meteo_cell%Tsemi
        meteo%dQ(i) = meteo_cell%dQ
        meteo%z0_m(i) = meteo_cell%z0_m
    enddo
    close(32)


    !< @brief calling flux module
    if (model == model_esm) then
        call get_surface_fluxes_vec_esm(sfx, meteo, numerics_esm, num)
    else if (model == model_log) then
        call get_surface_fluxes_vec_log(sfx, meteo, numerics_log, num)
    else if (model == model_most) then
        call get_surface_fluxes_vec_most(sfx, meteo, numerics_most, num)
    else if (model == model_sheba) then
        call get_surface_fluxes_vec_sheba(sfx, meteo, numerics_sheba, num)
    end if


    !< @brief write output data
    call write_ascii_vec11(filename_out, &
        sfx%zeta, sfx%Rib, &
        sfx%Re, sfx%B, sfx%z0_m, sfx%z0_t, &
        sfx%Rib_conv_lim, &
        sfx%Cm,sfx%Ct, sfx%Km, sfx%Pr_t_inv, num, '(11(f10.4,3x))', status)


    !< @brief deallocate input & output data
    call deallocate_meteo_vec(meteo)
    call deallocate_sfx_vec(sfx)

end subroutine


program sfx_main

    ! modules used
    ! --------------------------------------------------------------------------------
#ifdef USE_CONFIG_PARSER
    use parser_sub_f
    use parser
    use sfx_surface
    use iso_c_binding, only: C_NULL_CHAR
#endif

    use sfx_common
    use sfx_config
    ! --------------------------------------------------------------------------------

    ! directives list
    ! --------------------------------------------------------------------------------
    implicit none
    ! --------------------------------------------------------------------------------


    type(sfxDatasetType) :: dataset
    integer :: model

    character(len = 256) :: filename_out

    ! command line arguments
    ! --------------------------------------------------------------------------------
    integer :: num_args
    character(len = 128) :: arg

    character(len = 128), parameter :: arg_key_model = '--model'
    character(len = 128), parameter :: arg_key_dataset = '--dataset'
    character(len = 128), parameter :: arg_key_output = '--output'
    character(len = 128), parameter :: arg_key_nmax = '--nmax'
    character(len = 128), parameter :: arg_key_help = '--help'

    integer :: is_output_set
    ! --------------------------------------------------------------------------------

    ! local variables
    ! --------------------------------------------------------------------------------
    integer :: i
    integer :: status
    integer :: id
    ! --------------------------------------------------------------------------------

#ifdef USE_CONFIG_PARSER
    character, allocatable :: config_field(:)
#endif


    !< @brief define default model & dataset
    model = model_esm                           !< default = ESM
    call set_dataset(dataset, dataset_mosaic)   !< default = MOSAiC
    is_output_set = 0                           !< default = auto define output filename

    !< @brief command line arguments processing
    num_args = command_argument_count()
    do i = 1, num_args
        call get_command_argument(i, arg)
        if (trim(arg) == trim(arg_key_help)) then
            write(*, *) ' sfx model, usage:'
            write(*, *) ' --help '
            write(*, *) '    print usage options '
            write(*, *) ' --model [key]'
            write(*, *) '    key = esm (default) || log || most || sheba'
            write(*, *) ' --dataset [key]'
            write(*, *) '    key = mosaic (default) || irgason || sheba'
            write(*, *) '        = lake || papa || toga || user [filename] [param]'
            write(*, *) '    param = [h] [z0(m)] [z0(h)]'
            write(*, *) ' --output [filename]'
            write(*, *) '    set output filename '
            write(*, *) ' --nmax [value]'
            write(*, *) '    max number of data points > 0 '
            stop
        end if
        if (trim(arg) == trim(arg_key_model)) then
            if (i == num_args) then
                write(*, *) ' FAILURE! > missing model [key] argument'
                stop
            end if
            
            call get_command_argument(i + 1, arg)
            model = get_model_id(arg)
            if (model == -1) then
                write(*, *) ' FAILURE! > unknown model [key]: ', trim(arg)
                stop
            end if
        else if (trim(arg) == trim(arg_key_dataset)) then
            if (i == num_args) then
                write(*, *) ' FAILURE! > missing dataset [key] argument'
                stop
            end if
            
            call get_command_argument(i + 1, arg)
            id = get_dataset_id(arg)
            if (id == -1) then
                write(*, *) ' FAILURE! > unknown dataset [key]: ', trim(arg)
                stop
            end if
            call set_dataset(dataset, id)

            if (dataset%id == dataset_user) then
                !< @brief user-defined dataset
                if (i + 5 > num_args) then
                    write(*, *) ' FAILURE! > incorrect arguments for [user] dataset'
                    stop
                end if

                !< @brief filename
                call get_command_argument(i + 2, dataset%filename)
                
                !< @brief reading 'h'
                call get_command_argument(i + 3, arg)
                call str2real(dataset%h, arg, status)
                if (status /= 0) then
                    write(*, *) ' FAILURE! > expecting real h [value]'
                    stop
                end if
                if (dataset%h <= 0) then
                    write(*, *) ' FAILURE! > h [value] should be positive'
                    stop
                end if

                !< @brief reading 'z0(m)'
                call get_command_argument(i + 4, arg)
                call str2real(dataset%z0_m, arg, status)
                if (status /= 0) then
                    write(*, *) ' FAILURE! > expecting real z0(m) [value]'
                    stop
                end if
                
                !< @brief reading 'z0(h)'
                call get_command_argument(i + 5, arg)
                call str2real(dataset%z0_h, arg, status)
                if (status /= 0) then
                    write(*, *) ' FAILURE! > expecting real z0(h) [value]'
                    stop
                end if
                
            end if 
        else if (trim(arg) == trim(arg_key_output)) then
            if (i == num_args) then
                write(*, *) ' FAILURE! > missing output [key] argument'
                stop
            end if
            call get_command_argument(i + 1, filename_out)
            is_output_set = 1
        else if (trim(arg) == trim(arg_key_nmax)) then
            if (i == num_args) then
                write(*, *) ' FAILURE! > missing nmax [key] argument'
                stop
            end if
            call get_command_argument(i + 1, arg)
            call str2int(dataset%nmax, arg, status)
            if (status /= 0) then
                write(*, *) ' FAILURE! > expecting int nmax [value]'
                stop
            end if
            if (dataset%nmax <= 0) then
                write(*, *) ' FAILURE! > nmax [value] should be positive'
                stop
            end if
        end if
    end do

#ifdef USE_CONFIG_PARSER
    call run("config.txt"//C_NULL_CHAR)

    call is_varname("model.id"//C_NULL_CHAR, status)
    if (status /= 0) then
        call get_charf("model.id"//C_NULL_CHAR, config_field)
        model = get_model_id(char_array2str(config_field))
        if (model == -1) then
            write(*, *) ' FAILURE! > unknown model [key]: ', trim(char_array2str(config_field))
            stop
        end if
    end if 

    call is_varname("dataset.id"//C_NULL_CHAR, status)
    if (status /= 0) then
        call get_charf("dataset.id"//C_NULL_CHAR, config_field)
        id = get_dataset_id(char_array2str(config_field))
        if (id == -1) then
            write(*, *) ' FAILURE! > unknown dataset [key]: ', trim(char_array2str(config_field))
            stop
        end if
        call set_dataset(dataset, id)

        if (dataset%id == dataset_user) then
            call get_charf("dataset.filename"//C_NULL_CHAR, config_field)
            dataset%filename = char_array2str(config_field)
    
            call get_float("dataset.h"//C_NULL_CHAR, dataset%h)
            call get_float("dataset.z0_m"//C_NULL_CHAR, dataset%z0_m)
            call get_float("dataset.z0_h"//C_NULL_CHAR, dataset%z0_h)
        end if
    end if 

    call is_varname("dataset.nmax"//C_NULL_CHAR, status)
    if (status /= 0) then
        call get_int("dataset.nmax"//C_NULL_CHAR, dataset%nmax)
    end if 

    call is_varname("output.filename"//C_NULL_CHAR, status)
    if (status /= 0) then
        call get_charf("output.filename"//C_NULL_CHAR, config_field)
        filename_out = char_array2str(config_field)

        is_output_set = 1
    end if

    deallocate(config_field)
#endif

    !< @bried set auto-defined output filename 
    if (is_output_set == 0) then 
        filename_out = 'output-' // trim(get_dataset_tag(dataset%id)) // '.txt'
    end if

    open(32, file = "data/MOSAiC_zh.txt", iostat = status, status = 'old')
    if (status /= 0) then
        stop
    end if
    read(32, *) dataset%h, dataset%z0_m
    close(32)

    !< @brief running main driver
    call run_dataset(filename_out, dataset, model)
    stop

end program