diff --git a/CMakeLists.txt b/CMakeLists.txt
index 95d98684d745aaabda13b2e82de34a06125a138d..df461be3016ff213e48b8ba55601d301bb048d0e 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -67,6 +67,7 @@ set(SOURCES_F
     srcF/sfx_esm_param.f90
     srcF/sfx_log.f90
     srcF/sfx_log_param.f90
+    srcF/sfx_run.f90
     srcF/sfx_main.f90
     srcF/sfx_phys_const.f90
     srcF/sfx_surface.f90
diff --git a/srcF/sfx_main.f90 b/srcF/sfx_main.f90
index c96cc2890879301a40c0ea162e94389b33a10dac..0b3599e482ece3413c10958f6addcef5ee09026f 100644
--- a/srcF/sfx_main.f90
+++ b/srcF/sfx_main.f90
@@ -1,161 +1,9 @@
-!> @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 :: io, 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
-
-
-    !< @brief define number of cells
-    open(newunit = io, 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 (io, *, iostat = status) meteo_cell%U, meteo_cell%dT, meteo_cell%Tsemi, meteo_cell%dQ
-        num = num + 1
-    enddo
-    num = num - 1
-
-    close(io)
-
-    ! --- 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(newunit = io, 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(io, *) 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(io)
-
-
-    !< @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 config_parser
-    use iso_c_binding, only: C_NULL_CHAR
-#endif
-
-    use sfx_common
     use sfx_config
+    use sfx_run
     ! --------------------------------------------------------------------------------
 
     ! directives list
@@ -167,241 +15,12 @@ program sfx_main
     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'
-    character(len = 128), parameter :: arg_key_config = "--config"
-
-    integer :: is_output_set
-    ! --------------------------------------------------------------------------------
-
-    ! local variables
-    ! --------------------------------------------------------------------------------
-    integer :: i
-    integer :: status
-    integer :: id, nmax
+    character(len=:), allocatable :: filename_out
     ! --------------------------------------------------------------------------------
 
-#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(*, *) ' --config [filename]'
-            write(*, *) '    use configuration file'
-            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
-            !< save nmax if previously set
-            nmax = dataset%nmax
-            call set_dataset(dataset, id)
-            dataset%nmax = nmax
-
-            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
-
-                !< filename
-                call get_command_argument(i + 2, dataset%filename)
-                
-                !< 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
-
-                !< 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
-                
-                !< 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
-        else if (trim(arg) == trim(arg_key_config)) then
-            if (i == num_args) then
-                write(*, *) ' FAILURE! > missing configuration file [key] argument'
-                stop
-            end if
-            
-            call get_command_argument(i + 1, arg)
-#ifdef USE_CONFIG_PARSER
-            call c_config_run(trim(arg)//C_NULL_CHAR, status)
-            if (status == 0) then
-                write(*, *) ' FAILURE! > unable to parse configuration file: ', trim(arg)
-                stop
-            end if 
-
-            call c_config_is_varname("model.id"//C_NULL_CHAR, status)
-            if (status /= 0) then
-                call c_config_get_string("model.id"//C_NULL_CHAR, config_field, status)
-                ! *: check status
-                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 c_config_is_varname("dataset.id"//C_NULL_CHAR, status)            
-            if (status /= 0) then
-                call c_config_get_string("dataset.id"//C_NULL_CHAR, config_field, status)
-                ! *: check status
-                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
-                !< save nmax if previously set
-                nmax = dataset%nmax
-                call set_dataset(dataset, id)
-                dataset%nmax = nmax
-
-                call c_config_is_varname("dataset.filename"//C_NULL_CHAR, status)
-                if ((status /= 0).or.(dataset%id == dataset_user)) then
-                    !< mandatory in user dataset
-                    call c_config_get_string("dataset.filename"//C_NULL_CHAR, config_field, status)
-                    ! *: check status
-                    dataset%filename = char_array2str(config_field)
-                end if
-
-                call c_config_is_varname("dataset.h"//C_NULL_CHAR, status)
-                if ((status /= 0).or.(dataset%id == dataset_user)) then
-                    !< mandatory in user dataset
-                    call c_config_get_float("dataset.h"//C_NULL_CHAR, dataset%h, status)
-                end if
-
-                call c_config_is_varname("dataset.z0_m"//C_NULL_CHAR, status)
-                if ((status /= 0).or.(dataset%id == dataset_user)) then
-                    !< mandatory in user dataset
-                    call c_config_get_float("dataset.z0_m"//C_NULL_CHAR, dataset%z0_m, status)
-                end if
-
-                call c_config_is_varname("dataset.z0_h"//C_NULL_CHAR, status)
-                if ((status /= 0).or.(dataset%id == dataset_user)) then
-                    !< mandatory in user dataset
-                    call c_config_get_float("dataset.z0_h"//C_NULL_CHAR, dataset%z0_h, status)
-                end if
-            end if
-
-            call c_config_is_varname("dataset.nmax"//C_NULL_CHAR, status)
-            if (status /= 0) then
-                call c_config_get_int("dataset.nmax"//C_NULL_CHAR, dataset%nmax, status)
-            end if 
-        
-            call c_config_is_varname("output.filename"//C_NULL_CHAR, status)
-            if (status /= 0) then
-                call c_config_get_string("output.filename"//C_NULL_CHAR, config_field, status)
-                ! *: check status
-                filename_out = char_array2str(config_field)
-            
-                is_output_set = 1
-            end if
-
-            if (allocated(config_field)) deallocate(config_field)
-#endif
-        end if        
-    end do
 
-    !< @bried set auto-defined output filename 
-    if (is_output_set == 0) then 
-        filename_out = 'output-' // trim(get_dataset_tag(dataset%id)) // '.txt'
-    end if
+    !< @brief run setup
+    call set_run(filename_out, dataset, model)
 
     !< @brief running main driver
     call run_dataset(filename_out, dataset, model)
diff --git a/srcF/sfx_run.f90 b/srcF/sfx_run.f90
new file mode 100644
index 0000000000000000000000000000000000000000..3a2eb33b341789a152e932bfac5369bfe3b91ebe
--- /dev/null
+++ b/srcF/sfx_run.f90
@@ -0,0 +1,426 @@
+module sfx_run
+
+    ! directives list
+    ! --------------------------------------------------------------------------------
+    implicit none
+    private
+    ! --------------------------------------------------------------------------------
+
+    ! public interface
+    ! --------------------------------------------------------------------------------
+    public :: set_run
+    public :: run_dataset
+    ! --------------------------------------------------------------------------------
+
+contains 
+
+    !> @brief main run sfx on dataset subroutine
+    ! ----------------------------------------------------------------------------
+    subroutine run_dataset(filename_out, dataset, model)
+
+        ! modules used
+        ! --------------------------------------------------------------------------------
+        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 :: io, 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
+    
+    
+        !< @brief define number of cells
+        open(newunit = io, 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 (io, *, iostat = status) meteo_cell%U, meteo_cell%dT, meteo_cell%Tsemi, meteo_cell%dQ
+            num = num + 1
+        enddo
+        num = num - 1
+    
+        close(io)
+    
+        ! --- 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(newunit = io, 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(io, *) 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(io)
+    
+    
+        !< @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
+    
+    !> @brief main setup sfx & dataset subroutine
+    ! ----------------------------------------------------------------------------
+    subroutine set_run(filename_out, dataset, model)
+    
+        ! modules used
+        ! --------------------------------------------------------------------------------
+#ifdef USE_CONFIG_PARSER
+        use config_parser
+        use iso_c_binding, only: C_NULL_CHAR
+#endif
+    
+        use sfx_common
+        use sfx_config
+        ! --------------------------------------------------------------------------------
+    
+        ! directives list
+        ! --------------------------------------------------------------------------------
+        implicit none
+        ! --------------------------------------------------------------------------------
+    
+        character(len=:), allocatable, intent(out) :: filename_out
+        type(sfxDatasetType), intent(out) :: dataset
+        integer, intent(out) :: model
+    
+        ! 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'
+        character(len = 128), parameter :: arg_key_config = "--config"
+    
+        integer :: is_output_set
+        ! --------------------------------------------------------------------------------
+    
+        ! local variables
+        ! --------------------------------------------------------------------------------
+        integer :: i
+        integer :: status
+        integer :: id, nmax
+    
+#ifdef USE_CONFIG_PARSER
+        character, allocatable :: config_field(:)
+#endif
+        ! --------------------------------------------------------------------------------
+    
+    
+        !< 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
+    
+        !< 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(*, *) ' --config [filename]'
+                write(*, *) '    use configuration file'
+                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
+                !< save nmax if previously set
+                nmax = dataset%nmax
+                call set_dataset(dataset, id)
+                dataset%nmax = nmax
+    
+                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
+    
+                    !< filename
+                    call get_command_argument(i + 2, dataset%filename)
+                    
+                    !< 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
+    
+                    !< 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
+                    
+                    !< 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, arg)
+                filename_out = trim(arg)
+                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
+            else if (trim(arg) == trim(arg_key_config)) then
+                if (i == num_args) then
+                    write(*, *) ' FAILURE! > missing configuration file [key] argument'
+                    stop
+                end if
+                
+                call get_command_argument(i + 1, arg)
+#ifdef USE_CONFIG_PARSER
+                call c_config_run(trim(arg)//C_NULL_CHAR, status)
+                if (status == 0) then
+                    write(*, *) ' FAILURE! > unable to parse configuration file: ', trim(arg)
+                    stop
+                end if 
+    
+                call c_config_is_varname("model.id"//C_NULL_CHAR, status)
+                if (status /= 0) then
+                    call c_config_get_string("model.id"//C_NULL_CHAR, config_field, status)
+                    ! *: check status
+                    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 c_config_is_varname("dataset.id"//C_NULL_CHAR, status)            
+                if (status /= 0) then
+                    call c_config_get_string("dataset.id"//C_NULL_CHAR, config_field, status)
+                    ! *: check status
+                    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
+                    !< save nmax if previously set
+                    nmax = dataset%nmax
+                    call set_dataset(dataset, id)
+                    dataset%nmax = nmax
+    
+                    call c_config_is_varname("dataset.filename"//C_NULL_CHAR, status)
+                    if ((status /= 0).or.(dataset%id == dataset_user)) then
+                        !< mandatory in user dataset
+                        call c_config_get_string("dataset.filename"//C_NULL_CHAR, config_field, status)
+                        ! *: check status
+                        dataset%filename = char_array2str(config_field)
+                    end if
+    
+                    call c_config_is_varname("dataset.h"//C_NULL_CHAR, status)
+                    if ((status /= 0).or.(dataset%id == dataset_user)) then
+                        !< mandatory in user dataset
+                        call c_config_get_float("dataset.h"//C_NULL_CHAR, dataset%h, status)
+                    end if
+    
+                    call c_config_is_varname("dataset.z0_m"//C_NULL_CHAR, status)
+                    if ((status /= 0).or.(dataset%id == dataset_user)) then
+                        !< mandatory in user dataset
+                        call c_config_get_float("dataset.z0_m"//C_NULL_CHAR, dataset%z0_m, status)
+                    end if
+    
+                    call c_config_is_varname("dataset.z0_h"//C_NULL_CHAR, status)
+                    if ((status /= 0).or.(dataset%id == dataset_user)) then
+                        !< mandatory in user dataset
+                        call c_config_get_float("dataset.z0_h"//C_NULL_CHAR, dataset%z0_h, status)
+                    end if
+                end if
+    
+                call c_config_is_varname("dataset.nmax"//C_NULL_CHAR, status)
+                if (status /= 0) then
+                    call c_config_get_int("dataset.nmax"//C_NULL_CHAR, dataset%nmax, status)
+                end if 
+            
+                call c_config_is_varname("output.filename"//C_NULL_CHAR, status)
+                if (status /= 0) then
+                    call c_config_get_string("output.filename"//C_NULL_CHAR, config_field, status)
+                    ! *: check status
+                    filename_out = char_array2str(config_field)
+                
+                    is_output_set = 1
+                end if
+    
+                if (allocated(config_field)) deallocate(config_field)
+#endif
+            end if        
+        end do
+    
+        !< set auto-defined output filename 
+        if (is_output_set == 0) then 
+            filename_out = 'output-' // trim(get_dataset_tag(dataset%id)) // '.txt'
+        end if
+    
+    end subroutine
+    
+
+end module sfx_run